Analysis of pacemaker activity in a two-component model of some brainstem neurons Henry C. Tuckwell{}^{1,2*}, Ying Zhou{}^{3,}, Nicholas J. Penington {}^{4,5}   {}^{1} School of Electrical and Electronic Engineering, University of Adelaide,Adelaide, South Australia 5005, Australia {}^{2} School of Mathematical Sciences, Monash University, Clayton, Victoria 3800, Australia 3 Department of Mathematics, Lafayette College, 1 Pardee Drive, Easton, PA 18042, USA {}^{4} Department of Physiology and Pharmacology,{}^{5} Program in Neural and Behavioral Science and Robert F. Furchgott Center for Neural and Behavioral Science State University of New York, Downstate Medical Center,Box 29, 450 Clarkson Avenue, Brooklyn, NY 11203-2098, USA {}^{*} Corresponding author: Henry C. Tuckwell, School of Electrical and Electronic Engineering, University of Adelaide, North Terrace, Adelaide, SA 5005, Australia;  Tel: 61-481192816; Fax: 61-8-83134360: email:  

Analysis of pacemaker activity in a two-component model of some brainstem neurons
Henry C. Tuckwell, Ying Zhou, Nicholas J. Penington
School of Electrical and Electronic Engineering, University of Adelaide,
Adelaide, South Australia 5005, Australia
School of Mathematical Sciences, Monash University, Clayton, Victoria 3800, Australia
Department of Mathematics, Lafayette College, 1 Pardee Drive, Easton, PA 18042, USA
Department of Physiology and Pharmacology,
Program in Neural and Behavioral Science and Robert F. Furchgott Center for Neural and Behavioral Science
State University of New York, Downstate Medical Center,
Box 29, 450 Clarkson Avenue, Brooklyn, NY 11203-2098, USA
Corresponding author: Henry C. Tuckwell, School of Electrical and Electronic Engineering, University of Adelaide, North Terrace, Adelaide, SA 5005, Australia;  
Tel: 61-481192816; Fax: 61-8-83134360: email:


Serotonergic, noradrenergic and dopaminergic brainstem (including midbrain) neurons, often exhibit spontaneous and fairly regular spiking with frequencies of order a few Hz, though dopaminergic and noradrenergic neurons only exhibit such pacemaker-type activity in vitro or in vivo under special conditions. A large number of ion channel types contribute to such spiking so that detailed modeling of spike generation leads to the requirement of solving very large systems of differential equations. It is useful to have simplified mathematical models of spiking in such neurons so that, for example, features of inputs and output spike trains can be incorporated including stochastic effects for possible use in network models. In this article we investigate a simple two-component conductance-based model of the Hodgkin-Huxley type. Solutions are computed numerically and with suitably chosen parameters mimic features of pacemaker-type spiking in the above types of neurons. The effects of varying parameters is investigated in detail, it being found that there is extreme sensitivity to eight of them. Transitions from non-spiking to spiking are examined for two of these, the half-activation potential for an activation variable and the added (depolarizing) current and contrasted with the behavior of the classical Hodgkin-Huxley system. The plateaux levels between spikes can be adjusted, by changing a set of voltage parameters, to agree with experimental observations. Experiment has shown that in, in vivo, dopaminergic and noradrenergic neurons’ pacemaker activity can be induced by the removal of excitatory inputs or the introduction of inhibitory ones. These properties are confirmed by mimicking opposite such changes in the model, which resulted in a change from pacemaker activity to bursting-type phenomena. The article concludes with a brief review of previous modeling of these types of neurons and a brief discussion of their involvement in some pathological pysychiatric conditions.

Keywords: Dorsal raphe nucleus, serotonergic neurons, locus coeruleus, noradrenergic neurons, dopaminergic neurons, reduced computational models, pacemaker activity.



5-HT, 5-hydroxytryptamine (serotonin); AC, adenylate cyclase; AHP, afterhyperpolarization; cAMP, cyclic adenosine monophosphate; CREB, cAMP response element binding protein; CRH, corticotropin-releasing hormone; CRN, caudal raphe nucleus; DA, dopamine or dopaminergic; DRN, dorsal raphe nucleus; GABA, gamma-aminobutyric acid; HH, Hodgkin-Huxley; Hz, hertz; HPA, hypothalamus-pituitary-adrenal; ISI, interspike interval; LC, locus coeruleus; MDD, major depressive disorder; NA, noradrenaline or noradrenergic; OCD, obsessive-compulsive disorder; PTSD, post-traumatic stress disorder; PVN, paraventricular nucleus (of hypothalamus); REM, rapid eye movement; SE, serotonin or serotonergic.

1 Introduction

Neurons which exhibit (approximately) periodic spiking in the presumed absence of synaptic input are called autonomous pacemakers, and include neurons found in the subthalamic nucleus, nucleus basalis, globus pallidus, raphe nuclei, cerebellum, locus ceruleus, ventral tegmental area, and substantia nigra (Ramirez et al., 2011). In order to sustain spiking, some pacemaking cells may require small amounts of depolarizing inputs, natural or laboratory. Thus, for example Pan et al. (1994) found that all of 42 rat LC neurons fired spontaneously, whereas in Williams et al. (1982) and Ishimatsu et al. (1996) it was reported that most cells did not require excitatory input to fire regularly. These latter three sets of results were obtained in vitro. Sanchez-Padilla et al. (2014) reported that spike rate in mouse LC neurons was not affected by blockers of glutamatergic or GABAergic synaptic input, supporting the idea that these cells were autonomous pacemakers. In Grace and Onn (1989), it was found that in vitro, rat midbrain dopamine neurons fired periodically with a frequency which could be increased or decreased by the application of depolarizing or hyperpolarizing currents without changing the pattern of spikes - see Figure 1 herein. The mechanisms of pacemaking activity of DA neurons in the rat ventral tegmental area (slice) was analyzed in Khaliq and Bean (2008).

The neurons with which we are mainly concerned are serotonergic neurons of the (dorsal) raphe nucleus, noradrenergic neurons of the locus coeruleus and midbrain dopaminergic neurons. The electrophysiological properties of these cells have been much investigated over the last several decades (Aghajanian and Vandermaelen, 1982; Aghajanian et al., 1983; Alreja and Aghajanian, 1991; Williams et al., 1984; Grace and Bunney, 1983, 1984; Harris et al., 1989). The roles of serotonergic and noradrenergic neurons in stress related disorders such as MDD and PTSD by means of reciprocal interactions with, inter alia, the HPA axis (especially through the PVN), hippocampus and prefrontal cortex are well documented (Lopez et al., 1999; Carrasco and Van de Kar, 2003; Lanfumey et al., 2008; Mahar et al., 2014). For example, CRH neurons of the PVN project directly to SE neurons of the DRN (Lowry et al., 2000) and NA neurons of the LC (Valentino et al., 1993) and stress, via an upregulation of the cAMP pathway, leads to an increase in the excitability of LC neurons (Nestler et al., 1999). In addition many of these brainstem neurons are endowed with glucocorticoid receptors which are activated by high levels of corticosterone (cortisol) (Joëls et al., 2007).

In most common experimentally employed animals except cat, the locus coeruleus is almost completely homogeneous, consisting of noradrenergic neurons which in rat number about 1500 (Swanson, 1976; Berridge and Waterhouse, 2003). The number of neurons in the rat DRN is between about 12000 and 15000 (Jacobs and Azmitia, 1992; Vertes and Crane, 1997) of which up to 50% are principal serotonergic cells (Vasudeva and Waterhouse, 2014) but there are also present  1000 dopaminergic cells (Lowry et al., 2008) and GABAergic cells, whose density varies throughout the divisions of the nucleus, as well as several other types of neuron. In rat, the number of dopamine neurons (bilateral counts) in the midbrain groups is approximately 71,000, with 40,000 in the VTA and 25,000 in substantia nigra (German and Manaye, 1993; Nair-Roberts et al., 2008). In the latter two cell groups there are also about 21,000 and 41,000 GABA-ergic cells, respectively.

