- 1Medical Big Data Research Center, Northwest University, Xi'an, China
- 2Department of Neurology, Massachusetts General Hospital, Boston, MA, United States
It has been suggested that cholinergic neurons shape the oscillatory activity of the thalamocortical (TC) network in behavioral and electrophysiological experiments. However, theoretical modeling demonstrating how cholinergic neuromodulation of thalamocortical rhythms during non-rapid eye movement (NREM) sleep might occur has been lacking. In this paper, we first develop a novel computational model (TC-ACH) by incorporating a cholinergic neuron population (CH) into the classical thalamo-cortical circuitry, where connections between populations are modeled in accordance with existing knowledge. The neurotransmitter acetylcholine (ACH) released by neurons in CH, which is able to change the discharge activity of thalamocortical neurons, is the primary focus of our work. Simulation results with our TC-ACH model reveal that the cholinergic projection activity is a key factor in modulating oscillation patterns in three ways: (1) transitions between different patterns of thalamocortical oscillations are dramatically modulated through diverse projection pathways; (2) the model expresses a stable spindle oscillation state with certain parameter settings for the cholinergic projection from CH to thalamus, and more spindles appear when the strength of cholinergic input from CH to thalamocortical neurons increases; (3) the duration of oscillation patterns during NREM sleep including K-complexes, spindles, and slow oscillations is longer when cholinergic input from CH to thalamocortical neurons becomes stronger. Our modeling results provide insights into the mechanisms by which the sleep state is controlled, and provide a theoretical basis for future experimental and clinical studies.
1. Introduction
Sleep plays a pivotal role in mental and physical health. During sleep, the brain alternates between two stages, rapid-eye movement (REM) and NREM. This alternation is reflected in the electrical rhythms generated by the thalamocortical system, evident in the electroencephalogram (EEG). Specifically, the thalamocortical rhythms in NREM sleep include K-complexes and slow oscillations (dominated by low-frequency [0.5, 2]Hz, high amplitude oscillations), as well as spindles (characterized by a waxing and waning waveform in the range [11, 16]Hz). In contrast, REM sleep exhibits low amplitude activity of higher frequency, which is similar with rhythms in wakefulness, including the alpha oscillation (Rasch and Born, 2013). Studies reveal that spindles and slow oscillations are very helpful in protecting and stabilizing sleep (Roth, 2009; Dang-Vu et al., 2010; Kim et al., 2012).
A number of results in behavioral and electrophysiological experiments suggest that cholinergic neurons, which are distributed in the basal forebrain and brainstem, help shape the oscillatory activity of the thalamocortical network (Bellingham and Funk, 2000; Kobayashi et al., 2003; Boutrel and Koob, 2004; Brown et al., 2012). Meanwhile, it has also been shown that thalamocortical processing is subject to the action of modulatory of acetylcholine (ACH) released by cholinergic neurons (McCormick, 1992). However, a theoretical framework for understanding the mechanisms of cholinergic neuromodulation on thalamocortical rhythms during NREM sleep is still lacking.
Neural mass models (NMMs), which describe the dynamics of neural populations (Wilson and Cowan, 1973; Lopes da Silva et al., 1974), have shown success in clarifying the mechanisms involved in the generation and regulation of oscillation patterns of the “wakefulness-sleep cycle.” For example, Robinson et al. (2002) developed a continuum model of large-scale brain electrical activity to understand the transition from the resting background state to the spike-wave state. Suffczynski et al. (2004) proposed a thalamocortical NMM to explain the relation between mechanisms that generate sleep spindles and those for that generate spike-wave activity. More recently, Weigenand et al. (2014) extended a cortical NMM to explore the mechanisms for generating K-complexes and slow oscillations. Cona et al. (2014) presented a new NMM to describe the sleeping thalamocortical system, where thalamic neurons exhibit two firing modes: bursting and tonic.
Although prior (NMM-based) studies have made progress in explaining sleep rhythms, the role of ACH modulation of thalamocortical activities, especially during sleep has not been thoroughly explored. In this work we model the details of cholinergic neuromodulation in the cortex and thalamus, and study the effects of cholinergic projections on thalamocortical rhythms during NREM sleep.
In this study, a new computational model (TC-ACH) is first proposed by integrating a cholinergic neuron population (CH) into classical thalamo-cortical circuitry, with connections between populations in TC-ACH constructed on the basis of the existing mechanisms. We concentrate on the neurotransmitter ACH released by cholinergic neurons, which impacts firing activity of thalamocortical neurons. The proposed model TC-ACH is used to investigate effects of ACH modulation on thalamocortical rhythms during NREM sleep in two ways: (1) measurement of effects of cholinergic modulation in the thalamus; (2) measurement of effects of cholinergic modulation in the thalamocortical system.
The rest of the paper is organized as follows. In section 2, the model framework and mathematical expression of TC-ACH are introduced. Section 3 presents the simulation results. A summary and discussion of our work are given in the last section.
2. Methods
2.1. Model Framework
This subsection introduces the topological structure and computational operators of TC-ACH systematically.
2.1.1. Topological Structure
Our computational model TC-ACH is composed of classical thalamo-cortical circuitry as well as a cholinergic neuron population, whose topological structure is illustrated in Figure 1.
The thalamo-cortical circuitry consists of two mutually interconnected modules: a cortical module and a thalamic module. Each module comprises two neural populations, where pyramidal cells (PY) and inhibitory interneurons (IN) are included in the cortical module; the thalamic module also includes thalamocortical cells (TC) and thalamic reticular cells (RE). Regarding connections between modules, excitatory neurons (PY or TC) in each module connect to both populations in the other module. Within each module, two populations project to each other, and there exist self-connections within PY, IN, and RE (Costa et al., 2016). Note that in the thalamo-cortical circuitry, the excitatory and inhibitory projections are mediated by AMPA and GABA receptors, respectively.
The cholinergic neuron population (CH) is located in the basal forebrain and brainstem. The neurons in CH release neurotransmitter (ACH) into the synaptic cleft, which changes the activity of thalamocortical neurons. In our proposed model, there exist three modulated cholinergic projections from CH to PY, RE, as well as the connection between TC and PY (represented by red dashed lines in Figure 1), and one mediated cholinergic projection from CH to RE (represented by blue the solid line in Figure 1). The coupling coefficients of these projections are denoted by KPY, , , KTP, and , respectively. Here, all projections involve by ACH.
2.1.2. Computational Operators
The dynamical evolution in each population is implemented in two computational blocks.
The first computational block transforms the average membrane potential V(t) into the average firing rate Q(t), formulated by a sigmoid function with the form (Jansen et al., 1993)
where Qmax, θ, σ represent the maximal firing rate, the firing threshold and the neural gain, respectively.
In the second computational block, the firing rate Q(t) is first transformed into the fraction of open channels rξ(t) by a convolution with an alpha function hξ(t), that is,
Then the membrane potential V(t) is solved by
which is similar to the classical conductance-based form of Hodgkin and Huxley (1952) with one leak and two synaptic currents (say, ξ ∈ {AMPA, GABA}) (Weigenand et al., 2014). Moreover, for the convenience of calculating the convolution ⊗, rξ(t) can be equivalently obtained by solving the following second-order differential equation
In (2)-(5), Eξ and EL denote the reversal potential of the synaptic current and leak current, respectively, represents the synaptic input rate that scales rξ and Eξ, gL is the maximal conductivity of the leak current conductance, τ is the membrane time constant, N stands for the connectivity constant, and γξ is the rate constant of synaptic response. The detailed dynamical evolution in each population is illustrated in Figure 2.
2.2. Mathematical Expression
This subsection formulates the proposed model TC-ACH mathematically in the logic of five populations included in the model.
2.2.1. Cholinergic Neuron Population (CH)
We first consider the cholinergic neuronal population, which is an essential part of TC-ACH in studying the effects of cholinergic modulation on thalamocortical rhythms during NREM sleep. Due to the special role of ACH, here, we apply the concentration of ACH ([ACH]) rather than the firing rate applied in Equation (1) to complete the transformation from the membrane potential, with Destexhe et al. (1994)
where [ACH]max is the maximum concentration of ACH.
Moreover, CH receives an external input current Iext. Therefore, the membrane potential Vc(t) in CH obeys
here, Iext has the form (Rudolph et al., 2004)
where g0 is the average conductance, τext is the time constant, σext is the noise standard deviation (SD) value, and ϕC(t) denotes the independent Gaussian white noise process of unit SD and zero mean.
2.2.2. Populations PY and IN in the Cortical Module
According to Clearwater et al. (2008), we know that ACH acts on PY both pre-synaptically and post-synaptically over time scales from milliseconds to minutes.
In the post-synaptic case, ACH acts on PY through muscarinic receptors, which can cause a certain suppression of voltage-gated potassium channels (McCormick, 1992). Hence, the M-type potassium current generated by the inactivation of these channels will change the intrinsic excitability of PY. Motivated by this observation, the membrane potential Vp(t) in PY is then formulated by
where are calculated by Equations (2)–(4). The newly added current is described as follows (Rich et al., 2018)
where z(t) represents the unitless gating variable of the ionic current conductance satisfying
Here, stands for the voltage-sensitive steady-state activation function and τz is the time constant. Note that gM in Equation (11) is not a constant anymore, but varies in term of [ACH], that is,
where is a constant with the nominal value listed in Table 1.
In the pre-synaptic case, ACH acts on PY to modulate the properties of thalamocortical synapses via nicotinic receptor, which can result in changing the thalamocortical synaptic connection strength from TC to PY (denoted by Ntp) (Gil et al., 1997). Furthermore, as summarized by Kimura (2000), the thalamocortical connection strength will increase when the value of [ACH] becomes larger. Hence, we apply the way in Clearwater et al. (2008) to reformulate the connection strength N in Equation (5) as
Obviously, there have
where and N* represent the connectivity values under two extreme situations, and μ is a concentration constant.
Because in the cortical module IN contributes less than PY to the effects of cholinergic modulation on thalamocortical rhythms (Picciotto et al., 2012), in this paper, the dynamic evolution in IN is assumed to be not affected by [ACH]. That is to say, the mathematical formulation of IN obeys Equations (1)–(5).
2.2.3. Populations RE and TC in Thalamic Module
We see in Figure 1 that population RE appears to be a central hub in connecting the cortex and population CH. The electrical activity occurring in RE is important for regulating information transmission in the thalamocortical system and the shape of thalamocortical rhythms. Therefore, we first model the discharge activity of RE in this section.
On one hand, recent research indicates that thalamic reticular neurons can be directly excited by activating α7-containing nicotinic ACH receptor (Ni et al., 2016). Consequently, we introduce a new current IACH−n as an input of RE from CH (corresponding the blue solid line in Figure 1).
On the other hand, it has been shown that the potassium leak current and T-type calcium current are essential for the generation of spindle oscillation (Langdon et al., 2012). Based on this fact, two more currents ILK and IT are considered in the formulation of the membrane potential Vr(t) in RE simultaneously.
Hence, the final equation is formulated by
In Equation (15), , and are defined, respectively, as
where ω(t) is solved by
Here, ω denotes the proportion of open ion-channels caused by binding of ACH. gACH−n and EACH−n stand for the conductance and reversal potential, respectively. Pd and h represent the apparent dissociation constant and Hill coefficient of ACH binding to the receptors, and ωmax is the maximum conductance. m(t) and h(t) are activation and inactivation functions of T-type current. Details about Equations (16)–(19) can be seen in Baran et al. (2010), Sethuramanujam et al. (2016), and Destexhe et al. (1993). In addition, the calcium conductance in Equation (19) is recognized to be important for generating bursting oscillations in RE, whose value increases with increasing [ACH] (Fisher and Johnston, 1990). Specifically, it is expressed by
where and a are positive constants (Omori and Horiguchi, 2004).
In addition, it has been found that the GABAergic projection within RE is decreased by increasing [ACH] (Fisher and Johnston, 1990). This mechanism is then absorbed into defining the self-feedback connectivity of RE in our work and formulated by
where , b and c are positive constants (Omori and Horiguchi, 2004).
Next, we turn to modeling the discharge activity of TC. Different from RE, there exists one more h-type current in TC, which is responsible for the waxing and waning structure of spindle rhythms in thalamus. Hence, the membrane potential Vt(t) in TC is formulated by
where are calculated by Equations (18)–(19), and Ih is described as Destexhe et al. (1996)
Here, gh and Eh represent the conductance and reversal potential, ginc is the conductivity scaling. The details of functions mh1(t) and mh2(t) can be found in Destexhe et al. (1993).
2.2.4. Full Mathematical Expression of TC-ACH
There is, to emphasize, one key point in our model. On account of the long range afferent, there exist conduction delays between thalamus and cortex module. According to Costa et al. (2016), this delay is approximated as a convolution with the alpha function h(t). In this case, Equation (5) is then reformulated by
during four transmissions “PY → RE,” “PY → TC,” “TC → PY,” and “TC → IN.”
The full mathematical expression of TC-ACH is then formulated by Equations (22)–(36) (see page 5). Note that the model output is Vp, which can be viewed as the simulated EEG signals.
Remark 1. The description and nominal values of all parameters included in the model TC-ACH are listed in Table 1, whose values are the nominal ones reported in Costa et al. (2016), Omori and Horiguchi (2004), Clearwater et al. (2008), Rich et al. (2018), Baran et al. (2010), Sethuramanujam et al. (2016), Bhattacharya et al. (2012), and Rudolph et al. (2004).
3. Results
In this section, we apply the proposed model TC-ACH to verify the effects of cholinergic modulation on thalamocortical rhythms during NREM in two ways: (1) measurement of effects of cholinergic modulation in the thalamus; (2) measurement of effects of cholinergic modulation in the thalamocortical system.
All of the numerical simulations are performed in MATLAB R2017b (MathWorks,USA), using a stochastic Runge-Kutta method of 4th order (Rößler, 2010) with a step size of 0.1ms.
3.1. The Oscillatory Phenomena Due to Cholinergic Modulation in Thalamus
The thalamus is believed to act as a “pacemaker” for thalamocortical rhythms, and is able to independently generate multiple brain rhythms during sleep (Hughes and Crunelli, 2005; Li et al., 2017). Hence, oscillatory phenomena due to cholinergic modulation in thalamus are first explored, where only the thalamic module and projections from CH to RE are considered. Here, the model output is Vt.
We first show different thalamic oscillations and transitions between them caused by cholinergic projections. To this end, we apply one-dimensional bifurcation analysis for several key parameters, represented by gLK, , , and . The bifurcation diagram is obtained by plotting the stable local minimum, as well as top three maximum values of Vt over changes in each of four parameters. All simulations are executed for 40 s and those minimum and maximum values are obtained from the latter stable 30 s of the time series.
Figures 3A–D shows bifurcation diagrams of Vt over changes in gLK, , , and , respectively. Seen from Figure 3A, it reveals that there are four different dynamical states including low firing (I), spindle (II), fast oscillation (III) and slow oscillation (IV) as gLK varies in [0.018, 0.055]. Specifically, when the value of gLK is extremely small, the model exhibits a low firing state and no oscillation behavior can be observed. As gLK becomes a little larger, the model mainly experiences a spindle oscillation pattern for a period, in which multiple pairs of maximum and minimum values are found within each periodic complex. With further growth of gLK, the model exhibits a fast oscillation pattern, in which only one pair of maximum and minimum values emerges within each periodic complex. When strong hyperpolarization through increasing gLK attains, the model moves from a high frequency oscillation pattern into a slow oscillation pattern, in which multiple pairs of maximum and minimum values can also be observed within each periodic complex, but the distances between the top maximum and other maximum values become larger. Figure 3E illustrates the obtained thalamic oscillations, which correspond to four dynamical states (gLK = 0.018, 0.03, 0.045, 0.052), respectively.
![www.frontiersin.org](https://www.frontiersin.org/files/Articles/503660/fncom-13-00100-HTML-r1/image_m/fncom-13-00100-g003.jpg)
Figure 3. The bifurcation diagrams of Vt over changes in gLK, , , and . (A) Bifurcation diagram (varying gLK). (B) Bifurcation diagram (varying ). (C) Bifurcation diagram (varying ). (D) Bifurcation diagram (varying ). (E) Four thalamic oscillations with respect to different values of gLK.
A similar explanation can be given for the cases of , , and , while only two dynamical states (i.e., II and III) are obtained. Moreover, the state transition of is different from that of and . It is apparent that the model has a pattern transition from state III to an increasingly stable state II as increases (see Figure 3B). By contrast, a transition from state II to state III with the increase of (or ) can be observed in Figure 3C (or Figure 3D). That is to say, the model appears to demonstrate a stabler spindle oscillation state with the increasing value of , while the growth of (or ) leads to reduced spindle rhythms. Prior work has established that spindle rhythms are helpful to protect sleep (Dang-Vu et al., 2010; Kim et al., 2012). Therefore, the obtained results implicate that cholinergic modulation in RE may help promote sleep (or arousal) states.
Next, we check whether our results can be generalized within a certain range of parameters. The above bifurcation analysis allows us to further distinguish different dynamical state regions in the two-parameter space (for example, see Figure 5A). Moreover, the power spectral analysis is applied to estimate the dominant frequency (Df) and second dominant frequency (Sec-Df) of neural oscillations from the time series of Vt. Figure 4 illustrates the box-plots of Df and Sec-Df with respect to three oscillations (fast, slow and spindle). It is clear from Figure 4 that the Sec-Df has better differentiation capability than Df. Therefore, the corresponding Sec-Df regions are drawn for each pair of parameters (for example, see Figure 5B) in the following analysis.
![www.frontiersin.org](https://www.frontiersin.org/files/Articles/503660/fncom-13-00100-HTML-r1/image_m/fncom-13-00100-g004.jpg)
Figure 4. The box-plots of (A) Df and (B) Sec-Df with respect to three oscillations (fast, slow, and spindle).
Figures 5A,B illustrate the dynamical state regions and Sec-Df regions with 41 × 41 grids in the space of . As shown in Figure 5A, four different state regions are displayed, whose identification is same as above [i.e., low firing (I), spindle (II), fast oscillation (III) and slow oscillation (IV)]. It can be observed that along with rising, the model generates more spindle oscillation pattern after lasting a certain fast oscillation period (see the black arrow). By combining the results of frequency analysis shown in Figure 5B, we can outline the spindle oscillation region that falls into the 11–14 Hz frequency range. It should be noted that, the obtained result is consistent with the frequency range [11, 16]Hz characteristic of sleep spindles.
Similar results are obtained in the spaces of and , which are illustrated in Figures 5C–F. We observe that increasing promotes the generation of spindle oscillations in our model (see the black solid arrow in Figure 5A), while increasing (or ) suppresses generation of spindle oscillations (see the black dash arrow in Figures 5C,E). These results are in line with the results mentioned above.
Additionally, in order to demonstrate the interactions among parameters , and , the dynamical state regions and Sec-Df regions for each pairwise combination of them are shown in Figure 6. As expected, both state and Sec-Df analysis in three panels (, ), (, ), and (, ) provide the same evidence as the above description, that is, the spindle oscillation is generated as increases (or as , decreases).
![www.frontiersin.org](https://www.frontiersin.org/files/Articles/503660/fncom-13-00100-HTML-r1/image_m/fncom-13-00100-g006.jpg)
Figure 6. The dynamical state regions (A,C,E) and Sec-Df regions (B,D,F) in three panels (, ), (, ), and (, .
3.2. The Oscillatory Phenomena Due to Cholinergic Modulation in the Thalamocortical System
Here we concentrate on the effects of cholinergic modulation on thalamocortical rhythms during NREM sleep, including K-complexes, spindles, and slow oscillations. We fix the values of the four parameters to be gLK = 0.034, , , and . Based on the results obtained in the proceeding subsection, we know that spindle oscillations can be generated in the thalamus module with such settings.
We first show the dynamical behavior of our TC-ACH model with variation of KPY, which is the strength of the cholinergic projection from CH to PY. In this case, we assume that there is no cholinergic projection from CH to “PY-TC connection” (KTP = 0), and we apply linear ramps to increase or decrease the parameter value and observe the effects of continuous changes in KPY on the model. Figure 7B illustrates the model output Vp as KPY varies in [2.3, 5] linearly (see Figure 7A). In the beginning, when KPY = 2.3, the model exhibits K-complex oscillations. Then we can find that a small ascent of KPY drives a transition from K-complexes to slow oscillations. Furthermore, with the gradual increase of KPY, the slow oscillation is kept for a short period until KPY reaches its higher threshold. At this time, the slow oscillation is replaced by the α-like oscillation. As KPY ramps down, the model transitions back into its original state gradually (that is, from α-like to slow oscillations, and then from slow oscillations to K-complexes). The enlarged graphs of three different oscillations in 20 s are displayed in Figures 7E–G.
![www.frontiersin.org](https://www.frontiersin.org/files/Articles/503660/fncom-13-00100-HTML-r1/image_m/fncom-13-00100-g007.jpg)
Figure 7. The dynamical behavior of model TC-ACH with the variation of KPY: (A) the variation of KPY over time t; (B) the model output Vp; (C) the spectrogram of Vp; (D) the sample entropy of Vp; (E) the enlarged model output in [10,30]s (including K-complex); (F) the enlarged model output in [30,50]s (including slow oscillation); (G) the enlarged model output in [65,85]s (including α-like oscillation).
On the basis of the above observations, we conclude that K-complexes and slow oscillations can be triggered by weakening the strength of the cholinergic projection from CH to cortex, while the emergence of α-like activity, characteristic of REM, requires stronger cholinergic input. Here, if we roughly consider that K-complexes and slow oscillations mainly emerge during NREM sleep and α-like rhythms mostly appear during REM sleep, the obtained results in our work are consistent with the conclusion in Lena et al. (2005), that is, the concentration of ACH is lowest during NREM sleep and highest during wake and REM sleep.
The spectral analysis and the sample entropy extracted from the model output Vp also support this conclusion from another point of view. Figure 7C illustrates the spectrogram of Vp. During the first and last 40 s where K-complexes and slow oscillations are emerging, the spectral power of Vp mainly falls within the [0, 5]Hz range; while in the middle stage where α-like rhythms begin to appear, it increases to [0, 10]Hz. These results obviously conform with the frequency ranges of certain sleep rhythms. Besides, the sample entropy extracted from Vp also shows a similar tendency. In detail, it is larger in the α-like stage, but smaller in the K-complex and slow wave oscillation stage (see Figure 7D). Sample entropy measures signal complexity, hence this observation is consistent with the conclusion that the complexity of cortical rhythms are decreased as the deeper sleep state, but increased during REM and wake (Bruce et al., 2009; See and Liang, 2011).
Next, we show the dynamical behavior of constructed model TC-ACH with variation of KTP, which is the strength of the cholinergic projection from CH to the connection between PY and TC. Here, KPY shows the trend as in Figure 7, and KTP is set to be 1.5, 3.5, 5.5, respectively, in simulations. Figure 8 illustrates the model output Vp under three different cases.
![www.frontiersin.org](https://www.frontiersin.org/files/Articles/503660/fncom-13-00100-HTML-r1/image_m/fncom-13-00100-g008.jpg)
Figure 8. The dynamical behavior of model TC-ACH in the cases where KTP = 1.5, KTP = 3.5, KTP = 5.5. (A1,B1,C1) The model output Vp; (A2,B2,C2) the enlarged model output in [20,40]s; (A3,B3,C3) the enlarged model output in [40,60]s; (A4,B4,C4) the enlarged model output in [65,85]s.
In the first case where the connectivity between TC and PY is regulated by a lower strength of cholinergic projection from CH (i.e., KTP = 1.5), we observe in Figure 8A1 that some weaker spindle rhythms emerge when the value of KPY is relatively small. It can be observed more clearly in Figure 8A2 where the enlarged graph exhibits spindle rhythms (with red circle) and K-complex oscillations in [20, 40] s. Meanwhile, in Figures 8A3,A4, we see that the duration of slow oscillation is much longer, and α-like activity is much less by comparing the results shown in Figure 7.
This conclusion can be further verified under other two cases. It can be seen in Figures 8B1–B4 that more spindle rhythms appear as KPY gradually increases and the α-like activities have been largely replaced by slow waves in the case where KTP = 3.5. Furthermore, when KTP is set to be 5.5, the spindle rhythms run through the model output from the beginning to end. We can see from Figure 8C4 that there are still spindle rhythms when KPY attains its maximum value. These observations reveal that the strength of cholinergic projections from CN to “PY-TC connection” play an important role in promoting the spindle rhythms and prolonging the duration of the NREM state in the thalamocortical system.
4. Summary and Discussion
In this paper, we first proposed a novel computational model (TC-ACH) by integrating a neuron population CH into classical thalamo-cortical circuitry (including populations PY, IN, TC, and RE). The connections between five populations are built in accordance with the established mechanisms. Our model considers the neurotransmitter ACH released by neurons in CH, which alters discharge activities of thalamocortical neurons. For simplicity, we represent these additions to the classic model by four cholinergic projections, where coupling coefficients represent different projections from CH to thalamus, and KPY and KTP represent the projections from CH to cortex and thalamocortical system, respectively. On the basis of established model framework of TC-ACH, the corresponding mathematical expression has been formulated in the logic of five populations systematically, where the average membrane potential V(t) is solved to simulate the rhythms generated by each population.
Next, we applied the developed model TC-ACH to study the effects of ACH modulation on thalamocortical rhythms during NREM sleep in two ways:
(1) Measurement of effects of cholinergic modulation in the thalamus. In this case, only the thalamic module and projections from CH to RE were considered. Simulation results suggest that cholinergic projection activity is a key factor in modulating oscillation patterns in the thalamic module. Specifically, the model appears to be a stabler spindle oscillation state with the increasing value of , while the growth of (or ) leads to reduced spindle rhythms. Moreover, with variation of the potassium leak conductance gLK, which is dramatically modulated by the concentration of ACH (McCormick, 1989, 1992), four different dynamical states including the low firing, spindle, fast oscillation, and slow oscillation can be obtained and transited.
(2) Measurement of effects of cholinergic modulation in the thalamocortical system. In this case, the dynamical behavior of our TC-ACH model was studied by varying KPY and KTP, respectively. Simulation results show that the K-complex and slow oscillations can be triggered by weakening strength of KPY, while emergence of α-like activity requires stronger input. Furthermore, when there exists cholinergic input from CH to “PY-TC connection” (that is, KTP ≠ 0), we found that the duration of oscillation patterns during NREM sleep including K-complexes, spindles and slow oscillations is longer. Additionally, when KTP is relatively larger, more spindle rhythms appear and α-like activities are largely replaced by slow waves.
It should be noted that a number of electrophysiological experiments have investigated cholinergic modulation of cortex, thalamus or thalamocortical system, respectively (McCormick and Prince, 1986; McCormick, 1989, 1992; Clarke, 2004; Mesulam, 2004; Hasselmo and Giocomo, 2006; Beierlein, 2014). However, few studies correlate ACH modulation with rhythmic activities (Steriade et al., 1993; Steriade, 2004), where the model-based work to study in theory the effects of cholinergic modulation on thalamocortical rhythms during sleep is far less. More than that, what few existing studies are not comprehensive enough in studying the effects of ACH modulation on thalamocortical rhythms, especially during sleep. For example, the effects of ACH modulation are considered only on cortex or thalamus separately, but not on the whole thalamocortical system (Omori and Horiguchi, 2004; Li et al., 2017); the mechanism regarding ACH modulation is only modeled as a certain parameter, but not as a whole neuron population (Li et al., 2017); the model is constructed at the microscopic level, which cannot relate directly to thalamocortical rhythms at the mesoscopic level (Omori and Horiguchi, 2004). Therefore, in order to overcome such limitations, we have constructed a novel computational model (TC-ACH) by incorporating a cholinergic neuron population into the classical thalamo-cortical circuitry at the mesoscopic level. By thus doing, a deeper understanding of the role of cholinergic modulation on thalamocortical system will be got, and further, a critical insight into the mechanisms controlling sleep state may be found.
Besides five projections applied in the modeling (represented by , KPY, and KTP), there still exist cholinergic modulators on other projections, such as corticoreticular and corticothalamic projections. Castro-Alamancos and Calcagnotto (2001) demonstrated that the corticothalamic activity can be filtered by cholinergic activation during arousal high-pass according to the experiments performed in vitro and in vivo. In addition, as indicated in Itier and Bertrand (2002), the cholinergic modulation on corticoreticular projection may induce the generation of more complex brain rhythms (such as spike and wave) during sleep. However, the quantitative description relating such cholinergic mechanisms to the generation and transition of typical rhythms during NREM sleep is very limited. Therefore, we hope these mechanisms could be further considered with more attempts, to study the cholinergic modulation of thalamocortical rhythms during sleep in different cases.
Another point to stress here is the synaptic connection ways in intra-RE. Previous electrophysiological studies indicated that the thalamic reticular neurons are functionally connected through chemical/electrical synapses (Sanchez-Vives et al., 1997; Landisman et al., 2002; Shu and McCormick, 2002; Long et al., 2004; Deleuze and Huguenard, 2006; Lam et al., 2006). A potentially paradigm-shifting question has been presented as to whether chemical synapses between thalamic reticular neurons are altogether absent in certain mammals or degenerate as a function of increasing age (Landisman et al., 2002; Cruikshank et al., 2010; Hou et al., 2016). However, there still have other works supporting the existence of GABAergic intrareticular synapses. For example, one of the most recent model-based study shows that intrareticular synapses, both chemical and electrical, manifest certain effects on the signal propagation and oscillation (Brown et al., 2019). In our modeling work, by means of considering the chemical synapses in intra-RE population, the measurement of effects of cholinergic modulation in thalamocortical system has been completed well. In contrast, the electrical synapses do not work effectively under the same situation. Consequently, we roughly hypothesize that GABAergic synapses may contribute much more than electrical synapses to the effects of cholinergic modulation in thalamocortical rhythms, at least at the mesoscopic level. But on the other hand, we certainly acknowledge that the electrical intrareticular synapses play important roles vis-à-vis thalamic signaling. Therefore, an important extension of the current work would be the modeling of electrical synapses between reticular neurons from the microscopic point of view, such as the model-based works in Pham and Haas (2018) and Brown et al. (2019).
An important limitation of our TC-ACH model is that it considers only ACH modulation. However, it is known that thalamocortical rhythms during sleep are also directly affected by other neuromodulators, such as noradrenalin (NE), serotonin (5-HT), histamine (HA), and dopamine (DA) from the hypothalamus and brainstem, whose concentrations vary over the night (Lena et al., 2005). Therefore, some apparent questions need to be answered: whether our method can be further developed to shed light on other sleep-related neuromodulators? If so, how to model the corresponding mechanisms and explore their effects on thalamocortical rhythms during sleep? Fortunately, based on the progressive mathematical description of sleep regulatory networks (Kumar et al., 2012; Booth et al., 2017), it may well be possible to carry our model further by constructing a new thalamocortical NMM which contains various sleep-related neuron populations (releasing NE, 5-HT, HA, DA); this is a topic we plan to pursue in future work. In addition, such model-based research can make a contribution to understanding sleep related pathological conditions, such as sleep-related epilepsy.
Data Availability Statement
The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation, to any qualified researcher.
Author Contributions
QL, J-LS, S-HL, MW, and RZ designed and performed the research as well as wrote the paper.
Funding
This work was supported by the National Natural Science Foundation of China under Grant 61473223, the Innovative Talents Promotion Plan of Shaanxi Province under Grant 2018TD-016, the Key Industry Innovation Chain (Group) of Shaanxi Province under Grant 2019ZDLSF02-09-02. MW was supported during this research by the Glenn Foundation for Medical Research and the American Federation for Aging Research through a Breakthroughs in Gerontology Grant; through the American Academy of Sleep Medicine through an AASM Foundation Strategic Research Award; and by grants from the NIH (1R01NS102190, 1R01NS102574, 1R01NS107291, 1RF1AG064312).
Conflict of Interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
References
Baran, I., Iftime, A., and Popescu, A. (2010). Diffusion–convection effects on drug distribution at the cell membrane level in a patch-clamp setup. BioSystems 102, 134–147. doi: 10.1016/j.biosystems.2010.09.003
Beierlein, M. (2014). Synaptic mechanisms underlying cholinergic control of thalamic reticular nucleus neurons. J. Physiol. 592, 4137–4145. doi: 10.1113/jphysiol.2014.277376
Bellingham, M. C., and Funk, G. D. (2000). Cholinergic modulation of respiratory brain-stem neurons and its function in sleep–wake state determination. Clin. Exp. Pharmacol. Physiol. 27, 132–137. doi: 10.1046/j.1440-1681.2000.03192.x
Bhattacharya, B. S., Coyle, D., Maguire, L. P., and Stewart, J. (2012). “Kinetic modelling of synaptic functions in the alpha rhythm neural mass model,” in International Conference on Artificial Neural Networks, eds A. E. P. Villa (Berlin; Heidelberg: Springer), 645–652. doi: 10.1007/978-3-642-33269-2_81
Booth, V., Xique, I., and Diniz Behn, C. G. (2017). One-dimensional map for the circadian modulation of sleep in a sleep-wake regulatory network model for human sleep. SIAM J. Appl. Dyn. Syst. 16, 1089–1112. doi: 10.1137/16M1071328
Boutrel, B., and Koob, G. F. (2004). What keeps us awake: the neuropharmacology of stimulants and wakefulness promoting medications. Sleep 27, 1181–1194. doi: 10.1093/sleep/27.6.1181
Brown, J. W., Taheri, A., Kenyon, R. V., Berger-Wolf, T., and Llano, D. A. (2019). A computational model of intrathalamic signaling via open-loop thalamo-reticular-thalamic architectures. bioRXiv. 574178. doi: 10.1101/574178
Brown, R. E., Basheer, R., McKenna, J. T., Strecker, R. E., and McCarley, R. W. (2012). Control of sleep and wakefulness. Physiol. Rev. 92, 1087–1187. doi: 10.1152/physrev.00032.2011
Bruce, E. N., Bruce, M. C., and Vennelaganti, S. (2009). Sample entropy tracks changes in EEG power spectrum with sleep state and aging. J. Clin. Neurophysiol. 26:257. doi: 10.1097/WNP.0b013e3181b2f1e3
Castro-Alamancos, M. A., and Calcagnotto, M. E. (2001). High-pass filtering of corticothalamic activity by neuromodulators released in the thalamus during arousal: in vitro and in vivo. J. Neurophysiol. 85, 1489–1497. doi: 10.1152/jn.2001.85.4.1489
Clarke, P. B. (2004). Nicotinic modulation of thalamocortical neurotransmission. Prog. Brain Res. 145, 253–260. doi: 10.1016/S0079-6123(03)45017-6
Clearwater, J., Rennie, C., and Robinson, P. (2008). Mean field model of acetylcholine mediated dynamics in the thalamocortical system. J. Theor. Biol. 255, 287–298. doi: 10.1016/j.jtbi.2008.08.010
Cona, F., Lacanna, M., and Ursino, M. (2014). A thalamo-cortical neural mass model for the simulation of brain rhythms during sleep. J. Comput. Neurosci. 37, 125–148. doi: 10.1007/s10827-013-0493-1
Costa, M. S., Weigenand, A., Ngo, H.-V. V., Marshall, L., Born, J., Martinetz, T., et al. (2016). A thalamocortical neural mass model of the EEG during NREM sleep and its response to auditory stimulation. PLoS Comput. Biol. 12:e1005022. doi: 10.1371/journal.pcbi.1005022
Cruikshank, S. J., Urabe, H., Nurmikko, A. V., and Connors, B. W. (2010). Pathway-specific feedforward circuits between thalamus and neocortex revealed by selective optical stimulation of axons. Neuron 65, 230–245. doi: 10.1016/j.neuron.2009.12.025
Dang-Vu, T. T., McKinney, S. M., Buxton, O. M., Solet, J. M., and Ellenbogen, J. M. (2010). Spontaneous brain rhythms predict sleep stability in the face of noise. Curr. Biol. 20, R626–R627. doi: 10.1016/j.cub.2010.06.032
Deleuze, C., and Huguenard, J. R. (2006). Distinct electrical and chemical connectivity maps in the thalamic reticular nucleus: potential roles in synchronization and sensation. J. Neurosci. 26, 8633–8645. doi: 10.1523/JNEUROSCI.2333-06.2006
Destexhe, A., Babloyantz, A., and Sejnowski, T. J. (1993). Ionic mechanisms for intrinsic slow oscillations in thalamic relay neurons. Biophys. J. 65, 1538–1552. doi: 10.1016/S0006-3495(93)81190-1
Destexhe, A., Bal, T., McCormick, D. A., and Sejnowski, T. J. (1996). Ionic mechanisms underlying synchronized oscillations and propagating waves in a model of ferret thalamic slices. J. Neurophysiol. 76, 2049–2070. doi: 10.1152/jn.1996.76.3.2049
Destexhe, A., Mainen, Z. F., and Sejnowski, T. J. (1994). Synthesis of models for excitable membranes, synaptic transmission and neuromodulation using a common kinetic formalism. J. Comput. Neurosci. 1, 195–230. doi: 10.1007/BF00961734
Fisher, R., and Johnston, D. (1990). Differential modulation of single voltage-gated calcium channels by cholinergic and adrenergic agonists in adult hippocampal neurons. J. Neurophysiol. 64, 1291–1302. doi: 10.1152/jn.1990.64.4.1291
Gil, Z., Connors, B. W., and Amitai, Y. (1997). Differential regulation of neocortical synapses by neuromodulators and activity. Neuron 19, 679–686. doi: 10.1016/S0896-6273(00)80380-3
Hasselmo, M., and Giocomo, L. (2006). Cholinergic modulation of cortical function. J. Mol. Neurosci. 30, 133–135. doi: 10.1385/JMN:30:1:133
Hodgkin, A. L., and Huxley, A. F. (1952). A quantitative description of membrane current and its application to conduction and excitation in nerve. J. Physiol. 117, 500–544. doi: 10.1113/jphysiol.1952.sp004764
Hou, G., Smith, A. G., and Zhang, Z.-W. (2016). Lack of intrinsic gabaergic connections in the thalamic reticular nucleus of the mouse. J. Neurosci. 36, 7246–7252. doi: 10.1523/JNEUROSCI.0607-16.2016
Hughes, S. W., and Crunelli, V. (2005). Thalamic mechanisms of EEG alpha rhythms and their pathological implications. Neuroscientist 11, 357–372. doi: 10.1177/1073858405277450
Itier, V., and Bertrand, D. (2002). Mutations of the neuronal nicotinic acetylcholine receptors and their association with adnfle. Neurophysiol. Clin. 32, 99–107. doi: 10.1016/S0987-7053(02)00294-0
Jansen, B. H., Zouridakis, G., and Brandt, M. E. (1993). A neurophysiologically-based mathematical model of flash visual evoked potentials. Biol. Cybernet. 68, 275–283. doi: 10.1007/BF00224863
Kim, A., Latchoumane, C., Lee, S., Kim, G. B., Cheong, E., Augustine, G. J., et al. (2012). Optogenetically induced sleep spindle rhythms alter sleep architectures in mice. Proc. Natl. Acad. Sci. U.S.A. 109, 20673–20678. doi: 10.1073/pnas.1217897109
Kimura, F. (2000). Cholinergic modulation of cortical function: a hypothetical role in shifting the dynamics in cortical network. Neurosci. Res. 38, 19–26. doi: 10.1016/S0168-0102(00)00151-6
Kobayashi, Y., Inoue, Y., and Isa, T. (2003). Cognitive function of the brainstem cholinergic system. Int. Congr. Ser. 1250, 397–406. doi: 10.1016/S0531-5131(03)00967-1
Kumar, R., Bose, A., and Mallick, B. N. (2012). A mathematical model towards understanding the mechanism of neuronal regulation of wake-nrems-rems states. PLoS ONE 7:e42059. doi: 10.1371/journal.pone.0042059
Lam, Y.-W., Nelson, C. S., and Sherman, S. M. (2006). Mapping of the functional interconnections between thalamic reticular neurons using photostimulation. J. Neurophysiol. 96, 2593–2600. doi: 10.1152/jn.00555.2006
Landisman, C. E., Long, M. A., Beierlein, M., Deans, M. R., Paul, D. L., and Connors, B. W. (2002). Electrical synapses in the thalamic reticular nucleus. J. Neurosci. 22, 1002–1009. doi: 10.1523/JNEUROSCI.22-03-01002.2002
Langdon, A. J., Breakspear, M., and Coombes, S. (2012). Phase-locked cluster oscillations in periodically forced integrate-and-fire-or-burst neuronal populations. Phys. Rev. E 86:061903. doi: 10.1103/PhysRevE.86.061903
Lena, I., Parrot, S., Deschaux, O., Muffat-Joly, S., Sauvinet, V., Renaud, B., et al. (2005). Variations in extracellular levels of dopamine, noradrenaline, glutamate, and aspartate across the sleep–wake cycle in the medial prefrontal cortex and nucleus accumbens of freely moving rats. J. Neurosci. Res. 81, 891–899. doi: 10.1002/jnr.20602
Li, G., Henriquez, C. S., and Fröhlich, F. (2017). Unified thalamic model generates multiple distinct oscillations with state-dependent entrainment by stimulation. PLoS Comput. Biol. 13:e1005797. doi: 10.1371/journal.pcbi.1005797
Long, M. A., Landisman, C. E., and Connors, B. W. (2004). Small clusters of electrically coupled neurons generate synchronous rhythms in the thalamic reticular nucleus. J. Neurosci. 24, 341–349. doi: 10.1523/JNEUROSCI.3358-03.2004
Lopes da Silva, F., Hoeks, A., Smits, H., and Zetterberg, L. (1974). Model of brain rhythmic activity. Biol. Cybernet. 15, 27–37. doi: 10.1007/BF00270757
McCormick, D. A. (1989). Cholinergic and noradrenergic modulation of thalamocortical processing. Trends Neurosci. 12, 215–221. doi: 10.1016/0166-2236(89)90125-2
McCormick, D. A. (1992). Neurotransmitter actions in the thalamus and cerebral cortex and their role in neuromodulation of thalamocortical activity. Prog. Neurobiol. 39, 337–388. doi: 10.1016/0301-0082(92)90012-4
McCormick, D. A., and Prince, D. A. (1986). Acetylcholine induces burst firing in thalamic reticular neurones by activating a potassium conductance. Nature 319:402. doi: 10.1038/319402a0
Mesulam, M.-M. (2004). The cholinergic innervation of the human cerebral cortex. Prog. Brain Res. 145, 67–78. doi: 10.1016/S0079-6123(03)45004-8
Ni, K.-M., Hou, X.-J., Yang, C.-H., Dong, P., Li, Y., Zhang, Y., et al. (2016). Selectively driving cholinergic fibers optically in the thalamic reticular nucleus promotes sleep. Elife 5:e10382. doi: 10.7554/eLife.10382
Omori, T., and Horiguchi, T. (2004). Dynamical state transition by neuromodulation due to acetylcholine in neural network model for oscillatory phenomena in thalamus. J. Phys. Soc. Japan 73, 3489–3494. doi: 10.1143/JPSJ.73.3489
Pham, T., and Haas, J. S. (2018). Electrical synapses between inhibitory neurons shape the responses of principal neurons to transient inputs in the thalamus: a modeling study. Sci. Rep. 8:7763. doi: 10.1038/s41598-018-25956-x
Picciotto, M. R., Higley, M. J., and Mineur, Y. S. (2012). Acetylcholine as a neuromodulator: cholinergic signaling shapes nervous system function and behavior. Neuron 76, 116–129. doi: 10.1016/j.neuron.2012.08.036
Rasch, B., and Born, J. (2013). About sleep's role in memory. Physiol. Rev. 93, 681–766. doi: 10.1152/physrev.00032.2012
Rich, S., Zochowski, M., and Booth, V. (2018). Effects of neuromodulation on excitatory–inhibitory neural network dynamics depend on network connectivity structure. J. Nonlinear Sci. 8, 1–24. doi: 10.1007/s00332-017-9438-6
Robinson, P., Rennie, C., and Rowe, D. (2002). Dynamics of large-scale brain activity in normal arousal states and epileptic seizures. Phys. Rev. E 65:041924. doi: 10.1103/PhysRevE.65.041924
Rößler, A. (2010). Runge–kutta methods for the strong approximation of solutions of stochastic differential equations. SIAM J. Num. Anal. 48, 922–952. doi: 10.1137/09076636X
Rudolph, M., Piwkowska, Z., Badoual, M., Bal, T., and Destexhe, A. (2004). A method to estimate synaptic conductances from membrane potential fluctuations. J. Neurophysiol. 91, 2884–2896. doi: 10.1152/jn.01223.2003
Sanchez-Vives, M. V., Bal, T., and McCormick, D. A. (1997). Inhibitory interactions between perigeniculate gabaergic neurons. J. Neurosci. 17, 8894–8908. doi: 10.1523/JNEUROSCI.17-22-08894.1997
See, A. R., and Liang, C.-K. (2011). “A study on sleep eeg using sample entropy and power spectrum analysis,” in 2011 Defense Science Research Conference and Expo (DSR) (Singapore: IEEE), 1–4. doi: 10.1109/DSR.2011.6026802
Sethuramanujam, S., McLaughlin, A. J., deRosenroll, G., Hoggarth, A., Schwab, D. J., and Awatramani, G. B. (2016). A central role for mixed acetylcholine/gaba transmission in direction coding in the retina. Neuron 90, 1243–1256. doi: 10.1016/j.neuron.2016.04.041
Shu, Y., and McCormick, D. A. (2002). Inhibitory interactions between ferret thalamic reticular neurons. J. Neurophysiol. 87, 2571–2576. doi: 10.1152/jn.00850.2001
Steriade, M. (2004). Acetylcholine systems and rhythmic activities during the waking–sleep cycle. Prog. Brain Res. 145, 179–196. doi: 10.1016/S0079-6123(03)45013-9
Steriade, M., Contreras, D., Dossi, R. C., and Nunez, A. (1993). The slow (<1 hz) oscillation in reticular thalamic and thalamocortical neurons: scenario of sleep rhythm generation in interacting thalamic and neocortical networks. J. Neurosci. 13, 3284–3299. doi: 10.1523/JNEUROSCI.13-08-03284.1993
Suffczynski, P., Kalitzin, S., and Lopes Da Silva, F. (2004). Dynamics of non-convulsive epileptic phenomena modeled by a bistable neuronal network. Neuroscience 126, 467–484. doi: 10.1016/j.neuroscience.2004.03.014
Weigenand, A., Costa, M. S., Ngo, H.-V. V., Claussen, J. C., and Martinetz, T. (2014). Characterization of k-complexes and slow wave activity in a neural mass model. PLoS Comput. Biol. 10:e1003923. doi: 10.1371/journal.pcbi.1003923
Keywords: thalamocortical neural mass model, cholinergic projection, acetylcholine (ACH), thalamocortical rhythm, NREM sleep
Citation: Li Q, Song J-L, Li S-H, Westover MB and Zhang R (2020) Effects of Cholinergic Neuromodulation on Thalamocortical Rhythms During NREM Sleep: A Model Study. Front. Comput. Neurosci. 13:100. doi: 10.3389/fncom.2019.00100
Received: 09 September 2019; Accepted: 30 December 2019;
Published: 23 January 2020.
Edited by:
Julie S. Haas, Lehigh University, United StatesReviewed by:
Daniel Llano, University of Illinois at Urbana-Champaign, United StatesQing Yun Wang, Beihang University, China
Copyright © 2020 Li, Song, Li, Westover and Zhang. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Rui Zhang, rzhang@nwu.edu.cn