Input-output behaviour of a model neuron with alternating drift
aa r X i v : . [ q - b i o . N C ] J a n Input-output behaviour of a model neuronwith alternating drift
A.BUONOCORE (1)
A. DI CRESCENZO (2) ∗ , E. DI NARDO (3) ,(1) Dipartimento di Matematica e Applicazioni, Universit`a di Napoli Federico IIVia Cintia, I-80126 Napoli, ItalyE-mail: [email protected] (2) Dipartimento di Matematica e Informatica, Universit`a di SalernoVia Ponte don Melillo, I-84084 Fisciano (SA), ItalyE-mail: [email protected] (3) Dipartimento di Matematica, Universit`a degli Studi della BasilicataC.da Macchia Romana, I-85100 Potenza, ItalyE-mail: [email protected] Abstract
The input-output behaviour of the Wiener neuronal model subject to alternating inputis studied under the assumption that the effect of such an input is to make the drift itselfof an alternating type. Firing densities and related statistics are obtained via simulations ofthe sample-paths of the process in the following three cases: the drift changes occur duringrandom periods characterized by (i) exponential distribution, (ii)
Erlang distribution with apreassigned shape parameter, and (iii) deterministic distribution. The obtained results arecompared with those holding for the Wiener neuronal model subject to sinusoidal input.
Keywords:
Wiener neuronal model; firing densities; alternating drift.
During the last four decades numerous efforts have been devoted to the construction of mathe-matical models aiming at the description of single neuron’s firing processes. A customary featureof widespread existing models is the assumption that the neurons’ firing is described by the firstpassage of the membrane potential through a threshold, the membrane potential being viewed asa continuous stochastic process. Great attention has being put on the identification of suitableMarkov processes and firing thresholds, especially with a view to include in the model certainrelevant features, such as the effects of external input on the neuron’s response (see Ricciardi,1995, and Ricciardi et al. , 1999, for a description of neuronal models and methods to face thefirst-passage-time problem). ∗ corresponding author (i) exponential distribution, (ii) Erlang distribution with shape parameter 2, and (iii) de-terministic distribution. The computational results obtained for the considered model are thencompared with those holding for the Wiener neuronal model subject to sinusoidal input.Let us now point out some specific features of this paper with a view to relate them toother works in the field. First of all, as mentioned in L´ansk´y et al. (2001), the interpretation ofthe neuronal firing times in terms of alternating periods appears to be natural for some typesof neurons. It should be also mentioned that a similar problem was analyzed in Bulsara et al.(1994), where the first passage time problem of a Wiener process with sinusoidal time-varying driftthrough a constant threshold is treated. Our proposed model, however, offers the advantagesthat the epochs of high and low stimulation can be of different lengths; moreover, due to itsgreater flexibility, it allows to describe time-varying inputs characterized by more general shapesthan sinusoidal ones.
As customary, we assume that the neuron’s membrane potential is described by a continuous one-dimensional stochastic process { X ( t ); t ≥ } representing the changes in the membrane potentialbetween two consecutive neuronal firings, the firing threshold being a constant β . Assuming that X (0) = x < β , the FPT random variable T = inf { t ≥ X ( t ) > β } describes the neuronal interspike intervals. An essential problem is to determine the p.d.f. of T ,i.e. the firing density g ( β, t | x ) = ∂∂t Pr( T ≤ t ) , t > . A first well-known neuronal model was proposed by Gerstein and Mandelbrot (1964), whoassumed that the neuron’s membrane potential undergoes a simple random walk under the effectof excitatory and inhibitory synaptic actions. This assumption, upon a suitable limit procedure,leads to the Wiener process X ( t ) = x + ξ t + σ B ( t ) , t > , with x ∈ R , ξ ∈ R and σ >
0, and where B ( t ) denotes the standard Brownian motion. Eventhough the assumptions of this model are oversimplified as they do not take into account certainproperties such as the spontaneous exponential decay of the membrane potential in the absence ofinput, the Wiener model can be usefully treated as a starting point for investigations on neuronalstochastic models. For instance, the Wiener model with random jumps has been consideredrecently by Giraudo et al. (2001) as a model for the description of changes in the membranepotential depending on the distance between the trigger zone and the synaptic ending.2 Wiener model with alternating drift
Our present aim is to analyse the input-output behaviour of the Wiener model subject to alter-nating input. To such a purpose, we assume that the effects of alternating stimuli reflect into thedrift of the Wiener process. This leads to a Wiener process with alternating drift, described by X ( t ) = x + Z t ξ ( s ) d s + σ B ( t ) , t ≥ , (1)where x ∈ R , σ > { ξ ( t ); t ≥ } is a dichotomous stochastic process on states c and − v ,with c, v > ξ (0) = c . Eq. (1) describes a Brownian motion that starts from x at time0, with positive drift c and infinitesimal variance σ . Denoting by Y , Y , . . . the independentrandom times separating consecutive changes of the drift value, we have ξ ( t ) = c during therandom periods of lengths Y , Y , Y , . . . and ξ ( t ) = − v during the remaining random periods.Formally, ξ ( t ) = c − v c + v − N ( t ) , t ≥ , where N ( t ) = ∞ X k =1 { Y + ··· + Y k ≤ t } , t ≥ . At the random times Y + · · · + Y k the value of ξ ( t ) thus changes from c to − v if k is odd, andfrom − v to c if k is even. We also assume that the non-negative random variables Y , Y , Y , . . . and Y , Y , Y , . . . have distribution functions F U ( t ) and F D ( t ), respectively, so that the randomperiods during which the drift is positive (resp. negative) are i.i.d. We remark that the probabilitylaw of (1) can be expressed as a mixture of Gaussian densities (see Di Crescenzo, 2000). In this Section we shall determine the firing density of the Wiener neuronal model with alternatingdrift in the presence of a constant firing threshold β . Making use of a direct analysis of the sample-paths of process (1), it is possible to obtain a series expression of the firing density. Unfortunately,the series involves integrals of progressively increasing dimensionality, so that such an analyticresult is unsuitable for practical purposes. Nevertheless, by truncation to the first few terms thefollowing lower bound for the firing density is obtained: For all t > β > x we have g ( β, t | x ) ≥ [1 − F U ( t )] e g c ( β, t | x ) + Z t Z β −∞ [1 − F D ( t − u )] e g − v ( β, t − u | x ) α ( x, u | x ) d x d F U ( u ) , where e g η ( β, t | x ) is the first-passage-time density of a Wiener process with drift η throught β , e g η ( β, t | x ) = β − x √ πσ t exp ( − ( β − x − ηt ) σ t ) , while α ( x, t | x ) is the β -avoiding density of a Wiener process with drift c , α ( x, t | x ) = 1 √ πσ t exp ( − ( x − x − c t ) σ t ) (cid:20) − exp (cid:26) − σ t ( β − x ) ( β − x ) (cid:27)(cid:21) . In order to obtain histograms as estimates of the firing densities, we resort to simulationsof the sample paths of process (1). This is performed by making use of a procedure that is3trictly based on the properties of the underlying Wiener process. Such procedure simulates thesample-paths of process (1) at the random times where the drift alternates.A sketch of the procedure for the simulation of T is given below: Begin Procedure x := x , t := 0 , ξ := c ; generate an inversion instant τ ; p := Pr { the first passage occurs before τ } ; generate an uniform pseudo-random number u in (0 , ; if u > p then goto Step 7 ; generate a pseudo-random number θ in (0 , τ ) from p.d.f. f ( θ ) ; f pt := θ + t ; output( f pt ); stop; generate a pseudo-random number z in ( −∞ , β ) from p.d.f. f ( z ) ; x := z , t := t + τ ; update ξ ; goto Step 1 ;End Procedure
For
Step 1 we assume that the random inversion times having distributions F U and F D can benumerically simulated.Let us informally discuss the underlying ideas of this procedure. After the first drift inversiontime τ has been generated, a Bernoulli trial with success probability p is simulated ( Step 2 ),where p is the FPT probability p := P ( T ≤ τ ) = Φ (cid:18) − β − x − ξτσ √ τ (cid:19) + e ξ ( β − x ) /σ Φ (cid:18) − β − x + ξτσ √ τ (cid:19) , Φ denoting the standard normal distribution function. A success of the Bernoulli trial means thatthe first passage has occurred before time τ . In such a case the FPT is simulated via f , which isthe FPT density of a Wiener process through β conditional on { T ≤ τ } , and the simulation ends.If a failure occurs in the Bernoulli trial, then the value attained by the Wiener process at time τ is simulated via density f , which is the density of X ( τ ) conditional on { T > τ } . Densities f and f , whose analytic expressions are known, are simulated by making use of a standard VonNeumann acceptance-rejection method (see for instance Rubinstein, 1981). After X ( τ ) has beengenerated, the current drift η is updated to the new value, i.e. − v if it was c , and vice-versa. Thesimulation than proceeds, restarting afresh by generating a new inversion time of the drift, andso on until the first passage through β occurs.We point out that, according to the assumptions of our model, the case of endogenous pe-riodicity is developed in the given procedure. In other words, the input is reset after that thethreshold is reached. However, we point out that the case of exhogenous periodicity could beconsidered upon an easy modification of the procedure (useful references on the effects of en-dogenous and exhogenous periodicity are the papers by L´ansk´y (1997) and Plesser and Geisel(2001)).In the following sections we discuss the results obtained in three special cases for the driftinversion times: (i) exponentially distributed times, (ii) Erlang distributed times, and (iii) deter-ministic times. The shown estimates of the firing density g ( β, t | x ) are obtained via simulationof 10 sample-paths of X ( t ). 4
20 40 60 80 10000.010.020.030.040.050.060.070.080.090.1
Figure 1:
Estimated firing density of the Wiener model with exponentially alternating drift for σ = 1, λ = 0 . c = v = 1, with the choices µ = 0 . , . , . , . Let us assume that the inversion times of the Wiener neuronal model (1) have distributionfunctions F U ( t ) = 1 − e − λt and F D ( t ) = 1 − e − µt , t ≥
0. Thus, λ − is the mean of the randomperiods during which the excitatory stimuli prevail on inhibitory stimuli, while µ − is the meanof the random periods with opposite behaviour. We have considered some cases of interest in thephysiological contexts, with µ ≥ λ , with threshold β = 10, and with x = 0.The most evident feature of the estimated firing densities is their being unimodal (see Figure1). In particular, if µ increases the mode m increases while the median q decreases, as well asthe quartiles q and q , and the inter-quartile range IQR . This is in agreement with the factthat if µ increases, the random periods with prevailing inhibitory stimuli becomes ‘smaller’ instochastic sense. A quite similar behaviour is observed if λ decreases (see Figure 2). Moreover,as σ increases the firing density becomes less peaked, and q , q , q and m decrease (see Figure3). Finally, when c and v increase and c = v , then q , q , q and m decrease (see Figure 4), as itshould be expected since for high values of c the sample-paths of X ( t ) are more likely to reachthe firing threshold. In this Section we consider the case when the inversion times of the Wiener neuronal model (1)have distribution functions of Erlang type with shape parameter 2: F U ( t ) = 1 − λt e − λt and F D ( t ) = 1 − µt e − µt , t ≥
0. The mean values of the times intervals separating consecutiveinversions of the drift are now given by 2 /λ and 2 /µ . The analysis has been performed on somecases analogous to those presented in Section 5, with µ ≥ λ , with threshold β = 10, and with x = 0.As for the case of exponential inversion times, the estimated firing densities are unimodal(see Figures 5 ÷ µ increases the median q decreases, as well as the quartiles q and q , and the inter-quartile range IQR (see Figure 5). A similar behaviour is again observed5
20 40 60 80 10000.010.020.030.040.050.060.070.080.090.1
Figure 2:
As in Figure 1, with σ = 1, µ = 0 . c = v = 1, for λ = 0 . , . , . , .
05 (bottom to top nearthe origin).
Figure 3:
As in Figure 1, with λ = 0 . µ = 0 . c = v = 1, for σ = 1 . , , .
20 40 60 80 10000.010.020.030.040.050.060.070.080.090.1
Figure 4:
As in Figure 1, with λ = 0 . µ = 0 . σ = 1, for c = v = 0 . , , . if λ decreases (see Figure 6). We emphasize that also in this case the firing density becomesless peaked as σ increases (see Figure 7). Finally, when c and v increase, with c = v , then q , q , q decrease (see Figure 8). A similar behaviour has also been exhibited by the estimatedfiring densities in the case of exponentially distributed inversion times (compare Figures 5 ÷ ÷ With reference to the model (1) let us assume that c = v and that the drift inversion timesare deterministic, with F U ( t ) = F D ( t ) = { t ≥ τ } . This means that the second term at the right-hand-side of (1) increases and decreases alternately in a deterministic fashion every τ instants.In other words, u ( t ) := Z t ξ ( s ) d s (2)for t > τ .As in the previous cases, we have performed 10 simulations to obtain estimates of the firingdensity g ( β, t | x ), with threshold β = 10 and initial state x = 0.Firing densities now appear to be quite different from the previous cases. Indeed, the esti-mated firing densities are multimodal. In Figure 9 two examples are shown where the peaks arelocated near the maxima of (2), i.e. near τ, τ, τ, . . . . In this Section we shall analyse the effect of a smoother input for the Wiener neuronal modelof Section 7 by changing function (2) to a smoother one. More specifically, we assume that themembrane depolarization is now described by the Gauss-Markov diffusion process { ˆ X ( t ); t ≥ }
10 20 30 40 5000.020.040.060.080.10.120.140.16
Figure 5:
Estimated firing density of the Wiener model with Erlang alternating drift for σ = 1, λ = 0 . c = v = 1, for µ = 0 . , . , . , .
25 (bottom to top near abscissa 20).
Figure 6:
As in Figure 5, with σ = 1, µ = 0 . c = v = 1, for λ = 0 . , . , . , .
025 (bottom to topnear the origin).
10 20 30 40 5000.050.10.150.20.250.30.350.40.45
Figure 7:
As in Figure 5, with λ = 0 . µ = 0 . c = v = 1, for σ = 1 . , . , . , . , . Figure 8:
As in Figure 5, with λ = 0 . µ = 0 . σ = 1, for c = v = 0 . , , . Figure 9:
Estimated firing density of the Wiener model with drift alternating at deterministic times for τ = 5, c = v = 2, σ = 1 (lower peaks) and σ = 0 . given by ˆ X ( t ) = x + Z t η ( s ) d s + σ B ( t ) , t ≥ , (3)with x ∈ R , σ >
0, and η ( t ) = α π τ sin (cid:18) πτ t (cid:19) , with α and τ positive. Note that mean and covariance of ˆ X ( t ) are given by E [ ˆ X ( t )] = α (cid:20) − cos (cid:18) πτ t (cid:19)(cid:21) , t ≥ , (4)and Cov[ ˆ X ( s ) , ˆ X ( t )] = σ min { s, t } , s, t ≥ , respectively. Hence, α is the maximum value attained by the mean of ˆ X ( t ). To compare thepresent model with that considered in the previous section, we choose the involved parameters insuch a way that the functions (2) and (4) possess equal maxima and minima. Hence, accordingto the parameters’ values chosen in Figure 9 we set α = c τ , with c = v , β = 10 and τ = 5.To determine the firing density for the model (3) we resort to a numerical procedure due toDi Nardo et al. (2001). The results obtained show that the firing densities are still multimodal(see Figure 10), as for the model described in Section 7, and that the peaks are located aroundthe maxima of (4). The shapes of the firing densities are similar to those obtained for the Wienermodel with alternating drift, though being now characterized by less sharp peaks. This is alsoevident by comparing Figure 10 with Figure 9. The Wiener neuronal model has been considered by assuming that the effects of alternatingstimuli are included into the drift. Making use of a numerical procedure to simulate sample-10
Figure 10:
The firing density of the Wiener model with sinusoidal drift for τ = 5, α = 10 σ = 1 (lowerpeaks) and σ = 0 . paths of the resulting process, we have obtained estimates of the firing densities as FPT densitiesthrough a constant threshold. Three cases have been considered: the time periods separatingconsecutive changes of the drift (i) are exponentially distributed, (ii) they have an Erlang-type distribution, and (iii) they possess constant length. For sensible choices of the involvedparameters, the results obtained in these three cases are quite different. Indeed, while thepresence of randomness in the drift of the process produces unimodal firing densities in cases (i) and (ii) , the deterministic behaviour of the drift yields multimodal densities in case (iii) . Thecomputational results found in case (iii) have been finally compared with those holding for theWiener neuronal model subject to oscillating input of sinusoidal type. In this case the firingdensities are still multimodal, though less peaked than in case (iii) . Acknowledgements
This work has been performed within a joint cooperation agreement between Japan Science andTechnology Coorporation (JST) and Universit`a di Napoli Federico II, under partial support byMIUR (PRIN 2000).
References [1] Bulsara A.R., Lowen S.B., Rees C.D., 1994. Cooperative behavior in the periodically mod-ulated Wiener process: Noise-induced complexity in a model neuron. Phys. Rev. E, 49,4989–5000.[2] Di Crescenzo A., 2000. On Brownian motions with alternating drifts. In: Cybernetics andSystems 2000 (Trappl R ed.), 324–329. Austrian Society for Cybernetics Studies, Vienna.[3] Di Nardo E., Nobile A.G., Pirozzi E., Ricciardi L.M., 2001. A computational approach tofirst-passage-time problems for Gauss-Markov processes. Adv. Appl. Prob., 33, 453–482.114] Gerstein G.L., Mandelbrot B., 1964. Random walk models for the spike activity of a singleneuron. Biophys. J., 4, 41–68.[5] Giraudo M.T., Sacerdote L., Sirovich R., 2001. Effects of random jumps on a very simpleneuronal diffusion model. In: Abstract book of the 4th International Workshop “NeuronalCoding 2001”, 145. 10-14 September 2001, Plymouth, UK.[6] L´ansk´y P., 1997. Sources of periodical force in noisy integrate-and-fire models of neuronaldynamics. Phys. Rev. E, 55, 2040–2043.[7] L´ansk´y P., Krivan V., Rospars J.P., 2001. Ligand-receptor interaction under periodic stim-ulation: a modeling study of concentration chemoreceptors. European Biophys. J., 30, 110–120.[8] Plesser H.E., Geisel T., 2001. Stochastic resonance in neuron models: Endogenous stimula-tion revisited. Phys. Rev. E, 63, 31916–31921.[9] Ricciardi L.M., 1995. Diffusion models of neuron activity. In: The Handbook of Brain Theoryand Neural Networks (Arbib MA ed.), 299–304. The MIT Press, Cambridge.[10] Ricciardi L.M., Di Crescenzo A., Giorno V., Nobile A.G., 1999. An outline of theoreticaland algorithmic approaches to first passage time problems with applications to biologicalmodeling. Math. Japonica, 50, 247–322.[11] Rubinstein R.Y., 1981.