The SE neurons of the DRN and NA neurons of the LC, often exhibit a slow regular pattern of firing with frequencies of order 0.5 to 2 Hz in slice (Milnar et al., 2016) and sometimes higher in vivo. Midbrain dopaminergic neurons exhibit similar firing patterns in vitro. The origins of pacemaker firing differ amongst various neuronal types. Thus, brainstem dopaminergic neurons may fire regularly without excitatory synaptic input (Grace and Bunney 1983; Harris et al., 1989). Underlying the rhythmic activity are subthreshold oscillations which were demonstrated with a mathematical model to possibly reflect an interplay between an L-type calcium current and a calcium-activated potassium current (Amini et al, 1999) - but see Kuznetsova et al. (2010). The mechanisms of pacemaker firing in LC neurons are not fully understood, although there is the possibility that it is sustained by a TTX-insensitive persistent sodium current (De Carvalho et al., 2000; Alvarez et al., 2002). For serotonergic neurons of the DRN, there have been no reports of a persistent sodium current and L-type calcium currents are relatively small or absent (Penington et al., 1991) so the main candidate for depolarization underlying pacemaking is a combination of T-type calcium current and the classical fast TTX-sensitive sodium current which dominates the pre-spike interval (Tuckwell and Penington, 2014). In some cells the hyperpolarization-activated cation current may also play a role.

In Figure 1 are shown portions of spike trains in NA neurons of mouse and rat LC, SE neurons in rat DRN and CRN and DA neurons in rat midbrain. It is noteworthy that the frequency of firing of LC and DRN principal neurons depends on sleep stage. Thus for example, in rats, waking, slow-wave sleep and REM sleep are accompanied by LC firing rates of about 2.2 Hz, 0.7 Hz and 0.02 Hz, respectively (Foote et al., 1980; Aston-Jones and Bloom, 1981; Luppi et al., 2012).

Figure 1: Some representative spikes from rat and mouse raphe nuclei and LC neurons. A. Part of a train of spikes in a rat LC neuron in slice. Markers 20 mV and 200 ms. (Andrade and Aghajanian, 1984). B. Detail of the course of the average membrane potential in a mouse LC neuron during an interspike interval. (De Oliveira et al., 2010). C. Action potentials in infant (7 to 12 days) and adult (8 to 12 weeks) mice. (De Oliveira et al., 2011). D. Whole-cell current-clamp recording of spikes in a mouse (21 to 32 days) LC neuron. (Sanchez-Padilla et al., 2014). E. A few spikes from rat dorsal raphe nucleus (slice), (Vandermaelen and Aghajanian, 1983). F. Portion of a spike train from rat caudal raphe nucleus. (Bayliss et al., 1997). G. Train of spontaneous spikes at a mean frequency of 0.85 Hz for a rat LC neuron in vitro. (Jedema and Grace, 2004). H. In vitro spiking in a rat dopaminergic neuron with two levels of applied current. (Grace and Onn, 1989).

Figure 2 shows computed spikes for the detailed model of rat DRN SE neurons of Tuckwell and Penington (2014) with 4 different parameter sets. In this model the main variables are membrane potential and intracellular calcium ion concentration which satisfy ordinary differential equations. There are, however, 11 membrane currents which drive the model, which results in a system with 18 components and over 120 parameters. In order to study quantities like interspike interval distributions with various sources of random synaptic input, it is helpful to have a simpler system of differential equations which might yield insight into the properties of the complex model whose execution with random inputs over hundreds of trials would be overly time consuming. It is also pointed out that LC and DRN principal neurons are each responsive to activation of about 20 different receptor types, making computational tasks even more cumbersome with an 18-component model (Kubista and Boehm, 2006; Maejima et al., 2013). A simplified model of spike generation would be useful in modeling the dynamics of serotonin release and uptake (Flower and Wong-Lin, 2014) and in simplified models of brainstem neural networks (Joshi et al., 2017). This approach has been proven useful by several authors, some of whom reduced the number of component currents (for example Fitzhugh, 1961; Nagumo et al., 1962; Destexhe et al., 2001) whereas others have simplified the geometry of the dendritic tree (Rall, 1962 ; Walsh and Tuckwell, 1985; Bush and Sejnowski, 1993).

Figure 2: Examples of computed spikes in a model for serotonergic neurons of the rat dorsal raphe nucleus from the model of Tuckwell and Penington (2014). Illustrated are the prolonged afterhyperpolarizations after a spike. The subsequent climb to threshold is plateau-like, sometimes being almost horizontal. Membrane potentials in mV are plotted against time in ms.

2 A reduced physiological model with two component currents

With the multi-component neuronal model for DRN SE neurons (Tuckwell and Penington, 2014) there are 11 component currents and a variety of solution behaviors for different choices of parameters. It is often difficult to see which parameters are responsible for various spiking properties. Some parameter sets lead to satisfactory spike trains with properties similar to experiment, but sometimes spike durations are unacceptably long. Other sets give rise to spontaneous activity but with spikes which have uncharacteristic large notches on the repolarization phase. Doublets and triplets are also often observed, and although these are sometimes found in experimental spike trains, it is desirable to find solutions depicting the regular singlet spiking and relatively smooth voltage trajectories usually observed.

Here we describe a simplified two-component model whose solutions exhibit membrane potential trajectories resembling those for pacemaker spiking in many brainstem neurons, as depicted in Figure 1. The several currents which act to mainly depolarize the neurons are lumped together in a single current denoted by , and similarly the currents which act to repolarize the cells are lumped into . Additionally there is an applied current which may be required to maintain pacemaker activity and which may also contain excitatory and inhibitory synaptic input. and also subsume a leak current. The basic differential equation for the voltage at time is written


where is capacitance. The initial value of is usually taken to be the resting potential, . All potentials are in mV, all times are in ms, all conductances are in S and capacitances are in nF.

The form of the current is based on that of the usual fast sodium transient currrent in the 11-component model for DRN SE cells and in the original Hodgkin-Huxley (1952) model. The current is based on the form of the delayed rectifier potassium current in the 11-component model. The inclusion of just the currents , , and if necessary a small constant driving applied current, is found, with appropriate sets of parameters, to be sufficient to give spikes whose trajectories resemble many found in pacemaker spiking of brainstem neurons as exemplified by those illustrated in Figures 1 and 2. These solutions have properties contrasting with the usual ones for the Hodgkin-Huxley system - see for example Hodgkin and Huxley (1952), Tuckwell and Ditlevsen (2016) and Figure 10 herein.

The depolarizing current, is given by the classical form


with activation variable , inactivation variable , maximal conductance and reversal potential . For the steady state activation we put


If the corresponding time constant is assumed to be voltage-dependent, we write it as


The steady state inactivation is given by


and if its time constant is voltage dependent, we put


We take the hyperpolarizing current to be


where (the traditional symbol) is the activation variable, is a usually positive integer-valued index, is the maximal conductance and is the reversal potential. The steady state activation for is written


and the time constant is, if voltage dependent,


Activation and inactivation variables for are determined by the first order equations


The activation variable for satisfies an equation like (10).

In some instances, the time constants are chosen to be constant, independent of V, in which case they are denoted by , and for , and , respectively.

2.1 Examples of pacemaker-like firing

We illustrate the pacemaker-like solutions with two choices of parameter sets. Set 1 is based on the parameters used for fast sodium and delayed rectifier potassium in the 11-current component model (Tuckwell and Penington, 2014). Set 2 is adopted from the parameters in a model of a rat sympathetic neuron with five conductances (Belluzzi and Sacchi ,1991). For the second set we also use the resting potential and the cell capacitance from Kirby et al. (2003) for rat DRN SE cells. The two sets of parameters are summarized in Table 1.

Parameter Set 1 Set 2
-33.1 -36
8 7.2
-50.3 -53.2
6.5 6.5
-60 -67.8
0.2 0.1
1.0 2.0
0.04 0.08861
-15 -6.1
7.0 8.0
1 1
- 3.5
1 -
4 -
-20 -
7 -
2.00 1.5
0.5 0.5
45 45
-93 -93
Table 1: Two basic parameter sets

Figure 3: Top. Repetitive spiking in the two-current () model for the parameter set 1 at the approximate threshold for spiking. Bottom left. The currents and during spikes. Bottom right. The function defined in Equ. (12) without added current (thin blue curve), whereby spiking does not occur, and with sufficient depolarizing current of 0.0342 nA (red curve) to give rise to pacemaker activity.

Both of these parameter sets led to repetitive spiking with the addition of a small depolarizing current (see Table 2). Typical spike trains are shown in Figures 3 and 4, and Table 2 contains lists of some of the details of the spike and spike train properties. Spiking for the second parameter set has a lower threshold for (repetitive) spiking, a longer ISI at threshold, a longer spike duration and a larger spike amplitude. For both parameter sets, most of these spike properties are in the ranges observed for DRN SE neurons. Figures 3 and 4 show well defined spikes with abruptly falling repolarization phases to a pronounced level of hyperpolarization followed by a steady increase in depolarization until an apparent spike threshold is reached.

Figure 4: Top. Repetitive spiking in the two-component model for the parameter set 2 at the approximate threshold for spiking. Bottom left. The currents and during spikes. Bottom right. The function defined in Equ. (12) without added current (thin blue curve), whereby spiking does not occur, and with sufficient depolarizing current of 0.018 nA (red curve) to give rise to pacemaker activity

In the lower left-hand panels of Figures 3 and 4 are shown, on an expanded time scale, the excitatory current and the inhibitory current during spikes. In the lower right-hand panels are shown plots (thin blue curves) of the function defined as the sum of the steady state () values of the quantities of Equ. (2) and Equ. (7),


The behavior of this function near the resting potential has been found to provide a heuristic indicator for the occurrence of spiking (Tuckwell, 2013; Tuckwell and Penington, 2014). It can be seen in both Figures 3 and 4 that, for both parameter sets, the function is negative for in an interval of considerable size around the resting potential , which is indicated by the vertical at -60 mV in Figure 3 and at - 67.8 mV in Figure 4. This means that around the derivative of with respect to time tends to be negative so that spontaneous spiking is unlikely. The magnitude of the smallest depolarizing current required to enable spiking is approximately the amount which must be added to make positive at and around . The resulting curves, shown in red in Figures 3 and 4, are approximately tangential to the V-axis, being obtained with for set 1 and for set 2. Hence these values of are estimates of the threshold depolarizing current required for (pacemaker) spiking.

Many brainstem (particularly DRN SE neurons) often have characteristically long plateau-like phases in the latter part of the ISI and this is apparent in Figures 3 and 4 for spikes elicited near the threshold for spiking for both sets of parameters. The plateau for the second set is nearly three times as long as that for the first set. For both parameter sets, (not shown) at a particular value of the frequency jumps from zero to a positive value, being 3.0 Hz for set 1 and 1.1 Hz for set 2, so that this model with the chosen parameters would be classified as one with Hodgkin (1948) type 2 neuron properties.

Figure 5 shows, for both parameter sets, the computed spike trajectories and ISIs for levels of excitation not much above threshold, being from 1 to 1.05 times the threshold values of . In each case the spike trajectory displays a typical pronounced and prolonged post-spike afterhyperpolarization followed by a plateau-like phase before the next spike.

The spiking properties for the two-component model with the parameter sets in Table 1 have several features in common with the experimental spike trains of many brainstem neurons including DRN SE and LC NA cells. The frequency of action potentials is in good agreement with experimental values near threshold. However, as the level of excitation is increased to much greater values, the frequency becomes somewhat high (not shown) relative to the most commonly reported values for sustained firing in these brainstem neurons, although there are exceptions; for example, in rat LC neurons, Korf et al. (1974) found in unanesthetized in vivo preparations, frequencies up to 30 Hz, and Sugiyama et al. (2012) reported in vivo frequencies over 7 Hz. In midbrain serotonergic raphe neurons, Kocsis et al. (2006) obtained in vivo rates with mean 5.4 Hz. For experiments with depolarizing current injection, Li and Bayliss (1998) obtained initial firing rates of around 8 Hz for caudal raphe with an injected current of 60 pA; Li et al. (2001) reported firing of rat DRN SE cells at frequencies as high as 35 Hz with current injection, and Ohliger-Frerking et al. (2003) found rates as high as 8 Hz and 11 Hz with 100 pA current injection in lean and obese Zucker rats, respectively. See also subsection 4.1 on firing rates.

Figure 5: Spike trajectories and ISIs obtained in the reduced model for the parameter sets of Table 1. Shown are results for Set 1 (left part) and Set 2 (right part), when applied currents are increased by relatively small amounts of 1%. 2% and 5% above threshold (top, middle, bottom curves).
Property Set 1 Basic set 2
Spike threshold -0.0342 -0.018
ISI at threshold 331 ms 948 ms
Spike duration (-40 mV) 1.6 ms 2.9 ms
Max V +8 +19.4
Min V -90.0 -91.2
Table 2: Spike train properties

The simplified model is expected to be useful in predicting approximate responses to random synaptic inputs in contradistinction to a sustained depolarizing current as employed in some experiments, which can lead to very large firing rates. However, with other parameter sets the model may not exhibit such high frequencies as depolarization level increases.

As Figure 1 shows, not all regular spiking in brainstem neurons entails a long afterhyperpolarization and a long plateau. Figure 6 shows that spikes similar to those in Figures 1B, 1D and 1F are generated in the second model with increased depolarizing current.

Figure 6: Spiking in the reduced model for the parameter sets of Table 3 where in both cases , a value considerably above threshold and which results in short ISIs, characterized by relatively short plateaux and short-lasting AHPs. Top part, set 1; bottom part, set 2.

3 Effects of varying parameters

The simplified model has 19 parameters for set 1, which compares with over 120 for the original model. The dependence of spiking on changes in the 19 parameters was investigated by computing solutions with each parameter above and below its standard value as given in Table 1. The focus in this subsection is the effect of parameter changes on the ISI, athough 6 other properties were routinely computed on each run, being spike duration, upslope on leading edge of spike, duration of AHP, amplitude of AHP, duration of plateau phase and approximate voltage at the half-way point of the plateau. The MATLAB computer programs employed to determine these properties are available from the authors on request.

The % change in the ISI for increases and decreases in parameters by 0.1 % is shown for 19 parameters in Table 3. The entries are ranked according to the degree to which they change the ISI, with rank 1 for the parameter with the greatest effect. Examination of these results reveals that small changes in 6 parameters, ranks 1 to 6, can result in a very large change in the ISI, making this quantity infinite as spiking is abolished altogether. Five of these 6 parameters are involved in the activation and inactivation of and the activation of , the other one being the reversal potential for . The inhibitory conductance, is the 7th ranked parameter making changes of 47.0 % and -17.8 % for increases and decreases of 0.1 %, respectively. Small changes in the parameters, ranked 8 to 11, make moderate changes in the ISI of order 10 to 20 %, being the excitatory conductance, the depolarizing drive current and the excitatory reversal potential. Changing the remaining parameters by 0.1% has minor to insignificant effects on the ISI. The parameter whose such relatively small changes have the least effect on the ISI is the resting potential, .

Rank Parameter +0.1 % -0.1 %
1 -48
3 -40.2
6 -26.6
11 2.7 -2.4
5 -28.7
2 -44.4
19 0.0 0.0
4 -30.3
15 0.009 -0.009
16 -0.007 0.010
13 -0.03 0.03
18 -0.002 0.002
17 0.007 -0.007
14 -0.01 0.01
8 -12.2 21.1
7 47.0 -17.8
9 12.1 -8.6
12 0.14 -0.14
10 -6.3 8.0
Table 3: Percentage changes in ISI for changes in parameter

More details of the effects of changing the parameters are given in Table 4. Here the parameters are grouped to reflect their different strengths in changing the ISI. The results point to the extreme sensitivity of the ISI to small relative changes in the following 8 parameters.

  • , the maximal hyperpolarizing conductance. Increases of as little as 0.005% are sufficient to prevent spiking. Decreases of 0.05% lead to a 43% decrease in the ISI, or about a 50% increase in firing rate.

  • , the half-activation potential for the depolarizing current, . When this parameter is increased by only 0.05%, spiking is extinguished and when it is decreased by 0.05% the ISI is shorter by 37%.

  • , the slope factor for the activation function of . Decreasing this quantity by 0.05% was sufficient to extinguish spiking, whereas increasing it by the same amount led to a 28% decrease in the ISI.

  • , the slope factor for the activation of . Decreasing this parameter by 0.05% also extinguished spiking and an increase by 0.05% gave a 33% reduction in ISI, or a 50% increase in firing frequency.

  • , the maximal depolarizing conductance. A reduction of as little as 0.05% ln this quantity was sufficient to abolish spiking, whereas a 0.05% increase resulted in a 72% decrease in the ISI.

  • , the reversal potential for the hyperpolarizing current, A decrease of 0.05% of this parameter gave a 65% increase in ISI, equivalent to a substantial reduction in firing rate, whereas a 0.05% increase led to a 20% decrease in the ISI.

  • , the whole cell capacitance. Reducing this quantity by the relatively small amount of 0.05% resulted in about a 50% increase in firing rate whereas a 0.5% increase in led to a 50% decrease in firing rate.

  • , the added current. Increasing this by 0.05% (that is making it less negative and hence less depolarizing) made the ISI infinite, whereas decreasing by 0.05% gave a 27% reduction in the ISI.

% Changes
0.05 0.01 0.005 -0.005 -0.01 -0.05
25.7 10.4 -7.6 -13.6 -37.0
-16.9 -4.2 -2.3 2.4 5.1 41.5
-18.6 -4.8 -2.5 2.8 6.4 53.7
-28.6 9.1 -4.9 5.8 13.0
1.3 0.3 0.1 -0.1 -0.3 -1.2
18.1 7.8 -6.1 -11.1 -32.8
0.0 0.0 0.0 0.0 0.0 0.0
-19.8 -5.3 -2.8 3.1 6.5 65.1
% Changes 50 25 10 -10 -25 -50
4.7 2.5 1.0 -1.0 -2.6 -5.3
-4.2 -2.1 -0.8 0.7 1.8 1.8
% Changes 50 20 10 -10 -20 -50
-17.2 -9.3 -3.4 3.3 6.6 15.4
-1.3 -0.5 -0.2 0.3 0.5 1.3
1.6 1.0 0.6 -0.8 1.9 -7.8
-9.7 -3.1 -1.4 1.1 1.8 2.5
% Changes 10 5 1 -1 -5 -10
-72.4 -66.7 -46.1 8.9 21.1
-10.4 -17.8 -43.1
% Changes 0.5 0.1 0.05 -0.05 -0.1 -0.5
12.1 5.5 -4.6 -8.6 -27.9
% Changes 50 25 10 -10 -25 -50
66.2 33.3 13.4 -13.5 -33.7 -67.4
% Changes 1 0.5 0.1 -0.1 -0.5 -1
-33.1 -22.3 -6.3 8.0 111.8
Table 4: % changes in ISI for various % changes in set 1 parameters

Figure 7: Voltage versus time plots for various values of the half-activation potential of the excitatory current. A. The standard set of parameters including the value of at -33.1 mV. Pacemaker activity is sustained with an ISI of 331 ms. B. The value of is slightly above the critical value for spiking at 0.02480 % above the standard value. The ISI is about 13000 ms. C. With the value of slightly higher at 0.02481 % above the standard value, there is no spiking and approaches a limiting value of about -53 mV after about 200 ms. D. An even larger value of , 2.5% above -33.1 mV, results in an asymptotically smaller departure from rest to about -56 mV which is achieved in about 20 ms.

Figure 8: Plots of versus for the corresponding plots A-D of versus in Figure 7. Note the similarity of the trajectories in A and B, despite enormously different ISIs. In C and D, corresponding to only subthreshold responses, trajectories are shown approaching asymptotic values at .

Changes in the remaining 11 parameters have moderate to minor effects on the ISI - see also Table 3.

3.1 More detailed analysis for and

In this subsection the transition from non-spiking to spiking for two of the parameters is examined more closely near the critical values.

3.1.1 Half-activation potential for

The parameter ranked 1 in Table 3 is which is the half-activation potential for the activation variable of the excitatory current. Its standard value was -33.1 mV which gave rise to pacemaker-type spiking with an ISI of 331 ms. Since only a small increase in this parameter abolished spiking, it is of interest to study the behavior of solutions near a critical value above -33.1 at which no spiking occurs. To this end, solutions were computed with set 1 parameters with varying values of . Results are shown in Figures 7 and 8. In Figure 7 are shown computed spike trains, i.e., versus , whereas in Figure 8 are shown plots of versus for the indicated time periods.

In panel A, the value of is the standard value (-33.1) and the spike train is as described above with a constant ISI (Figure 7A). The trajectories of versus for each spike are practically identical, motion being clockwise on these orbits (Figure 8A).

In the panels B, the value of is just 0.024800 % above the standard value and this is very close to the critical value (call it ) beyond which spiking does not occur because it is too far from the resting potential. As seen in Figure 7B, the first spike occurs at 12980 ms and the second at 26010 ms, but somewhat remarkably, the orbits (for the two spikes)as seen in Figure 8B are visually indistinguishable from those in Figure 8A where the ISI is only about 1/40 as long at 331 ms.

The versus curve for a value of which is slightly larger at 0.024810 % above the standard value, and presumably greater than the critical value , is shown in Figure 7C. Here there are no spikes (up to 30000 ms) and the values of have limiting values of about (-53.0038, 0.004368) as seen in Figure 8C. Again, somewhat remarkably, the trajectory from rest to lies almost exactly on the spiking orbit as it approaches the limiting value.

With the considerably higher value of , at 2.5 % above the standard value, there are, as expected, no spikes as seen in Figure 7D. The equilibrium point is closer to rest at =(-56.4600, 0.002670). The trajectory in space still adheres to the early part of the spike orbit.

3.1.2 Magnitude of the added current: comparison with usual Hodgkin-Huxley parameters

The parameter ranked 9 in Table 3 is which is the magnitude of the added depolarizing current which may be required to induce pacemaker-type spiking. Its set 1 value was which was near the threshold value for spiking when all other parameters had their set 1 values. The solutions are not as sensitive to changes in this parameter as to those in . Figure 9 shows the plots of versus for four values of .

In the left panel orbits are shown for the the standard value -0.0342 of , which gives an ISI of 331 (red curves), and for a value 0.44 % above (more positive and hence more hyperpolarizing) the standard value (green curves) which yields an extremely small spike rate with the first spike not occurring until 13,700 ms. The orbits for these two values of are almost coincident, but one can on close inspection distinguish parts of the red and green curves. The critical value of , called , for spiking is apparently not far above the second value.

In the right panel of Figure 9 are shown orbits, effectively from to , for two values of above , at which there is no evidence of spiking. The red curve is for a value of just 0.5 % above the standard value. The asymptotic value of is with and . Again the curve closely adheres to the early part of the spike trajectory. This is also the case for the much higher value of at 20 % above the standard value. The curve for this case is blue and terminates a shorter distance from the rest point, with and .

Usual HH parameters The manner in which the periodic orbit emerges as passes through a critical value in the above calculations is seemingly quite different from that in which periodic spiking arises in the HH model with the usual parameters as given, for example, in Tuckwell and Ditlevsen (2016). The latter case is illustrated in Figure 10 where there are shown in A, B, C, respectively the responses to input currents well below threshold for periodic spiking (), not far below () and above threshold (). Here the subthreshold responses are oscillatory and decay is into the rest point. This is in contradistinction to the behavior of solutions with the parameters of the present model - see Figure 9. The nature of the bifurcation to periodic solutions that occurs at the critical values for the original Hodgkin-Huxley parameters was discussed in Tuckwell et al. (2009). Further analysis is required to see if the bifurcation to periodic solutions in the present model has a different mathematical structure from that in the original model.

Figure 9: Plots of versus as the depolarizing current is changed above and below the critical value for spiking. In the left panel are shown spiking orbits for the standard value of =-0.0342 (red curves) which gives an ISI of 331 ms, and a value which is 0.44 % above the standard value (green curve) which results in an ISI of about 14,000 ms. The phase diagrams for these widely disparate values of the ISI are almost identical. In the right hand panel are orbits for values of which are too large (hyperpolarizing) to induce spiking. The orbits again asymptote to subthreshold values at . These subthreshold orbits coincide with the early parts of the spiking orbits shown in the left panel.

Figure 10: Solutions versus (left column) and phase diagrams versus (right column) for the Hodgkin-Huxley equations with the original parameters. In A, the input current is well below the spiking threshold, in B it is not far below and in C, it is above the threshold for repetitive spiking. The subthreshold behaviors are quite different from those illustrated in Figure 9, right hand diagram.

3.2 Changing plateau level and the minimum of the AHP

Examination of the membrane potential trajectories during pacemaker-type firing often shows a relatively long plateau between successive spikes during which the membrane potential increases slowly, sometimes very slowly before climbing past threshold to give rise to the subsequent spike. Such plateaux are exemplified by the recordings shown in Figure 11, where the dopamine neuron spikes are from Grace and Bunney (1995), those for noradrenaline from Andrade and Aghajanian (1984), and those for serotonin from Vandermaelen and Aghajanian (1983). In the case of the in vitro recording shown for the DA cell, the plateau phase is not so flat, but an example with a flat plateau phase is shown in Figure 1H (from Grace and Onn, 1989) for a cell with an applied hyperpolarizing current of 0.07 nA. The resting potential for the DA cell shown in Figure 11 (middle dotted line) is -50 mV, for rat NA cells the average resting potential is -58.2 mV based on 11 sets of experimental results, and the average resting potential for SE neurons of the DRN is about -60 mV.

Figure 11: Pacemaker-type spiking in a DA neuron (top) from Grace and Bunney (1995) (also in Grace and Onn, 1989), from an NA neuron (middle) from Andrade and Aghajanian (1984), and an SE neuron (bottom) from Vandermaelen and Aghajanian, (1983). In each trajectory there is evidence of an extended plateau phase during which the membrane potential undergoes relatively minor changes before ascending at the end of the ISI to give rise to a subsequent spike.

The DA and NA recordings are in vitro recordings from rat substantia nigra and locus coeruleus, respectively. It is likely that regular pacemaker-type spiking in these two types of cells usually only occurs in vitro. Grace and Onn (1989) stated that DA cells recorded in vitro fired exclusively in a highly regular pacemaker-like pattern. This is supported by results reported by Grenhoff et al. (1988) for DA neurons and by Svensson et al. (1989) on NA neurons. In these two sets of in vivo experiments, application of the excitatory amino acid antagonist kynurenate (or kynurenic acid) changed the irregular (bursting) nature of spiking of DA and NA neurons to regular pacemaker type activity.

Since in vitro recordings usually have diminished and possibly no synaptic input, the findings in these latter two references offer a plausible explanation why regular pacemaking often occurs in vitro but rarely in vivo (see for example Sugiyama et al., 2011). Findings related to that of Grenhoff et al. (1988) were that the GABA agonist baclofen made firing in DA neurons slower and more regular (Engberg et al.,1993), and that similar effects resulted from the application of NMDA receptor blockers (Overton and Clark, 1992; Cherugi et al., 1993). In another more recent study, Blythe et al. (2009) found that in DA neurons of the substantia nigra, bursting could be induced by the application of glutamate at somatic or dendritic sites. In a related interesting investigation in DA neurons, Putzier et al. (2009) found that the voltage characteristics of L-type channels rather than selectivity were important factors in pacemaking and bursting activity.

For SE neurons, regular pacemaking has been found in both in vivo and in vitro experiments, suggesting less fine tuning by excitatory synaptic inputs with a more robust pacemaker mechanism. Levine and Jacobs (1992) applied kynurenic acid to cat DRN serotonergic neurons and found that their spontaneous firing rates were only reduced by about 4 %.

It can be seen that the voltage levels of the plateaux between the NA and SE spikes in Figure 11 are often not much above resting potential. However for the model with set 1 parameters, as shown in Figures 3 and 5, the plateau level is about 6.8 mV above resting value. It turns out that the plateau level and minimum can be easily adjusted with a change of parameters which are voltages, being in particular . These were all changed together by an amount . With mV applied to set 1 parameters no spikes emerged.

With all set 1 parameters unchanged with the exception that the reversal potential for the inhibitory current was at a more depolarized level of mV, spiking was very rapid with an ISI of only 67.8 ms. If the maximal inhibitory conductance was then multiplied by 1.1427, the ISI was back to 338.9 ms, close to the value for set 1 parameters. Calling this set 3, with altered and altered , the plateau level was 7 mV above resting level. When the set was uniformly decreased by 5mV, (set 3’), the ISI was virtually unchanged, but the plateau level was reduced to only 2mV above resting potential. The spike trains for sets 3 and 3’ are shown in Figure 12, with red and blue curves, respectively. The minimum value of during the AHP for set 3 is -85.7 mV and -90.7 for set 3’. Some of the remaining spike train properties for set 3 parameters (with corresponding properties for set 3’ in brackets) are: spike duration 1.6 ms (1.27), upslope of spike 100 V/s (111), duration of AHP 38.8 ms (64), AHP amplitude 25.7 mV (30.7) and plateau duration 270 ms (247). Several of these properties are quite different for sets 3 and 3’ and it is remarkable that the ISIs are almost identical. However, it has been shown that the plateau level can be easily adjusted to various values above resting level.

Figure 12: Adjustment of plateau level and voltage minimum during AHP in pacemaker-type spiking in the 2-component model with parameters given in text (red curve) as well as with changes of mV in the 6 key voltages (blue curve).

4 Firing frequency and bursting

4.1 Firing rates

Examination of the experimental literature on firing rates and their responses to applied currents for the three main types of neuron, DA, NA, and SE, considered above, reveals a wide range of often disparate results. The origins of such variability lie in several factors including animal types and ages and experimental conditions and techniques. Generally, in pacemaking mode, firing rates of locus coeruleus (NA) neurons (Tuckwell, 2017, contains a summary) are less than 5 Hz with a few in vivo experiments reporting higher or much higher rates. For serotonergic neurons of the DRN (summarized in Tuckwell, 2013), firing rates are usually less than a few Hz, rarely around 4-5 Hz and with applied current up to 10 Hz (Ohliger-Frerking et al., 2003). Rates as high as 20 Hz in cat DRN were reported under glutamate application by Levine and Jacobs (1992). Bursting is also sometimes reported in these cells (Hájós et al., 1996; Hájós et al., 2007). As remarked above, DA cells in vivo generally exhibit burst firing and in vitro pacemaking was found to have an average value of 5.4 Hz in one experiment (Grenhoff et al., 1988). Firing rates for set 1 parameters in the simplified model were computed for both more hyperpolarizing applied currents (Figure 13A) and more depolarizing currents (Figure 13B) relative to the standard value of 34.2 pA. Frequencies ranged from 0 to 3 Hz in Figure 9A and from 3 to about 8 Hz in Figure 9B, the rather high frequencies in the latter resulting from large applied depolarizing currents. The frequencies for set 2 parameters were considerably less (see Figures 4 and 5).

Figure 13: A. Frequency of firing for set 1 parameters with values of which are more hyperpolarizing than the standard value of nA. B. As in part A but that applied currents are more depolarizing than nA.

4.2 Bursting

The results on the mechanisms of bursting in DA and LC neurons mentioned above (Grenhoff et al.,1988; Svensson et al., 1989; Overton and Clark, 1992; Cherugi et al., 1993; Engberg et al.,1993; Blythe et al., 2009) indicate that in several instances regular pacemaking in these cells may be replaced by bursting if a (synaptic) excitatory drive is activated or an inhibitory drive inactivated. This was explored in the simplified model by briefly perturbing the constant excitatory conditions that led to regular pacemaking. The perturbations mimic the effects of either increased excitatory input or decreased inhibitory input. The results are shown in Figure 14. In all cases the perturbations lead to rapid spiking which mimics bursting over the duration of the perturbation with a subsequent return to pacemaker spiking when the perturbation is terminated. In Figures 14A and 14B, the perturbations consist of increases in the maximal excitatory conductance for 200 ms (14A) or decreased inhibitory maximal conductance (14B). In Figures 14C and 14D, the depolarizing current parameter is increased for either 200 ms (14C) or 100 ms (14D). In all four cases the evidence is convincing that such perturbations may lead to a disturbance of pacemaker firing. Thus, a steady but intermittent flow of excitatory synaptic input or a reduction of inhibitory synaptic input could lead to non-pacemaker firing and occasional or frequent bursting.

Figure 14: Induction of bursting-type spiking by means of brief changes in excitatory or inhibitory parameters. Parameters are as in set 1 except with variations in a single parameter as described. A. Increasing the maximal excitatory conductance from the standard value of 2 to 2.5 for gives rise to a burst of three spikes, followed by a return to pacemaker activity. B. Reducing the maximal inhibitory conductance from the standard value of 0.5 to 0.3 for the same time period as in A produces a burst of 4 spikes with a subsequent return to pacemaker-type spiking. C. A 10% increase in the depolarizing current for the same time-period as in A and B gives rise to rapid spiking followed by the resumption of regular pacemaking activity. D. As in C but now a 50% increase in depolarizing current applied for the shorter time period gives rise to more rapid bursting-type spiking followed regular pacemaking activity.

5 Discussion

Principal brainstem neurons, particularly serotonergic cells of the dorsal raphe nucleus and noradrenergic cells of the locus coeruleus are of great importance in the functioning of many neuronal populations throughout cortical and subcortical structures. Of note is the modulatory role the firing of neurons in these brainstem nuclei have on neurons of the prefrontal cortex, including the orbitofrontal cortex, and hippocampus. These latter structures have been strongly implicated in various pathologies, including depression, and OCD. Lanfumey et al. (2008) contains a comprehensive summary of many of the biological processes which are influenced by serotonin including those originating from serotonergic neurons of the DRN. Modeling networks involving both serotonergic and noradrenergic afferents requires plausible models for the spiking activity of the principal SE and NA cells. Whereas detailed models of such activity are now available, their application to many thousands of cells has the disadvantage of leading to very large computation time and large memory requirements, so that the simplified models described in the present article may provide useful approximating components for such complex computing tasks. Accurate models of brainstem neurons (see below) involve many component currents and very large numbers of parameters, several of which are possibly uncertain. Hence realistic simplified mathematical models of brainstem neurons, beyond that provided by extremely simplified models such as the leaky integrate and fire (or Lapicque) model (Tuckwell, 1988), are useful in order to investigate approximately the responses of these cells to their complex array of synaptic and other input and to construct and analyze complex networks involving these cells and those in other centers such as hippocampus, frontal cortex and hypothalamus.

The simplified model we have used is a modified Hodgkin-Huxley model with currents relabeled generically. With certain parameter sets such a model predicts spike trajectories which are similar to those of brainstem neurons in pacemaker mode. The sensitivity of solutions to changes in parameters was investigated in detail and very small changes in some parameters were found to produce dramatic changes in the spiking behavior. Details of such features as plateaux and AHPs could be adjusted to match experimental results with judicious choices of voltage parameters. A preliminary investigation was able to convert regular pacemaking to bursting behavior in accordance with experimental findings on the application of kynurenic acid to DA and NA neurons. We end with a brief discussion of previous computational modeling of DA, NA and SE brainstem neurons

5.1 LC NA neurons

Thus far there have been several mathematical models of locus coeruleus neurons per se, which include a few ionic channels (Putnam et al., 2014; Contreras et al., 2015) or many ionic channels including the usual sodium and potassium, high and low threshold calcium currents, transient potassium , persistent sodium, leak and hyperpolarization activated cation current (De Carvalho et al., 2000; Alvarez et al. 2002; Carter et al., 2012). Noteworthy is the omission of in the model of Alvarez et al. (2002) and its inclusion in De Carvalho et al. (2000) and Carter et al. (2012). Also, a persistent sodium current is included in De Carvalho et al. (2000) and Alvarez et al. (2002) but not in Carter et al. (2012).

Despite such uncertainties in the mechanisms involved in pacemaker activity in LC neurons, some of these works have included synaptic input and gap-junction inputs from neighboring LC neurons. The pioneering article of De Carvalho et al. (2000) addressed the mechanisms of morphine addiction and included several biochemical reactions involving cAMP, -opioid receptors, morphine, G-protein, AC, CREB and Fos. Tuckwell (2017) contains a summary of previous LC modeling as well as a review of LC neuron anatomy and physiology. Brown et al. (2004) did employ Rose-Hindmarsh model neurons to study a network of LC neurons but there have not appeared any plausible simplified models of these cells per se. Thus the two-component model considered in this article provides a good starting point for investigating, for example, the effects of synaptic input on LC firing which will be performed in future articles.

5.2 DRN SE neurons

For serotonergic neurons of the dorsal raphe there has been only one detailed model as described in the introduction (Tuckwell and Penington, 2014). Some authors have addressed quantitatively serotonin release and included the effects of antidepressants but without an explicit model for SE cell spiking (Geldof et al., 2008). Wong-Lin et al. (2012) use a quadratic integrate and fire model for spiking DRN SE neurons in a network of such cells along with inhibitory neurons. The model has a reset mechanism after spikes which is artificial. In a related work, Cano-Colino et al. (2013, 2014) have modeled the influence of serotonin on networks of excitatory and inhibitory cells in spatial working memory. More recently, in a similar vein, Maia and Cano-Colino (2015) have made an interesting study of serotonergic modulation of the strength of attractors in orbitofrontal cortex and related this to the occurrence of OCD.

5.3 Midbrain DA neurons

Midbrain dopaminergic neuons have been implicated, inter alia, in reward pathways and in several pathological conditions including Parkinson’s disease (Kalia and Lang, 2015), psychiatric conditions such as schizophrenia (Weinstein et al., 2017) and addiction to drugs and natural rewards (Nestler, 2004). There have been several comprehensive computational modeling studies of these neurons, many of which were carried out by Canavier and coworkers since 1999. These include the model of Amini et al. (1999) with four calcium currents, two potassium currents, a hyperpolarization activated current and various pump currents which was used to explore the nature and role of a slow oscillatory potential which was observed when spiking was blocked by tetrodotoxin. Compartmental models were considered subsequently by Komendantov et al. (2004) and Kuznetsova et al. (2010), with soma and two dendritic compartments. In the second of these works it was found that spontaneous firing rates were determined mainly by an L-type calcium current and the A-type potassium current. In a more recent experimental and computational study using a compartmental model, Tucker et al. (2012) found that pacemaker frequency was dependent on the density of Na channels in the soma and their spatial distribution over the soma-dendritic regions.

6 Acknowledgements

This research was partially supported by a Mathematical Biosciences Institute postdoctoral fellowship to YZ from the National Science Foundation under Agreement No. 0931642.

7 References

Aghajanian, G.K., Vandermaelen, C.P., 1982. Intracellular recordings from serotonergic dorsal raphe neurons: pacemaker potentials and the effect of LSD. Brain Res. 238, 463-469.

Aghajanian, G.K., Vandermaelen, C.P., Andrade, R., 1983. Intracellular studies on the role of calcium in regulating the activity and reactivity of locus coeruleus neurons in vivo. Brain Res. 273, 237-243.

Alreja, M., Aghajanian, G.K., 1991. Pacemaker activity of locus coeruleus neurons: whole-cell recordings in brain slices show dependence on cAMP and protein kinase A. Brain Res. 556, 339-343.

Alvarez, V.A., Chow, C.C., Van Bockstaele, E.J., Williams, J.T., 2002. Frequency-dependent synchrony in locus ceruleus: role of electrotonic coupling. PNAS 99, 4032-4036.

Amini, B., Clark, J.W. Jr, Canavier, C.C., 1999. Calcium dynamics underlying pacemaker-like and burst firing oscillations in midbrain dopaminergic neurons: A computational study. J. Neurophysiol. 82, 2249-2261.

Andrade. R., Aghajanian, G.K., 1984. Locus coeruleus activity in vitro: intrinsic regulation by a calcium-dependent potassium conductance but not -adrenoceptors J. Neurosci 4, .161-170.

Aston-Jones, G., Bloom, F.E., 1981. Activity of norepinephrine-containing locus coeruleus neurons in behaving rats anticipates fluctuations in the sleep-waking cycle. J. Neurosci. 1, 876-886.

Bayliss, D.A., Li, Y.-W., Talley, E.M., 1997a. Effects of serotonin on caudal raphe neurons: activation of an inwardly rectifying potassium conductance. J. Neurophysiol. 77, 1349–1361.

Belluzzi, O., Sacchi, O.,1991. A five-conductance model of the action potential in the rat sympathetic neurone. Prog. Biophys. Molec. Biol. 55: 1-30

Berridge, C.W., Waterhouse, B.D., 2003. The locus coeruleus–noradrenergic system: modulation of behavioral state and state-dependent cognitive processes. Brain Res. Rev. 42, 33-84.

Blythe, S.N., Wokosin, D., Atherton, J.F., Bevan, M.D., 2009. Cellular mechanisms underlying burst firing in substantia nigra dopamine neurons. J. Neurosci. 29, 15531-15541.

Brown, E., Moehlis, J., Holmes, P. et al., 2004. The influence of spike rate and stimulus duration on noradrenergic neurons. J. Comp. Neurosci. 17, 13-29.

Bush, P.C., Sejnowski, T.J., 1993. Reduced compartmental models of neocortical pyramidal cells. Journal of Neuroscience Methods 46,159-166.

Cano-Colino, M., Almeida, R., Compte, A., 2013. Serotonergic modulation of spatial working memory: predictions from a computational network model. Frontiers in Integrative Neuroscience 7, 71.

Cano-Colino, M., Almeida, R., Gomez-Cabrero, D. et al., 2014. Serotonin regulates performance nonmonotonically in a spatial working memory network. Cerebral Cortex 24, 2449-2463.

Carrasco, G.A., Van de Kar, L.D., 2003. Neuroendocrine pharmacology of stress. Eur. J. Pharm. 463, 235-272.

Carter, M.E., Brill, J., Bonnavion, P. et al., 2012. Mechanism for hypocretin-mediated sleep-to-wake transitions. PNAS 109, E2635–E2644.

Cherugi, K., Charléty, P.J., Akaoka, H. et al., 1993. Tonic activation of NMDA receptors causes spontaneous burst discharge of rat midbrain dopamine neurons in vivo. Eur. J. Neurosci. 5,137-144.

Contreras, S., Quintero, M., Putnam, R., Santin, J., Hartzler, L. and Cordovez, J., 2015. A Computational Model of Temperature-Dependent Intracellular pH Regulation. The FASEB Journal, 29(1 Supplement), pp.860-10

De Carvalho, L.A.V., De Azevedo, L.O., 2000. A model for the cellular mechanisms of morphine tolerance and dependence. Math. Comp. Mod. 32, 933-953.

De Oliveira, R.B., Howlett, M.C.H., Gravina, F.S. et al., 2010. Pacemaker currents in mouse locus coeruleus neurons. Neurosci 170, 166-177.

De Oliveira, R.B., Gravina, F.S., Lim, R. et al., 2011. Developmental changes in pacemaker currents in mouse locus coeruleus neurons. Brain Res. 1425, 27-36.

Engberg, G., Kling-Peterson, T., Nissbrandt, H., 1993. GABA-receptor activation alters the firing pattern of dopamine neurons in the rat substantia nigra. Synapse 15, 229-238.

Fitzhugh, R.,1961. Impulses and physiological states in theoretical models of nerve membrane. Biophys. J. 1, 445-466.

Flower, G., Wong-Lin, K., 2014. Reduced computational models of serotonin synthesis, release and reuptake. IEEE Trans. Biomed. Eng. 61, 1054-1061.

Foote, S.L., Aston-Jones, G., Bloom, F.E., 1980. Impulse activity of locus coeruleus neurons in awake rats and monkeys is a function of sensory stimulation and arousal. Proc. Nati. Acad. Sci. USA 77, 3033-3037.

Geldof, M., Freijer, J.I., Peletier, L.A. et al., 2008. Mechanistic model for the acute effect of fluvoxamine on 5-HT and 5-HIAA concentrations in rat frontal cortex. Eur. J. Pharm. Sci. 33, 217-219.

German, D.C. and Manaye, K.F., 1993. Midbrain dopaminergic neurons (nuclei A8, A9, and A10): Three-dimensional reconstruction in the rat. J. Comp. Neurol. 331, 297-309.

Grace, A. A., Bunney, B.S., 1983. Intracellular and extracellular electrophysiology of nigral dopaminergic neurons. I. Identification and characterization. Neurosci. 10: 301-315.

Grace, A., Bunney, B., 1984. The control of firing pattern in nigral dopamine neurons: burst firing. J. Neurosci. 4, 2877-2890.

Grace, A.A., Onn, S-P., 1989. Morphology and electrophysiological properties of immunocytochemically identified rat dopamine neurons recorded in vitro. J. Neurosci. 9, 3463-3481.

Grace, AA., Bunney, B.S., 1995. Electrophysiological Properties of Midbrain Dopamine Neurons. Neuropsychopharmacology - 4th Generation of Progress. Eds. Bloom, F.E., Kupfer, D.J. Raven Press, New York.

Grenhoff, J., Tung, C-S., Svensson, T.H., 1988. The excitatory amino acid antagonist kynurenate induces pacemaker-like firing of dopamine neurons in rat ventral tegmental area in vivo. Acta Physiol. Scand. 134, 567-568.

Hájós, M., Sharp, T., Newberry, N.R., 1996. Intracellular recordings from burst-firing presumed serotonergic neurones the rat dorsal raphe nucleus in vivo. Brain Res. 737, 308-312.

Hájós, M., Allers, K.A., Jennings, K. et al., 2007. Neurochemical identification of stereotypic burst-firing neurons in the rat dorsal raphe nucleus using juxtacellular labelling methods. Eur. J. Neurosci. 25, 119-126.

Harris, N.C., Webb, C., Greenfield, S.A., 1989. A possible pacemaker mechanism in pars compacta neurons of the guinea-pig substantia nigra revealed by various ion channel blocking agents. Neurosci. 31, 355-362.

Hodgkin, A.L., 1948. The local changes associated with repetitive action in a non-medullated axon. J. Physiol. 107, 165-181.

Hodgkin, A.L., Huxley, A.F., 1952. A quantitative description of membrane current and its application to conduction and excitation in nerve. J. Physiol. 117, 500-544.

Jedema, H.P., Grace, A.A., 2004. Corticotropin-releasing hormone directly activates noradrenergic neurons of the locus ceruleus recorded in vitro. J. Neurosci. 24, 9703-9713.

Ishimatsu, M., Williams, J.T., 1996. Synchronous activity in locus coeruleus results from dendritic interactions in pericoerulear regions. J. Neurosci. 16, 5196-5204.

Jacobs, B.L., Azmitia, E.C., 1992. Structure and function of the brain serotonin system. Physiol. Rev. 72, 165-229.

Joëls, M, Karst, H., Krugers, H.J., Lucassen, P.J., 2007. Chronic stress: Implications for neuronal morphology, function and neurogenesis. Front. Neuroendocrin. 28, 72-96.

Joshi, A., Youssofzadeh, V., Vemana, V., McGinnity, T.M., Prasad, G., Wong-Lin, K., 2017. An integrated modelling framework for neural circuits with multiple neuromodulators. Journal of The Royal Society Interface 14, 20160902.

Kalia, L.V., Lang, A.E., 2015. Parkinson’s disease. Lancet 386, 896-912.

Khaliq, Z.M., Bean, B.P., 2008. Dynamic, nonlinear feedback regulation of slow pacemaking by A-type potassium current in ventral tegmental area neurons. J. Neurosci. 28, 10905-10917.

Kirby, L.G., Pernar, L., Valentino, R.J. et al., 2003. Distinguishing characteristics of serotonin and nonserotonin- containing cells in the dorsal raphe nucleus: electrophysiological and immunohistochemical studies. Neurosci. 116, 669-683.

Kocsis B, Varga V, Dahan L, Sik A (2006) Serotonergic neuron diversity: Identification of raphe neurons with discharges time-locked to the hippocampal theta rhythm. PNAS 103: 1059-1064.

Komendantov, A.O., Komendantova, O.G., Johnson, S.W., Canavier, C.C., 2004. A modeling study suggests complementary roles for GABA and NMDA receptors and the SK channel in regulating the firing pattern in midbrain dopamine neurons. J. Neurophysiol. 91, 346-357.

Korf, J., Bunney, B.S., Aghajanian, G.K., 1974. Noradrenergic neurons: morphine inhibition of spontaneous activity. Eur.J. Pharmacol. 25, 165-169.

Kubista, H., Boehm, S., 2006. Molecular mechanisms underlying the modulation of exocytotic noradrenaline release via presynaptic receptors. Pharmacol. Therap. 112, 213-242.

Kuznetsova, A.Y., Huertas, M.A., Kuznetsov, A.S., Paladini, C.A., Canavier, C.C., 2010. Regulation of firing frequency in a computational model of a midbrain dopaminergic neuron. Journal of computational neuroscience 28, 389-403.

Lanfumey ,L., Mongeau, R., Cohen-Salmon, C., Hamon, M., 2008. Corticosteroid-serotonin interactions in the neurobiological mechanisms of stress-related disorders. Neurosci. Biobehav. Rev. 32, 1174-1184.

Levine, E.S., Jacobs, B.L., 1992. Neurochemical afferents controlling the activity of serotonergic neurons in the dorsal raphe nucleus: microiontophoretic studies in the awake cat. J. Neurosci. 12, 4037-4044.

Li, Y-Q., Li, H., Kaneko, T., Mizuno, N., 2001. Morphological features and electrophysiological properties of serotonergic and non-serotonergic projection neurons in the dorsal raphe nucleus An intracellular recording and labeling study in rat brain slices. Brain Res. 900, 110-118.

Li, Y-W., Bayliss, D.A., 1998. Electrophysiological properties, synaptic transmission and neuromodulation in serotonergic caudal raphe neurons. Clin. Exp. Pharm. Physiol. 25, 468-473.

Lopez. J.F., Akil, H., Watson, S.J., 1999. Neural circuits mediating stress. Biol. Psychiatry 46:1461-1471.

Lowry, C.A., Evans, A.K., Gasser, P.J. et al., 2008. Topographic organization and chemoarchitecture of the dorsal raphe nucleus and the median raphe nucleus. In: Serotonin and sleep: molecular, functional and clinical aspects, p 25-68, Monti, J.M. et al., Eds. Basel: Birkhauser Verlag AG.

Lowry, C.A., Rodda, J.E., Lightman, S.L., Ingram, C.D., 2000. Corticotropin- releasing factor increases in vitro firing rates of serotonergic neurons in the rat dorsal raphe nucleus: evidence for activation of a topographically organized mesolimbocortical serotonergic system. J. Neurosci. 20, 7728-7736.

Luppi. P-H., Clement, O., Sapin, E. et al., 2012. Brainstem mechanisms of paradoxical (REM) sleep generation. Eur. J. Physiol. 463:43-52.

Maejima, T., Masseck, O.A., Mark, M.D., Herlitze, S., 2013. Modulation of firing and synaptic transmission of serotonergic neurons by intrinsic G protein-coupled receptors and ion channels. Frontiers in Integrative Neuroscience 7, 40.

Mahar, I., Bambico, F.R., Mechawar, N., Nobrega J.N., 2014. Stress, serotonin, and hippocampal neurogenesis in relation to depression and antidepressant effects. Neurosci. Biobehav. Rev. 38, 173-192.

Maia, T.V., Cano-Colino, M., 2015. The role of serotonin in orbitofrontal function and obsessive-compulsive disorder. Clinical Psychological Science 3, 460-482, and Supplemental-data.

Milnar, B., Montalbano, A., Piszczek et al., 2016. Firing properties of genetically identified dorsal raphe serotonergic neurons in brain slices. Front. Cell. Neurosci. 10, 195.

Nagumo, J. S., Arimoto, S., Yoshizawa, S., 1962. An active pulse transmission line simulating nerve axon. Proc. I.R.E. 50, 2061-2070.

Nair-Roberts, R.G., Chatelain-Badie, S.D., Benson, E., White-Cooper, H., Bolam, J.P., Ungless, M.A., 2008. Stereological estimates of dopaminergic, GABAergic and glutamatergic neurons in the ventral tegmental area, substantia nigra and retrorubral field in the rat. Neurosci. 152, , 1024-1031.

Nestler, E.J., 2004. Molecular mechanisms of drug addiction. Neuropharmacol. 47, 24-32.

Nestler, E. J., Alreja, M., Aghajanian, G.K., 1999. Molecular control of locus coeruleus neurotransmission. Biol. Psychiatry 46, 1131-1139.

Ohliger-Frerking, P., Horwitz, B.A., Horowitz, J.M., 2003. Serotonergic dorsal raphe neurons from obese zucker rats are hyperexcitable. Neurosci. 120, 627-634.

Overton, P., Clark., D.,1992. Iontophoretically administered drugs acting at the N-methyl-D-aspartate receptor modulate burst firing in A9 dopamine neurons in the rat. Synapse 10,131-140.

Pan, W.J., Osmanović, S.S., Shefner, S.A., 1994. Adenosine decreases action potential duration by modulation of A-current in rat locus coeruleus neurons. J. Neurosci. 14, 1114-1122.

Penington, N.J., Kelly, J.S., Fox, A.P., 1991. A study of the mechanism of Ca2+ current inhibition produced by serotonin in rat dorsal raphe neurons. J. Neurosci. I7, 3594-3609.

Putnam, R., Quintero, M., Santin, J. et al., 2014. Computational modeling of the effects of temperature on chemosensitive locus coeruleus neurons from bullfrogs. Faseb J. 28, Supp. 1128.3. (Abstract only).

Putzier, I., Kullmann, P.H.M., Horn, J.P., Levitan, E.S., 2009. Ca1.3 channel voltage dependence, not Ca selectivity, drives pacemaker activity and amplifies bursts in nigral dopamine neurons. J. Neurosci. 29, 15414-15419.

Rall, W., 1962. Theory of physiological properties of dendrites. Ann. NY Acad. Sci. 96: 1071-1092.

Ramirez, J-M., Koch, H., Garcia, A.J. III, et al., 2011. The role of spiking and bursting pacemakers in the neuronal control of breathing. J Biol. Phys. 37, 241-261.

Sanchez-Padilla, J., Guzman, J.N., Ilijic, E. et al., 2014. Mitochondrial oxidant stress in locus coeruleus is regulated by activity and nitric oxide synthase. Nat. Neurosci. 17, 832-842.

Svensson, T.H., Engberg, C-S., Tung, C-S., Grenhoff, J., 1989. Pacemaker-like firing of noradrenergic locus coeruleus neurons in vivo induced by the excitatory amino acid antagonist kynurenate in the rat. Acta Physiol. Scand. 135, 421-422.

Sugiyama, D., Hur, S.W., Pickering, A:E. et al. 2012. In vivo patch-clamp recording from locus coeruleus neurones in the rat brainstem. J Physiol. 590, 2225-2231.

Swanson, L.W., 1976. The locus coeruleus: a cytoarchitectonic, golgi and immunohistochemical study in the albino rat. Brain Res. 110, 39-56.

Tucker, K.R., Huertas, M.A., Horn, J.P., Canavier, C.C., Levitan, E.S., 2012. Pacemaker rate and depolarization block in nigral dopamine neurons: a somatic sodium channel balancing act. J. Neurosci 32, 14519-14531.

Tuckwell, H.C., 1988. Introduction to Theoretical Neurobiology. Cambridge University Press, Cambridge UK.

Tuckwell, H.C., 2013. Biophysical properties and computational modeling of calcium spikes in serotonergic neurons of the dorsal raphe nucleus. BioSystems 112, 204-213.

Tuckwell, H.C., 2017. Computational modeling of spike generation in locus coeruleus noradrenergic neurons. Preprint.

Tuckwell, H.C., Ditlevsen, S., 2016. The space-clamped Hodgkin-Huxley system with random synaptic input: inhibition of spiking by weak noise and analysis with moment equations. Neural Computation 28, 2129-2161.

Tuckwell, H.C., Jost, J., GutKin, B.S., 2009. Inhibition and modulation of rhythmic neuronal spiking by noise. Phys. Rev. E. 80, 031907.

Tuckwell, H.C., Penington, N.J., 2014. Computational modeling of spike generation in serotonergic neurons of the dorsal raphe nucleus. Prog. Neurobiol. 118, 59-101.

Valentino, R.J., Foote, S.L., Page, M.E., 1993. The locus coeruleus as a site for integrating corticotropin-releasing factor and noradrenergic mediation of stress responses. Ann. NY.Acad. Sci 697, 173-188.

Vandermaelen, C.P., Aghajanian, G.K., 1983. Electrophysiological and pharmacological characterization of serotonergic dorsal raphe neurons recorded extracellularly and intracellularly in rat brain slices. Brain Res. 289, 109–119.

Vasudeva, R.K., Waterhouse, B.D., 2014. Cellular profile of the dorsal raphe lateral wing sub-region: relationship to the lateral dorsal tegmental nucleus. J. Chem. Neuroanat. 57-58, 15-23.

Vertes, R.P., Crane, A.M., 1997. Distribution, quantification, and morphological characteristics of serotonin-immunoreactive cells of the supralemniscal nucleus (b9) and pontomesencephalic reticular formation in the rat. J. Comp. Neurol. 378, 411-424.

Walsh, J.B., Tuckwell, H.C., 1985. Determination of the electrical potential over dendritic trees by mapping onto a nerve cylinder. J. Theor. Neurobiol. 4, 27-46.

Weinstein, J.J., Chohan, M.O., Slifstein, M., Kegeles, L.S., Moore, H., Abi-Dargham, A., 2017. Pathway-specific dopamine abnormalities in schizophrenia. Biological psychiatry, 81(1), pp.31-42.

Williams, J.T., Egan, T.M., North, R.A., 1982. Enkephalin opens potassium channels on mammalian central neurons. Nature 299, 74-77.

Williams, J.T., North, R.A., Shefner S.A. et al., 1984. Membrane properties of rat locus coeruleus neurones. Neuroscience 13, 137-156.

Wong-Lin, K-F., Joshi, A., Prasad, G., McGinnity, T.M., 2012. Network properties of a computational model of the dorsal raphe nucleus. Neural Netw. 32, 15-25.

Comments 0
Request Comment
You are adding the first comment!
How to quickly get a good reply:
  • Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
  • Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
  • Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters
Add comment
Loading ...
This is a comment super asjknd jkasnjk adsnkj
The feedback must be of minumum 40 characters
The feedback must be of minumum 40 characters

You are asking your first question!
How to quickly get a good answer:
  • Keep your question short and to the point
  • Check for grammar or spelling errors.
  • Phrase it like a question
Test description