1. Neuroscience
Download icon

Cellular and neurochemical basis of sleep stages in the thalamocortical network

  1. Giri P Krishnan  Is a corresponding author
  2. Sylvain Chauvette
  3. Isaac Shamie
  4. Sara Soltani
  5. Igor Timofeev
  6. Sydney S Cash
  7. Eric Halgren
  8. Maxim Bazhenov
  1. University of California, San Diego, United States
  2. Université Laval, Canada
  3. Massachusetts General Hospital and Harvard Medical School, United States
Research Article
  • Cited 11
  • Views 1,745
  • Annotations
Cite this article as: eLife 2016;5:e18607 doi: 10.7554/eLife.18607

Abstract

The link between the combined action of neuromodulators in the brain and global brain states remains a mystery. In this study, using biophysically realistic models of the thalamocortical network, we identified the critical intrinsic and synaptic mechanisms, associated with the putative action of acetylcholine (ACh), GABA and monoamines, which lead to transitions between primary brain vigilance states (waking, non-rapid eye movement sleep [NREM] and REM sleep) within an ultradian cycle. Using ECoG recordings from humans and LFP recordings from cats and mice, we found that during NREM sleep the power of spindle and delta oscillations is negatively correlated in humans and positively correlated in animal recordings. We explained this discrepancy by the differences in the relative level of ACh. Overall, our study revealed the critical intrinsic and synaptic mechanisms through which different neuromodulators acting in combination result in characteristic brain EEG rhythms and transitions between sleep stages.

https://doi.org/10.7554/eLife.18607.001

eLife digest

There are several stages of sleep that cycle repeatedly through the night with each producing distinctive patterns of electrical activity in the brain. It is thought that these patterns may help us to remember things that have happened throughout the day. Cells in parts of the brain called the hypothalamus and the brainstem control transitions between sleep stages. They regulate the release of chemicals known as neuromodulators in many parts of the brain, including the cortex and thalamus, which play the roles in memory and learning. Researchers now know how the neuromodulators influence the properties of individual brain cells. However, it is not clear how coordinated action of many neuromodulators result in the patterns of electrical activity seen in the brain during each stage of sleep.

Krishnan et al. used a computer model to investigate how three of these neuromodulators – acetylcholine, histamine and GABA – shift electrical activity in the brain between sleep stages. The computer model was able to recreate the network of brain cells in the cortex and thalamus and how this network responds to the changes in the levels of neuromodulators. The study found that simultaneous and balanced changes of acetylcholine, histamine, and GABA work together to shift the brain between the stages of sleep and to initiate patterns of the brain electrical activity specific to the different sleep stages.

Krishnan et al. predict that the relative differences in the level of acetylcholine in the brains of humans, cats and mice may explain why different species have different patterns of electrical activity during sleep. The study also found that an anesthetic drug called propofol may induce sleep-like patterns of electrical activity in the human brain by affecting the levels of all three of the neuromodulators. More studies are needed to look at how the networks of cells in the cortex and thalamus communicate with the brainstem, and how changes in the levels of neuromodulators affect memory and learning.

https://doi.org/10.7554/eLife.18607.002

Introduction

During sleep, unique electrophysiological rhythms are observed in the EEG, intracellular recordings, EOG and EMG (Steriade et al., 1993a; Niedermeyer and da Silva, 2005). These features form the basis for classification of sleep into different stages: rapid eye movement sleep (REM), and stages N1 (Stage 1), N2 (Stage 2) and N3 (Stage 3) of non-REM sleep (Rechtschaffen and Kales, 1968; Iber and Medicine, 2007; Silber et al., 2007). A typical night of sleep consists of 4–5 sleep cycles of transitions across stages in healthy adults. Each sleep cycle shows progression of sleep stages in the following order, N1 N2 N3 N2 REM (Feinberg and Floyd, 1979). Approximately 75% of sleep is spent in NREM stages. There is also a greater amount of N3 (also termed Slow Wave Sleep, SWS) earlier in the night, whereas REM stage is relatively dominant during later cycles of sleep (Aeschbach and Borbely, 1993).

Each sleep stage is characterized by specific patterns of the brain electrical activity. The N1 stage consists of slow eye movements and low-amplitude low-frequency (4–7 Hz) EEG activity in humans. The sleep stage N2 is dominated by sleep spindle oscillations at 10–17 Hz (in humans) with waxing-and-waning field potentials. Spindles last about 0.5–2 s and recur every 2–20 s. Similar events have been studied in non-human mammals, where their frequency may be as low as 7 Hz. During the N3 delta band (0.2–4 Hz), rhythms are found in the EEG (Rechtschaffen and Kales, 1968; Amzica and Steriade, 1998). This stage is dominated by the slow oscillation, consisting of active (Up) and silent (Down) cortical states alternating at frequency 0.2–1 Hz, that are prominent and visible in the EEG, and in the extracellular and intracellular recordings (Steriade et al., 1991, 1993a; 2001, Werth et al., 1996; Timofeev et al., 2000, 2001). During Up states, most cells within the cerebral cortex are relatively depolarized and may generate action potentials; during Down states, most cortical neurons are hyperpolarized and do not fire (Steriade et al., 1993b; Contreras and Steriade, 1995; Timofeev et al., 2000). The slow oscillation may nest faster spindles, which are commonly found during Down-to-Up state transitions (Molle et al., 2002; Clemens et al., 2007). During REM sleep, electrical brain activity resembles that of the awake state, including mixed signals with bursts of alpha (7–12 Hz) activity (Cantero and Atienza, 2000). Precise coordination of different EEG rhythms during sleep is believed to be critical for memory consolidation which manifests in reactivation of specific neural activity patterns – sleep replay – that have been observed in different brain areas including the hippocampus, amygdala, neocortex and striatum (Nadasdy et al., 1999; Pennartz et al., 2004; Euston et al., 2007; Popa et al., 2010; Bendor and Wilson, 2012).

Neuromodulators, such as acetylcholine (ACh), GABA, histamine (HA), serotonin (5-HT) and norepinephrine (NE), are known to vary significantly during sleep and awake as well as across sleep stages. Compared to the awake state, in N2 and N3, the levels of ACh and monoamines such as HA and NE are reduced while the level of the inhibitory neurotransmitter GABA is increased (Aston-Jones and Bloom, 1981); (McCormick, 1992; Vazquez and Baghdoyan, 2001; Lena et al., 2005; Vanini et al., 2011). During REM sleep, the level of ACh is increased, but monoamines and GABA are reduced (Baghdoyan and Lydic, 2012). While intracellular and synaptic targets of specific neuromodulators are somewhat known, we still lack a clear understanding of how the orchestrated action of many neuromodulators leads to the very specific types of the brain electrical activity in awake and sleep.

In this study, we present a comprehensive computational model of the thalamocortical system implementing effects of neuromodulators and identify the critical intrinsic and synaptic neuronal mechanisms required to explain transitions between sleep stages within an ultradian cycle. We predict that the differences in the temporal dynamics of spindle and delta band oscillations observed in the LFP recordings of cats and mice and ECoG activity of human subjects during SWS (N3 in humans) could arise from the relative differences in the level of neuromodulators. We further apply the model to explain electrical activity observed in propofol-induced anesthesia. Our study predicts that even minor changes of neuromodulators may affect the properties of sleep spindles and sleep slow oscillation in NREM sleep, thus possibly affecting the dynamics of memory consolidation.

Results

Levels of the ACh, HA and GABA vary across sleep stages in vivo. During NREM stages (N2 and N3), the ACh and HA levels are reduced compared to awake; during REM, the ACh and HA levels are increased (Baghdoyan and Lydic, 2012) (Figure 1B). The level of the inhibitory neurotransmitter, GABA, is increased during N2 and N3 compared to awake, while reduced during REM sleep (McCormick, 1992; Vanini et al., 2011). To link the action of these neuromodulators to characteristic brain electrical activity in different sleep stages, we developed a thalamocortical network model implementing the basic effects of the ACh, HA and GABA (see 'Materials and methods'). The model included populations of pyramidal (PY) cells and inhibitory interneurons (IN) in the cortex, and thalamic relay (TC) and reticular nucleus (RE) neurons in the thalamus (Figure 1A); it was designed based on our previous studies of sleep slow oscillations and sleep spindles (Bazhenov et al., 1998, 1999; Timofeev et al., 2000; Bonjean et al., 2011; Chen et al., 2012; Wei et al., 2016).

Model design.

(A) Network structure of the model. PY - cortical principal neurons, IN - cortical inhibitory interneurons, TC - thalamocortical (relay) neurons, RE - thalamic reticular neurons. There were 500 PY, 100 IN, 100 RE and 100 TC cells in the network. (B) The activation function for the h-current under different HA levels, corresponding to the different awake and sleep states as indicated in the legend. ShiftHA = 0 for the REM state. (C) A sequence of changes to the maximal conductances for AMPA, GABA and K+ leak currents in the cortex to model sleep stage transitions. (D) A sequence of changes in K+ leak currents in thalamic cells, and the maximal conductance for GABA synaptic currents within the thalamus.

https://doi.org/10.7554/eLife.18607.003

Specific effects of neuromodulators were implemented by changing the strength of the intrinsic and synaptic currents in the cortical and thalamic neurons (see 'Materials and methods' for more details). While our model does not capture the entire spectrum of changes associated with known effects of neuromodulators, we identified a minimal set sufficient to account for characteristic changes of brain electrical activity across the sleep-wake cycle. Briefly, a reduction of ACh was implemented as an increase in the conductance of the K+ leak current in TC, PY and IN neurons and a reduction of the K+ leak current in RE cells (McCormick, 1992) (Figure 1C,D). The level of ACh also determined the maximal conductance of the excitatory AMPA connections (Figure 1C,D), so that reduction of ACh led to an increase in AMPA connection strength, in agreement with experimental findings (Gil et al., 1997; Kimura et al., 1999; Hsieh et al., 2000). The effect of HA was implemented as a shift in the activation curve of a hyperpolarization-activated cation current, Ih, in TC cells (Figure 1B) (McCormick and Williamson, 1991; McCormick, 1992a). Higher values of HA led to a positive shift in the Ih activation curve (Figure 1B) as seen in vitro (McCormick and Williamson, 1991). The level of HA also determined the strength of the K+ leak conductance in all neurons (see 'Materials and methods' section). The level of GABA determined the maximal conductance of the GABAergic synaptic connections. We also tested in the model the effect of extracellular GABA concentration on tonic inhibition.

Below we will first show that combined action of all these factors is sufficient to explain transitions between the awake state, the N2 and N3 stages and the REM sleep state. Next, we map the relationship between the space of electrical brain activities during sleep and the neurochemical space. We then discuss the impact of specific neuromodulators on brain electrical activity and sleep patterns. Next, we apply our model to explain discrepancies in the animal and human data regarding phase coupling between specific sleep rhythms. Finally, we apply our model to simulate the effect of propofol anesthesia on brain electrical activity.

Collective action of acetylcholine, histamine and GABA explains transitions between sleep and awake states

When the model implemented the intracellular and synaptic changes associated with the putative actions of ACh, HA and GABA, it produced the main features of the different sleep stages observed in human and animal studies. During the awake period, electrical activity in the cortical PY and IN and thalamic RE cells was sparse, that is mean firing was around 2–3 Hz and 4–5 Hz in PY and IN cells, respectively, while the TC neurons primarily remained silent and RE cells were busting intermittently (Figure 2A and Figure 3A). To test the stability of network dynamics, we applied brief thalamic stimulation in the form of a short (100 ms) DC pulse during the awake state. It triggered a transient increase of firing in cortical PY neurons and INs, after which the network returned to the baseline, as observed in experiments with a sensory input during alert awake conditions (Bruno and Sakmann, 2006; Hengen et al., 2016). The average activity (simulated LFP) during awake period (Figure 2C) showed no oscillatory activity or large deflections from baseline, reflecting desynchronized neuronal firing.

Network activity in awake, N2, N3 and REM sleep stages.

(A) Top: Activity of all classes of neurons in the thalamocortical network model (500 PY, 100 IN, 100 RE and 100 TC cells) with the neuromodulators (ACh, GABA and HA) varying as shown in the Figure 1B–D. X-axis is time and Y-axis is neuron index. Colors indicate membrane voltage. For awake, N2, N3 and REM stages, respectively, ACh was 100%, 80%, 50% and 115%; HA was 100%, 40%, 30% and 10%; and GABA was 100%, 115% 130% and 75%. (B) Membrane voltages of the representative selected neurons from the network. (C) LFP calculated from PY network (mean voltage of all 500 PY neurons). Note spindles and slow oscillations during N2 and N3, respectably. (D) Spectrogram (based on the sliding time window FFT) of the LFP shows activity in spindle frequency (8–15 Hz) during N2, rare alpha-burst events (8–15 Hz) during REM, and slow oscillation (0.5–1 Hz) during N3.

https://doi.org/10.7554/eLife.18607.004
Characteristic patterns of electrical activity during different sleep stages.

(A) Top: Network activity in PY, IN, TC and RE cells shown for a 5 s time window during awake. Y-axis is neuron index. Colors indicate membrane voltage. Bottom: membrane voltages of selected neurons from the network (neuron #250 in PY, #50 in IN, RE and TC populations). (BD) Network activity during N2 (B) N3 (C) and REM (D) sleep for a 10 s time window. Legend is same as for panel A. Note a single localized in space 'alpha-burst' event in D.

https://doi.org/10.7554/eLife.18607.005

The N2 stage was characterized in the model by the periods of spindles – thalamically-organized bursts of 7–15 Hz oscillations lasting 0.5 to 2 s each and recurring every 3–20 s (Steriade et al., 1993a). In the model, we observed a transition to N2-like activity following reduction of ACh and HA, and increase of GABA. No external stimulation was applied to initiate or maintain spindles; however, thalamic stimulation could trigger a spindle response. In this state of the model (see N2 in Figure 2), periods of 7–15 Hz oscillatory activity reappeared spontaneously and lasted for at least 2–5 s as revealed by the spectrogram (Figure 2D). Membrane voltages were hyperpolarized in the PY, IN and TC cells due to the reduction of the ACh (Figures 2B and 3B), similar to the intracellular recording (Steriade et al., 1993c). In agreement with the prior intracellular data (Contreras and Steriade, 1996) and computational models (Bazhenov et al., 2000; Bonjean et al., 2011), spindles recurred every 3–10 s.

When the levels of ACh and HA were reduced and GABA was increased compared to the simulations of the N2 state, N3-like activity appeared in the network (Figures 2 and 3C). In this state, a stereotypical slow oscillation consisting of transitions between Up and Down cortical states at a frequency of around 0.5–1 Hz was observed in the LFP, similar to experimental data (Steriade et al., 1993c). Up states were characterized by the higher spiking activity in all cortical pyramidal cells, interneurons and thalamic reticular neurons, while TC cells showed a depolarized state (similar to intracellular recordings from higher order thalamic nuclei [Sheroziya and Timofeev, 2014]) with a few spikes usually at transition from the silent to the active state; the Down states corresponded to the periods of the network silence (Timofeev et al., 2000; Bazhenov et al., 2002; Compte et al., 2003). As in vivo, faster spindle oscillations were found during Up states of the slow oscillation in the model.

To model REM sleep, the level of ACh was increased slightly compared to the awake period, whereas the levels of HA and GABA were reduced. In this network state, there was an overall increase in the cortical activity with very brief periods of 7–15 Hz oscillations that appeared on the background of desynchronized low-frequency cortical firing. The clusters of synchronized firing were somewhat similar to that observed during spindle oscillations. However, they appeared to be more localized in space. In the model, these brief periods of 7–15 Hz synchronized firing (Figure 3D) depended on the reduced level of monoamines in the thalamus and higher excitability in the cortex. We propose that these periods of synchronized oscillations correspond to the alpha or mu rhythms (7–13 Hz activity) that have been reported during REM sleep (Cantero and Atienza, 2000).

Effect of neuromodulators on spindle and slow oscillation power and network synchronization

In order to reveal the specific role of different neuromodulators on network dynamics, the levels of ACh, GABA and HA were varied across a wide range (ACh: 0 to 120%, GABA: 25–225%, HA: 34 to 100% of the awake stage). Power was measured in the spindle (7–15 Hz) and delta (0.5–4 Hz) frequency bands, during 10 s periods for each combination of the neuromodulators. In our model, during N3 the delta band activity was dominated by slow oscillations (0.5–1 Hz). Synchronization among network sites was measured using the phase locking value (PLV) in a broadband frequency (0.5–20 Hz) range, between the LFPs obtained by averaging membrane voltages within five groups of neurons, each comprising 100 PY neurons. Figure 4A shows the summary of the results from all simulations projected to the spindle vs. delta power plane. Several clusters became apparent. Using cluster analysis (Gaussian mixture model) with spindle and delta power and PLV as inputs, 10 different clusters were identified (Akaike information criteria saturated around 10 components). The mean values of the clusters are indicated by the different color spheres in Figure 4A and examples of characteristic activity in each cluster are shown in Figure 4C.

Cluster analysis of the electrophysiological and neuromodulatory spaces.

(A) Delta (0.5–4 Hz) power of the global LFP (mean voltage of the 500 PY neurons) is plotted against spindle (7–15 Hz) power for different combinations of ACh, GABA and HA. Different colors represent different clusters identified based on the cluster analysis (mixed Gaussian model-based clustering was computed with spindle and delta power and Phase Locking Value (PLV) as inputs; 10 different clusters were identified). The PLV was computed for all pairs taken from 5 LFPs (mean voltage time series that were obtained by averaging every 100 PY neurons). Color spheres indicate the mean values of the clusters. (B) Projections of the clusters identified in panel A to the neuromodulator space of the ACh, HA and GABA. Each ellipsoid represents the approximate location of the clusters identified in A. (C) Characteristic electrical activity of the thalamocortical network (top to bottom: PY, IN, TC and RE cells) for each cluster identified in panels A and B. Cluster index is indicated above each plot and corresponds to the index in the panels A and B.

https://doi.org/10.7554/eLife.18607.006

We found that identified cluster states corresponded to either well-defined characteristic types of sleep activity (e.g. sleep spindles or slow oscillations) or various mixed states. When spindle and delta power were both low (e.g. clusters 3 and 4), the activity resembled an awake or REM sleep. In many instances of activity in cluster 4 (Figure 4C), we observed very brief and spatially localized periods of 7–15 Hz oscillations, which were similar to the alpha or mu rhythms recorded during REM sleep. For higher levels of delta power, two groups of clusters were observed: clusters corresponding to low spindle power (e.g. cluster 2) and clusters where spindle power was positively correlated to the delta power (e.g. clusters 5, 6). Cluster 2 with low-spindle power but high-delta power corresponded to the stereotypic transitions between Up and Down states (Figure 4C) as observed during the sleep slow oscillation in vivo, while clusters with a correlated spindle and delta power (clusters 5 and 6) corresponded to the mixed states, where brief (cluster 6) or long (cluster 5) periods of spindles appeared between periods of slow oscillation. Finally, a distinct cluster 1 corresponded to low-delta power, but high-spindle power. This cluster represented quasi-continuous spindle activity (Figure 4C).

The clusters shown in the Figure 4A were identified based on the quantitative characteristics of the electrical activity in the network. In order to link them to specific levels of the neuromodulators, all detected clusters were projected to the neurochemical space. Such projection was performed by fitting an ellipsoid around the means of each cluster in the three-dimensional space of the neuromodulators (Figure 4B). Thus, location of each ellipsoid provides a guide to the range of neuromodulators that produced a specific electrographic pattern. A multivariate ANOVA on the location of the different clusters in the neuromodulator space was significant (p<10e-6), suggesting that the clusters occupy different regions. Clusters 3 and 4 (green ellipsoids in Figure 4B corresponding to the awake or REM sleep like electrical activity) were associated with relatively high ACh levels such as that observed during awake and REM periods. Cluster 2 (violet ellipsoid in Figure 4B corresponding to the slow oscillation) was observed in the low ACh conditions, as found during N3. Clusters 5 and 6 (sleep stage 2 based on electrical activity) corresponded to intermediate levels of the neuromodulators. Finally, cluster 1 was observed for lower HA and higher GABA levels, in agreement with the data showing that strong GABA promotes spindle activity (Parrino and Terzano, 1996). In sum, our analysis confirmed a specific role of neuromodulators in controlling the electrical activity of the thalamocortical network in different sleep and awake conditions, in agreement with past experimental studies (Baghdoyan and Lydic, 2012; McCormick, 1992).

Minimal models of sleep stage transitions

In order to further identify the critical mechanisms responsible for state transitions in the thalamocortical network, we tested the minimal changes of the neuromodulators that were sufficient to generate the essential electrographic characteristics of each sleep stage (Figure 5). Changes in the neuromodulators were restricted to either thalamus or cortex to isolate the role of different regions. Furthermore, we focused on manipulations of the ACh and HA and not GABA, since in our model GABA modulation (within the physiological range) alone was not able to induce any state transitions. We found that in order to induce a transition from awake to N2, a reduction of HA in the thalamus was necessary but produced very focal and scattered spindles (Figure 5C). Reduction of both ACh and HA in the thalamus (Figure 5D) resulted in spindle activity that was similar to that in a full model (Figure 5A, see also Figure 3B). Cortical changes alone did not lead to spindles (Figure 5B) but were required to obtain cortical membrane potentials in the physiological range during N2. We concluded that reduction of ACh and HA in the thalamus, along with reduction of ACh in the cortex, constitute a minimal model needed to simulate the main features of the transition from the awake state to the N2 sleep state.

Minimal models of sleep state transitions.

In all these simulations only specific indicated manipulations were performed. Other model parameters remained unchanged compared to the awake state. The time of the parameter changes is indicated by the bars at the top of each plot. Each plot shows activity of all classes of neurons in the thalamocortical network (500 PY, 100 IN, 100 RE and 100 TC cells) and the LFP. (A) All previously indicated manipulations corresponding to the awake to N2 transition (reduction of ACh and HA, and an increase in GABA) were performed in the thalamic network; no cortical changes. For awake stage: AChTC, AChRE, ShiftHA, LGABA (RE-RE, RE-TC) were 1.0, 1.0, −8.0 and 1.0 respectively. For N2 stage: AChTC, AChRE, ShiftHA, LGABA (RE-RE, RE-TC) were 1.25, 1.25, −3.0 and 1.15 respectively; LACh, AChPY and LGABA(IN-PY) were all fixed at 1.0. (B) All previously indicated manipulations corresponding to the awake to N2 transition (reduction of ACh and increase in GABA) were performed in the cortical network; no thalamic changes. For awake stage: LACh, AChPY, ShiftHA, LGABA(IN-CX) were 1.0, 1.0, -8.0 and 1.0 respectively. For N2 stage: LACh, AChPY, ShiftHAand LGABA(IN-CX) were 1.25, 1.25, -3.0 and 1.15 respectively; AChTC, AChREand LGABA (RE-RE, RE-TC) were fixed at 1.0. (C) HA level was reduced in the thalamic network; no any other changes were performed. For awake and N2: only  ShiftHA was changed from −8 to −3; all other variables were fixed at awake level. (D) HA and ACh levels were reduced in the thalamic network; no any other changes were performed. For awake stage: AChTC, AChRE, ShiftHA , LGABA (RE-RE, RE-TC) were 1.0, 1,0, −8.0 and 1.0 respectively. For N2 stage: AChTC, AChRE, ShiftHA , LGABA (RE-RE, RE-TC) were 1.25, 1.25, −3.0 and 1.15 respectively; LACh, AChPY, LGABA(IN-CX) were fixed at 1.0. (E) All specific manipulations corresponding to the N2 to SWS transition (ACh, HA and GABA levels modified) were performed in the thalamic network; no changes in the cortex. For N2 stage: AChTC, AChRE, ShiftHA, LGABA (RE-RE, RE-TC) were 1.25, 1.25, − 3.0 and 1.15 respectively. For N3 stage: AChTC, AChRE, ShiftHA, LGABA (RE-RE, RE-TC) were 2.0, 2.0, −8.0 and 1.3, respectively; LACh, AChPY and LGABA(IN-PY) were fixed at N2 stage level. (F) ACh level was reduced in the cortex; no any other changes were performed. For N2 and N3 stages, LACh was 1.0 and 1.25; AChPY was 1.0 and 1.25; all other variables were fixed at N2 level. (G) Only AMPA strength was increased in the cortex. For N2 and N3 stages, LACh was 1 and 1.25, correspondingly; all other variables were fixed at N2 level. (H) K+ leak current was decreased in the PY neurons. For N2 and N3 stages, AChPY was 1 and 1.25, correspodingly; all other variables were fixed at N2 level.

https://doi.org/10.7554/eLife.18607.007

Transition from N2 to N3 could be observed when ACh was reduced in the cortex alone (Figure 5F). In contrast, even when all the other proposed neuromodulatory changes (reduction of ACh and HA, and increase of GABA) were applied in the thalamus, the thalamocortical network did not display N3-like activity (Figure 5E), suggesting that appearance of the slow oscillation during N3 requires neuromodulatory action within the cortex itself. Since reduction of the ACh has two main actions in our model, an increase in excitatory connection strength and an increase of leak currents, we examined if either of these changes alone could result in transition to the N3. Increasing excitatory AMPA conductance alone led to an activated state (Figure 5G), and an increase in the K+ leak conductance alone led to a hyperpolarized quiescent state in the network (Figure 5H). Thus, both actions together were required for transition between the N2 and N3.

Comparison between sleep stages in cats, mice and humans

We next examined if our thalamocortical model implementing varying levels of neuromodulators can explain LFP data from regions of cat and mouse neocortex across sleep stages (Figure 6 A,B). Surprisingly, visual scoring of cat LFPs revealed that there were only relatively short periods of the S2 (equivalent to N2; to describe animal experiments we use here terminology common in the animal studies) compared to the SWS (equivalent to N3). This result was consistent across animals. The power in the spindle (8–15 Hz) and delta (0.2–4 Hz) frequency bands was measured across all channels and the average power was plotted in Figure 6C and D, separately for the cortical and thalamic (VPL) recordings. The amplitude of the slow oscillation and spindles were both high during the S2 and SWS, and were reduced in awake and REM. This suggests that the SWS in naturally sleeping cats is marked with the relatively high-spindle activity. Slow (0.2–1 Hz) oscillations provided the major contribution to delta band activity during SWS. To further examine the interaction between sleep spindles and slow oscillations, we plotted the power in the spindle vs delta frequency bands across different cortical LFP electrodes (Figure 6E). We found a strong correlation between slow oscillation and spindle activity. In mice (Figure 6F,G), a similar relationship between spindle and delta power was observed across recordings from five animals (3 channels per mouse); a positive correlation was observed in majority of the channels (9 channels out of 15 channels showed significant correlation of which 7 were positively correlated [r ranged from 0.17 to 0.43; p<0.05] and 2 were negatively correlated [r was −0.41 and –0.43; p<0.05]). Recordings from four out of five mice revealed significant correlations. In one mouse (Figure 6F–G), frontal and posterior somatosensory electrodes showed a significant positive correlation between slow oscillations and spindle activity.

LFP recordings from cats and mice show positively correlated spindle and delta power during NREM sleep.

(A). In vivo LFP recordings from different locations in the cortex (black lines) and the thalamus (red lines) of a non-anesthetized cat. (B) Locations of the recording electrodes. (C,D) Power at delta frequency (0.2–4 Hz) (C) and spindle frequency (8–15 Hz) (D) measured from 5 s time windows (1 s sliding window). (E) Spindle power is plotted against delta power for different channels. PCA was applied prior to computing the power and the principal components correlating strongly to EOG were removed; channel data were recomputed from remaining components to remove eye movement artifacts. Each dot represents the power in spindle and delta bands measured from a 30-s period of cat recordings. Note a positive correlation between delta and spindle power when combining S2 and SWS epochs. (F) Top: Hypnogram of recordings in mice. Bottom: In vivo LFP recordings from different cortical locations in mice. Note lack of S2 periods in mice. (G) Spindle power is plotted against delta power for different channels in mice. Each dot represents the power in spindle and delta bands measured from a 30-s window of data. There is a positive correlation between delta and spindle power during SWS epochs.

https://doi.org/10.7554/eLife.18607.008

Our analysis of the neuromodulator space (Figure 4) suggests that a moderate reduction of ACh would lead to a network state with correlated variations in spindle and delta power (clusters 4, 5 and 6 in Figure 4A). Thus, in order to replicate our animal data (cat and mouse LFP activity), a computational model with a relatively high level of ACh during SWS (reduction to 75% of the baseline awake level) was required (Figure 7A). In contrast, in our control simulations of SWS (N3), the ACh was reduced to 57% of the baseline (Figures 2 and 3). Using the new model with the higher level of the ACh during SWS, we found a relative increase in the spindle activity during the SWS (Figure 7A) and a strong positive correlation between slow oscillations and spindle activity (Figure 7B, right), as was observed in the animal LFP data (Figure 7B, left).

Model predicts temporal characteristics of the LFP activity in cats and mice.

(A) The network model implementing moderate ACh reduction (to 75% of the baseline awake level) reproduces sleep stage transitions with pattern of activity similar to that recorded from cats and mice. (B) Spindle power (7–15 Hz) averaged across all electrodes is plotted against delta band power (0.5–4 Hz) for different sleep stages. PCA was applied to remove eye movement artifacts in cat recordings. Right: Spindle power vs. delta power based on the model simulations. Average activity within groups of 100 neurons was used as an estimate of the LFP. Each dot represents the power in spindle and delta bands measured from 2 s windows of data obtained from the model simulations.

https://doi.org/10.7554/eLife.18607.009

Next, ECoG recordings from humans during different sleep stages were examined for the relationship between delta (0.01–2 Hz) and spindle (9–17 Hz) power (Figure 8A,B). (Note, these slightly different frequency bands are conventional in human studies.) In contrast to the pattern found in cats and mice, there was a significant negative correlation between the spindle and the delta power observed in successive 30 s epochs obtained during the N2 and N3 (Figure 8C and D) in humans. Generally, the slow oscillation was higher during N3 than during N2, whereas spindle activity was higher in N2 than in N3, resulting in a negative correlation across epochs. The average Pearson’s correlation coefficient across 42 channels was −0.2595 (sd = 0.2281) (Figure 8E). Fifteen out of 42 electrodes showed a significant negative correlation (p<0.001), and no electrodes showed a significant positive correlation. We also observed periods of increased power in the 7–15 Hz frequency band during REM and awake suggesting transient periods of alpha (or mu) rhythm. To simulate the negative correlation between sleep slow oscillation and spindle activities, we again applied results from our analysis of the neuromodulator space (Figure 4). When the ACh level was reduced to 45% of the baseline during the N3 (compared to 57% used in the control model, see Figures 2, 3), the network activity revealed a negative correlation between the delta and spindle frequency bands (Figure 8 F and G). These findings suggest that a relatively low level of the ACh during the N3 could possibly explain the negative correlation between the slow oscillation and spindle activities as observed in human ECoG recordings.

ECoG recordings from a human patient show negatively correlated spindle and delta power during NREM sleep.

(A) Patient MG29 had a 64-contact grid and 4 strips implanted over the left hemisphere. After rejection of electrodes that either had poor contact or showed significant epileptic interictal activity, 42 electrode contacts remained for analysis. (B) Recording of patient for 9 hr beginning in the evening, during which patient was both awake and asleep. Shown are LFP over four grid electrodes, along with the right EOG channel. Beneath that is the sleep staging done in 30 s epochs. Staging was scored using information from scalp, EOG and intracranial electrodes. Data from Awake, REM, N2 and N3 were used for further analysis. A significant decrease in overall LFP amplitude was seen in REM and occasionally during awake. (C) Delta and spindle power over the same period. A Fast Fourier Transform with a window size of 30 s was done for delta (band 0.01–2 Hz) and spindle (9–17 Hz) band frequencies to obtain their power. All 42 electrodes were averaged and are plotted here. (D) Delta vs. spindle bands for 30 s time epochs were plotted for individual electrodes and separated based on sleep staging. There is a negative correlation between delta (dominated by slow oscillation) and spindle activities for these electrodes when combining N2 and N3 epochs. (E) A histogram of Pearson's R correlation coefficients between delta and spindle power during N2 and N3 across all 42 electrodes. (F) The computational model implementing significant ACh reduction (to 45% of the baseline awake level) reproduces sleep stage transitions with pattern of activity similar to that recorded from humans. (G) Spectral analysis of model activity. Power at delta frequency (0.2–4 Hz, red) and spindle frequency (8–15 Hz, blue) were measured from FFT obtained with non-overlapping sliding 2 s windows, similar to the analysis shown in C of actual recordings. (H) Inverse correlation of delta and spindle power in the model during combined N2 and N3 activity, similar to human recordings in D.

https://doi.org/10.7554/eLife.18607.010

Modeling the effects of propofol anesthesia requires the combined action of ACh, HA and GABA

We next applied our model to explain the features of the electrical activity induced by propofol anesthesia. In humans, propofol has shown to increase 8–15 Hz oscillations frontally and to reduce alpha oscillations posteriorly (Murphy et al., 2011; Purdon et al., 2013). Propofol application also increased the slow oscillations in the delta (1–4 Hz) frequency range across all regions. The onset of the loss of consciousness was closely matched by the increase of slow (<1 Hz) and 8–15 Hz oscillations (Purdon et al., 2013). We examined if our model could capture these changes in the EEG due to the action of propofol.

Propofol is a GABA agonist and is known to increase the decay time constant of the GABA-A mediated inhibitory post-synaptic potential (Kitamura et al., 2003). Multiple lines of evidence suggest that propofol also acts by reducing the action of ACh and HA. Indeed, propofol is shown to decrease the level of ACh in frontal cortex (Kikuchi et al., 1998) and to attenuate ACh receptor responses (Flood et al., 1997; Murasaki et al., 2003). Increasing ACh transmission prevents the action of propofol in humans (Meuret et al., 2000). Propofol inhibits the activity of the tuberomammillary nucleus (Nelson et al., 2002), although its action on ventrolateral preoptic nucleus neurons (Liu et al., 2013) leading to reduction of HA.

When these effects of ACh and HA were implemented in our model, as described before, and the decay time constant of GABA was increased, we observed a large increase in slow oscillations as well as oscillatory activity in 8–15 Hz frequency band (Figure 9A). In comparison to the model of natural sleep (Figure 9C), there was elevated 8–15 Hz power in the propofol condition (Figure 9D). This is consistent with our recordings from cats (not shown), where we found that 8–15 Hz power under propofol was high compared to natural sleep or ketamine-xylazine anesthesia. These results are also consistent with the previous computational models of propofol application, where the increase in decay time constant of IPSPs led to increase in 8–15 Hz oscillations (Ching et al., 2010). However, in our model, when only the time constant of the GABA-A IPSP was increased compare to the awake condition, the model failed to generate slow oscillations or the 8–15 Hz synchronized activity (Figure 9B). Overall, our study suggests that the known action of propofol in vivo may require its indirect effect on both the ACh and HA neuromodulatory systems.

Mechanism of action of propofol induced anesthesia.

(A) Model of propofol anesthesia implementing inhibition of ACh and HA release and increase in the decay time constant of the GABAergic IPSPs. From top to bottom: Spatio-temporal pattern of activity in PY neurons, average activity in PY network (simulated LFP), band-filtered (7–15 Hz) PY LFP. The GABA decay time constant was increased by 150%; ACh and HA were reduced to the same level as in simulations of the natural sleep in cats. Note significant amount of spindle activity. (B) The model where only the GABA time constant was increased by 150%. Note an almost complete absence of spindle activity. (C) Network activity during simulated natural SWS in cat model. Note decreased spindle activity compared to the propofol simulations in panel A. (D) Power in 0.2–4 Hz and 8–15 Hz bands for all three conditions. The network was divided into 10 groups of 50 neurons. Membrane voltages were averaged within each group, then FFT was used to estimate power in each group. Bar height indicates average across 10 groups, error bar indicates standard deviation across 10 groups.

https://doi.org/10.7554/eLife.18607.011

Role of GABA in sleep stage transitions

Experimental studies suggest that the concentration of extracellular GABA in the neocortex is higher during NREM sleep and it is reduced during REM sleep compared to the awake state (Vanini et al., 2012). Increase of extracellular GABA may reflect an increase of synaptic inhibition. Indeed, during slow wave sleep, synaptic excitation and inhibition are balanced (Rudolph et al., 2007Haider et al., 2006). Since decrease of ACh during stage N2/N3 sleep is known to increase excitatory connections (Gil et al., 1997), synaptic inhibition and thus phasic GABA release may be also increased, which would then be reflected in elevated extracellular GABA. The effect of increase of synaptic inhibition was implemented in our baseline model as reported above.

Nevertheless, the origin of the change in the extracellular GABA is still not fully understood. Increase of extracellular GABA may increase tonic inhibition but not necessarily be associated with increase of the phasic GABA release. Further, several studies revealed co-release of GABA with ACh and HA (Saunders et al., 2015; Yu et al., 2015) which suggests a decrease of GABA during NREM sleep while the observations from microdialysis experiments (Vanini et al., 2012) suggest an increase of GABA during NREM sleep. Therefore, in our study, we also tested models with (a) no change in the level of the GABA release; (b) increase of tonic inhibition reflecting increase of the extracellular GABA, based on the observations from microdialysis experiments (Vanini et al., 2012); and (c) tonic inhibition proportional to the ACh and HA levels, based on the evidence of the co-release of GABA with ACh (Saunders et al., 2015) and HA (Yu et al., 2015). In all these conditions (Figure 1), we found transitions between sleep stages similar to those reported in the baseline model. However, the model implementing co-release of GABA with ACh had poorly formed spindle activity during N2 and elevated alpha activity during REM sleep. Future experiments are required to distinguish between changes in tonic and phasic inhibition to match specific features of the synchronized oscillatory activity across sleep stages.

Effect of GABA on the network state transitions.

(A) Both phasic and tonic GABA conductances were fixed for entire simulation. Note the state transitions similar to the baseline model (Figure 2). (B) Tonic GABA conductance (conductance of the miniature IPSPs) was varied based on the measured GABA levels in microdialysis experiments (mIPSPs were modified as 115%, 130% and 75% of awake stage for N2, N3 and REM stages, respectively), while phasic GABA conductance was remained fixed. (C) Assuming co-release of ACh with GABA, tonic GABA conductance was proportionally varied with ACh concentration (mIPSP for GABA conductance was scaled as 85%, 70% and 125% of awake stage for N2, N3 and REM stages respectively). Transition between sleep stages was observed in all these conditions, but spindles were less synchronous in C; in A and B alpha-bursts were less common.

https://doi.org/10.7554/eLife.18607.012

Discussion

The functional state and the patterns of electrical activity of human and animal brain are continuously influenced by a broad range of neuromodulators (McCormick, 1992; Steriade et al., 2001; Baghdoyan and Lydic, 2012). Changes in the level of neuromodulators have been correlated with sleep induction as well as transitions between characteristic EEG sleep patterns (McCormick, 1992; McCormick and Bal, 1997; Steriade et al., 2001). However, the mechanisms whereby specific cellular and synaptic changes triggered by the combined action of neuromodulators result in transitions between sleep stages and their precise electrical characteristics remain poorly understood. In this new study, we used biophysically realistic computational models of the thalamocortical system to identify the minimal set of cellular and network level changes, linked to the putative action of the neuromodulators, that is sufficient to explain the characteristic neuronal dynamics during sleep as well as the transitions between primary sleep stages and from sleep to awake state. Our study predicts a critical role of the neuromodulators in controlling the precise spatio-temporal characteristics of neuronal synchronization that manifest in major sleep rhythms – sleep spindles and slow oscillations. Our results are consistent with in vivo intracellular and LFP recordings from animals and ECoG data from human subjects, as reported in this study.

Previous studies identified intrinsic and synaptic mechanisms involved in generating specific types of sleep rhythms. The survival of the slow oscillation (<1 Hz global electrical activity found during deep (N3) sleep) after extensive thalamic lesions in vivo (Steriade et al., 1993b), its existence in cortical slice preparations (Sanchez-Vives and McCormick, 2000), and the absence of the slow oscillations in the thalamus of decorticated cats (Timofeev and Steriade, 1996) point to the primarily intracortical origin for this rhythm. (Note that thalamus can be actively involved in synchronization of the SO [Lemieux et al., 2014]). In vivo and in vitro studies suggest that the minimal substrate accounting for spindle oscillations (7–15 Hz recurrent oscillatory activity found during stage 2 (N2) of NREM sleep) is the interaction between thalamic reticular and relay cells (Steriade and Deschénes, 1984; Steriade et al., 1985; Steriade and Llinas, 1988, 1990; von Krosigk et al., 1993). Previous computational studies proposed a minimal set of the mechanisms sufficient to model sleep spindles and slow oscillations (Destexhe et al., 1996; Bazhenov et al., 1998, 1999; Timofeev et al., 2000; Bazhenov et al., 2002; Compte et al., 2003; Hill and Tononi, 2005; Bonjean et al., 2011). An increase in the K+ leak current was identified as a critical component for the transition between awake and slow wave sleep (Bazhenov et al., 2002; Hill and Tononi, 2005), and it was predicted that synchronization of the slow oscillation depends on cortico-cortical connections (Hill and Tononi, 2005). Several phenomenological and reduced mathematical models have been proposed to test the effects of neuromodulators on the sleep-wake cycle (McCarley and Hobson, 1975; Robinson et al., 2011; Schellenberger Costa et al., 2016). Nevertheless, none of these previous models examined transitions among all major sleep stages based on biophysical mechanisms.

In this study, we found that reduction in HA and ACh in the thalamus was required for induction of N2-like activity, characterized by recurrent spindles with characteristics consistent with animal and human data. Other neuromodulatory changes (such as change of GABA conductance) were not sufficient alone to induce N2-like activity, but played a role in determining the power and the synchrony of oscillations. Reduction in the ACh levels within the cortical network was identified to be a minimal requirement for the transition from the N2 (spindles) to the N3 (slow oscillation). Our study of the minimal models revealed the importance of the cell-type-specific action of the neuromodulators. Indeed, in order to model sleep stage transitions, the crucial changes were in HA, primarily acting on thalamic neurons, and in ACh, primarily acting within a cortex. All the network rhythms reported in this study resulted from intrinsic (autonomous) dynamics of the thalamocortical network; no external stimulation was applied to entrain any of the oscillatory activities in the model.

By parametrically varying the levels of the ACh, HA and GABA, we observed a wider range of activities in the thalamocortical network than have been previously reported in the computational models, but which do correspond to those that have been observed experimentally. While we found network states dominated by either spindles or slow oscillations, under certain conditions (such as only moderate reduction of the ACh level) we also observed mixed states combining these two major rhythms. Such mixed states may occur in vivo under normal conditions (Aeschbach and Borbely, 1993; Muller et al., 2006); our recent behavior study suggests that the strengths of the phase amplitude coupling between the spindle and the slow oscillation during NREM sleep correlates with memory consolidation (Niknazar et al., 2015). We speculate that the mixed states may become prevalent in some pathological conditions such as in Alzheimer’s disease where sleep is altered and changes to the neuromodulatory system are reported (Wulff et al., 2010). Furthermore, our results on the putative relationship between levels of neuromodulators and characteristic sleep electrical activities could be potentially relevant to understanding the changes observed with aging. Indeed, the level of the HA in cerebrospinal fluid is known to be elevated with aging (Prell et al., 1988), where spindle power is known to be reduced (Crowley et al., 2002). In our model, the HA level is one of the major predictors of spindle power. Thus, we can speculate that an increase of HA with age could contribute to the decline of the spindle power.

In our study, the effects of neuromodulators were limited by changing the following four model properties – maximal conductances of the K+-leak, Ih, GABA and AMPA currents. We found that specific combinations of changes of these core parameters were sufficient to cause transitions between sleep stages. We also tested a more detailed model, implementing the effects of ACh and HA on the other K+ currents as well as on the persistent Na+ and potassium M-currents (Constanti and Galvan, 1983; McCormick and Williamson, 1989; Mittmann and Alzheimer, 1998). We observed that this extended model qualitatively captured results similar to the minimal model. In contrast, after fixing the K+-leak, Ih, GABA and AMPA currents but leaving the action of neuromodulators on the other channels, the model failed to show transitions between sleep stages, suggesting that the cellular and synaptic properties identified in our study are critical in determining patterns of electrical activity in the thalamocortical network.

In developing our model, we assumed that the levels of ACh, GABA and HA were different between N2 (Stage 2) and N3 (Stage 3 or SWS). These assumptions are consistent with data on spiking activity in the Basal Forebrain region (Aston-Jones and Bloom, 1981; Szymusiak et al., 2000), suggesting a difference in ACh between the N2 and N3. Another assumption used in this study was that ACh modulates the strength of the excitatory AMPA synapses in the intracortical and thalamocortical circuits. Indeed, ACh suppresses the spread of excitation in vivo (Kimura et al., 1999) and reduces evoked responses (Chauvette et al., 2012). In slice experiments, ACh activation of muscarinic receptors suppressed both intracortical and thalamocortical excitatory connections, while nicotinic activation led to suppression of the intracortical connections and enhancement of the thalamocortical connections (Gil et al., 1997; Hsieh et al., 2000).

By comparing recordings from animals (cats and mice LFP) and human ECoG data across sleep stages, we found a significant difference in the interaction between major sleep rhythms (spindles and slow oscillation) during NREM sleep. In cats and mice, there was a strong positive correlation between spindle and delta power, while in humans there was a negative correlation between oscillations in these frequency bands. This difference may reflect the prominence of N2 (stage 2) in humans, in whom N2 has the longest duration of any sleep stage. In contrast, in cats, N2 is very rare compared to N3. In humans, N2 is mainly characterized by spindles with relatively short periods of slow oscillation or isolated K-complexes. (K-Complexes are large downward deflections in EEG often followed by spindle activity and is one of the electrophysiological markers of N2 (Loomis et al., 1937; Amzica and Steriade, 1997; Cash et al., 2009). In N3, spindles continue to occur in conjunction with the slow oscillation, but the EEG is dominated by slow activity. Our model predicts that this inter-species difference can be explained by a difference in the relative level of ACh during NREM sleep. This result is consistent with microdialysis measurements in cats where ACh was present during SWS albeit at low concentrations compared to wakefulness (Marrosu et al., 1995).

Surprisingly, we found that even minor changes of neuromodulators may lead to significant changes in the characteristics of brain electrical activity during sleep. Our model predicts that fluctuations of the level of ACh over time may explain the appearance of short periods of slow oscillations in N2 in humans (Achermann and Borbely, 1997) or modulate gamma activity during Up-states of the slow oscillation (Mena-Segovia et al., 2008). We hypothesize that transient activities such as K-Complexes during N2 sleep or PGO waves during REM stage may arise from fluctuations in ACh and/or other neuromodulators. Ultimately, our study predicts that transient changes in the level of neuromodulators may affect the interaction and characteristic properties of sleep rhythms and may thus affect functional outcomes of sleep, such as the impact of sleep on memory and learning.

Materials and methods

Model description

Neuromodulators and sleep stages

Our model implements a set of intrinsic and synaptic mechanisms that are known to be related to effects of ACh, HA and GABA. Specifically, ACh modulated the potassium leak currents in all neurons and the strength of excitatory AMPA connections in the cortex; HA modulated the strength of the hyperpolarization-activated cation current, Ih, in thalamic relay neurons; and GABA modulated the maximal conductance of the GABAergic synapses in cortex and thalamus. Based on experimental data, when compared to the awake stage, the levels of ACh and HA were reduced during NREM sleep (N2 and N3) and increased during REM sleep (Vanini et al., 2012). Concentration of the extracellular GABA in the neocortex is higher during NREM sleep and reduced during REM sleep compared to the awake state (Vanini et al., 2012). The origin of this change in the level of GABA is not fully understood. We examined the effect of changes in both phasic and tonic inhibitions. Thus, we tested three different models: (1) no change in the level of the GABA release; (2) miniature IPSPs were changed based on the measurements of extracellular GABA levels in microdialysis experiments (Vanini et al., 2012); and (3) miniature IPSPs were modified based on the assumption that GABA is co-released with ACh and HA, as shown in some experiments (Saunders et al., 2015; Yu et al., 2015).

Intrinsic currents

The cortical neuron model included dendritic and axo-somatic compartments, similar to the previous studies (Timofeev et al., 2000; Bazhenov et al., 2002; Chen et al., 2012) and was described by the following equations,

CmdVddt=[AChPY IdKleak]IdleakIdNaIdNapIdCaIdKCaIdKmIsyngcs(VdVs)=IsNaIsKIsNap

where the subscripts s and d correspond to axo-somatic and dendritic compartments, AChPY is the variable which is inversely proportional to the level of ACh in the cortex (AchPY=1 for awake state, AChPY= 1.25 for N2, AChPY=1.8 for N3 and AChPY= 0.85 during REM) and determines the strength of the potassium leak current (IdKleak), Idleak is the Cl leak currents, INa is fast Na+ channels, INap  is persistent sodium current, IK  is fast delayed rectifier K+ current, IKm is slow voltage-dependent non-inactivating K+ current, IKCa is slow Ca2+ dependent K+ current, ICa is high-threshold Ca2+ current, Ih is hyperpolarization-activated depolarizing current and Isyn is the sum of the synaptic currents to the neuron. All intrinsic currents were of the form (VE), where g is the conductance, V is the voltage of the corresponding compartment and E is the reversal potential. The detailed descriptions of individual currents are provided in previous publications (Bazhenov et al., 2002; Chen et al., 2012). In the extended model (Figure 11), the conductances of IdKCa, IdKm  and IdNap were also scaled by AChPY variable. The conductances of the leak currents were 0.007 mS/cm2 for IdKleak and 0.023 mS/cm2 for Idleak. The maximal conductances for the voltage and ion-gated intrinsic currents were, IdNap: 2.0 mS/cm2 IdNa: 0.8 mS/cm2, IdCa: 0.012 mS/cm2, IdKCa: 0.015 mS/cm2, IdKm: 0.012 mS/cm2, IsNa: 3000 mS/cm2, IsK: 200 mS/cm2 and IsNap: 15 mS/cm2. Cm was 0.075 μF/cm2.

Effect of additional parameter changes on the stage transitions.

(A) Network activity in an extended model with ACh and HA having effects on the all K+ channels (Ca2+-dependent K+ current, M-current, K+ leak current), persistent Na+ current and H-current. (B) Network activity in the model implementing only additional effects of ACh and HA (Ca2+-dependent K+ current, M-current, persistent Na+ current) but having core intrinsic and synaptic properties (K+-leak, Ih, GABA and AMPA) fixed as in the awake model. This model failed to simulate transitions between sleep stages. (C and D) Zoom in of the periods during spindle (C) and slow oscillation (D) from panel A.

https://doi.org/10.7554/eLife.18607.013

The following equation describes the cortical inhibitory (IN) neuron:

CmdVddt=[AChPY IdKleak]IdleakIdNaIdCaIdKCaIdKmIsyngcs(VdVs)=IsNaIsK

The conductances for the leak currents for IN neurons were 0.034 mS/cm2 for IdKleak and 0.006 for Idleak. Maximal conductances for the active currents were, IdNa: 0.8 mS/cm2, IdCa: 0.012 mS/cm2, IdKCa: 0.015 mS/cm2, IdKm: 0.012 mS/cm2, IsNa: 2500 mS/cm2 and IsK: 200 mS/cm2.

The thalamic relay (TC) neuron consisted of a single compartment,

dVdt=[AChTCIKleak]IleakINaIKILCaIhIsyn

where AChTC was 1, 1.25, and  0.85 during awake, N2, N3 and REM sleep correspondingly. The conductances of the leak currents were, Ileak: 0.01 mS/cm2, IKleak: 0.007 mS/cm2. The maximal conductance for the other currents were, fast Na+ (INa) current: 90 mS/cm2, fast K+ (IK) current: 10 mS/cm2, low threshold Ca2+( ILCa) current: 2.5 mS/cm2, hyperpolarization-activated cation current (Ih): 0.015 mS/cm2. The effect of HA on Ih was implemented as a shift in the activation curve based on experimental evidences as follows,

Ih=g([O]+k[OL])(VEh)Cα(V)OCβ(V)OP0+4Ca2+k1P1P0+4Ca2+k2P1O+P1k3OLO+P1k4OLα(V)=h(V)/τs(V)β(V)=(1h(V))/τs(V)h(V)=11+e(V+75+ShiftHA)/5.5

where ShiftHA was −8 mV, −3 mV, −2 mV and 0 mV for awake, stage 2, SWS and REM sleep correspondingly. The values of other constants were the same as in previous study (Bonjean et al., 2011).

The RE neuron was modeled as a single compartment as follows,

dVdt=[AChREIKleak]IleakINaIKICaIsyn

where AChRE was 1, 0.8, 0.5 and 1.15 for awake, N2, N3 and REM stage. Since it is known that ACh reduces excitability of the RE cells (McCormick, 1992), the ACh was reduced during NREM stage compared to awake. The conductances for leak current were, Ileak: 0.05 mS/cm2, IKleak: 0.016 mS/cm2. The maximal conductance for the other currents were, fast Na+ (INa) current: 100 mS/cm2, fast K+ (IK) current: 10 mS/cm2, low threshold Ca2+(ILCa) current: 2.2 mS/cm2.

Network structure

Network architecture of the thalamocortical model was based on the previous studies of the spindle and slow-wave activity (Bazhenov et al., 2002; Bonjean et al., 2011; Chen et al., 2012; Wei et al., 2016). The thalamocortical network consisted of five hundred cortical pyramidal (PY) neurons, one hundred cortical inhibitory (IN) neurons, one hundred thalamic relay (TC) neurons and one hundred thalamic reticular (RE) neurons. The network connectivity between different groups of neuron is shown in Figure 1A. One-dimensional network topology was used in each layer and connectivity was local and determined by the radius of connections. The PY neurons received AMPA and NMDA synapses from other PY neurons within radius of 5. The TC neurons projected to RE neurons through AMPA synapses (radius of 8), and connections from RE to TC neurons included GABA-A and GABA-B synapses (radius of 8). Thalamocortical connections were wider and mediated by AMPA synapses (TC to PY: radius of 10; TC to IN: radius of 2); corticothalamic connections included AMPA synapses (PY to TC: radius of 10; PY to RE: radius of 8) and were stronger to reticular neurons than to thalamic relay neurons.

Synaptic currents

GABA-A, NMDA and AMPA synaptic currents were described by first-order activation schemes (Destexhe et al., 1994). The equations for all synaptic currents used in this model are given in our previous publication (Bazhenov et al., 2002; Chen et al., 2012). Briefly, below we mention only the relevant equations that were used to model effect of the varying levels of ACh and GABA.

IsynAMPA=LAchgsyn[O](VEAMPA)
IsynGABA=LGABAgsyn[O](VEGABA)

LACh was 1, 1.25, 2.0 and 0.8 for awake, N2, N3 and REM stages respectively for cortical PY-PY, TC-PY and TC-IN connections; LGABA was 1, 1.15, 1.3 and 0.75 for awake, N2, N3 and REM stages respectively for cortical IN-PY, RE-RE and RE-TC connections. In these equations, gsyn is the maximal synaptic conductance and [O(t)] is the fraction of open channels, Eampa(= 0 mV) and Egaba(= −90 mV) are reversal potential for AMPA and GABA. The maximal conductances were (subscript indicates cell types and superscript indicates type of connection): gPYPYAMPA = 0.024 μS, gPYPYNMDA = 0.001 μS, gPYTCAMPA = 0.005 μS, gPYREAMPA = 0.015 μS,, gPYINAMPA = 0.012 μS, gPYINNMDA = 0.001 μS, gINPYGABAA = 0.024 μS, gTCCXGABAA = 0.02 μS, gTCINGABAA = 0.02 μS, gREREGABAA = 0.05 μS, gRETCGABAA = 0.05 μS, gRETCGABAB = 0.02 μS and gTCREAMPA = 0.05 μS. Detailed description of the O(t) is given in a previous publication (Bazhenov et al., 1998; Timofeev et al., 2000) and is based on first-order activation schemes. In addition, spontaneous miniature EPSPs and IPSPs were implemented for PY-PY, PY-IN and IN-PY connections. The synaptic dynamics followed the same equations as the regular PSPs, and their arrival times were modeled by Poisson processes with time-dependent mean rate, with next release time given by function: (2/(1 + exp(-(t-t0)/20))-1)/250, where t0 is a time instant of the last presynaptic spike. The maximal conductance for miniature PSP were: gmini(PYPY)AMPA = 0.033 μS, gmini(PYIN)AMPA = 0.02 μS and gmini(INPY)AMPA = 0.02 μS. Short-term depression of intracortical excitatory connections was also included, similar to previous publications (Bazhenov et al., 2002; Chen et al., 2012).

Stimulation protocol

A brief DC input lasting for 100 ms and amplitude of 1.5 nA was applied to TC neurons during awake state in some simulations to test stability of the network dynamics. The DC pulse was modeled using two Heaviside functions. No periodic stimulation was applied in any of the modeling experiments.

In vivo experiments in cats

Experiments were conducted on adult non-anesthetized cats. The cats were purchased from an established animal breeding supplier. Good health conditions of all animals were certified by the supplier and determined upon arrival to the animal house by physical examination, which was performed by animal facilities technicians and a veterinarian in accordance with requirements of the Canadian Council on Animal Care. The protocol for experiments involving cats was approved by CPAUL, Comité de protection des animaux de l'Université Laval (Protocol # 2012–174-4). The surgery was performed on animals 5–20 days after their arrival to the local animal house. We recorded field potentials from several cortical areas and from the VPL thalamic nucleus of cats during natural sleep/wake transitions.

Preparation for cat surgeries

Chronic experiments were conducted using an approach similar to that previously described (Steriade et al., 2001; Timofeev et al., 2001). For implantation of recording chamber and electrodes, cats were anesthetized with isoflurane (0.75–2%). Prior to surgery, the animal was given a dose of preanesthetic, which was composed of ketamine (15 mg/kg), buprenorphine (0.01 mg/kg) and acepromazine (0.3 mg/kg). After site shaving and cat intubation for gaseous anesthesia, the site of incision was washed with at least three alternating passages of a 4% chlorexidine solution and 70% alcohol. Lidocaine (0.5%) and/or marcaine (0.5%) was injected at the site of incision and at all pressure points. During surgery, electrodes for LFP recordings, EMG from neck muscle, and EOG were implanted and fixed with acrylic dental cement. Eight to ten screws were fixed to the cranium. To allow future head-restrained recordings without any pressure point, we covered four bolts in the dental cement that also covered bone-fixed screws and permanently implanted electrodes. Throughout the surgery, the body temperature was maintained at 37°C using a water-circulating thermoregulated blanket. Heart beat and oxygen saturation were continuously monitored using a pulse oximeter (Rad-8, MatVet) and the level of anesthesia was adjusted to maintain a heart beat at 110–120 per minute. A lactate ringer solution (10 ml/kg/hr, intravenously [i.v.]) was given during the surgery. After the surgery, cats were given buprenorphine (0.01 mg/ kg) or anafen (2 mg/kg) twice a day for 3 days and baytril (5 mg/kg) once a day for 7 days. About a week was allowed for animals to recover from the surgery before the first recording session occurred. Usually, 2–3 days of training were sufficient for cats to remain in head-restrained position for 2–4 hr and display several periods of quiet wakefulness, SWS, and REM sleep. The recordings were performed up to 40 days after the surgery.

In vivo recordings in cats

All in vivo recordings were done in a Faraday chamber. LFPs were recorded using tungsten electrodes (2 MΩ, band-pass filter 0.1 Hz to 10 kHz) and amplified with AM 3000 amplifiers (A-M systems) with custom modifications. We aimed to implant electrodes at 1 mm below the cortical surface. A silver wire was fixed either in the frontal bone over the sinus cavity or in the occipital bone over the cerebellum and was used as a reference electrode. All electrical signals were digitally sampled at 20 kHz on Powerlab (ADinstruments, Colorado Springs, USA) and stored for offline analysis. At the end of the experiments, the cats were euthanized with a lethal dose of pentobarbital (100 mg/kg, i.v.).

In vivo experiments in mice

Experiments were conducted on adult C57Bl6 mice. The implantation of recording electrodes was done under isoflurane anesthetized (0.75%–2%). After site shaving, the site of incision was washed with at least three alternating passages of a 4% chlorexidine solution and 70% alcohol. Lidocaine (0.5%) and/or marcaine (0.25%) was injected at the site of incision and at all pressure points. Prior to surgery, the animal was given a dose of buprenorphine (0.05–0.1 mg/kg). During surgery, electrodes for LFP recordings (custom-made, stainless steel) and EMG from neck muscle were implanted and fixed with acrylic dental cement. We aimed at implanting LFP electrodes at a depth of about 650 µm. Four or five anchoring screws were fixed to the cranium. A stainless steel screw implanted above the cerebellum was used as a reference. Throughout the surgery, the body temperature was maintained at 37°C using a water-circulating thermoregulated blanket. Saline (0.1–0.5 ml) and anafen (5 mg/kg) subcutaneous injection were performed at the end of the surgery. After the surgery, mice were given anafen (5 mg/kg) once a day for 2 days. About a week was allowed for animals to recover from the surgery before the first recording session occurred. LFPs were band-pass filtered (0.1–100 Hz). Freely moving mice were plugged to a custom-made headstage to reduce movement artifacts and then amplified with AM 3000 amplifiers (A-M systems). All electrical signals were digitally sampled at 1 kHz on Powerlab (ADinstruments, Colorado Springs, USA) and stored for offline analysis. At the end of the experiments, mice were euthanized with a lethal dose of euthanyl (120 mg/kg, i.v.). The protocol for mice experiments was approved by CPAUL, Comité de protection des animaux de l'Université Laval (Protocol # 2013-014-4).

Data analysis for animal recordings

Electrographic recordings in cats and mice (shown in Figure 6) were analyzed offline using custom-written routines in the IgorPro software. The delta power was calculated from 5 s time windows, using 1 s sliding time window, as the integral power between 0.2 and 4 Hz of the full spectrogram; the spindles power was calculated as the integral power between 8 and 15 Hz. We also computed spindle and delta power using longer time windows (2 s, 5 s, and 30 s). For scatter plots and correlation analysis in Figure 6, spindle and delta power was computed from 30 s windows in 0.2–4 Hz and 8–15 Hz frequency bands. In Figure 7B a shorter time window of 2 s was used. In all cases, we found significant positive correlation between spindle and delta power as reported in the Results section.

ECoG recordings and analysis

Recordings were obtained during natural sleep in a patient undergoing monitoring with intracranial electrodes in order to locate spontaneous seizure onset and guide surgical treatment. Recordings were obtained after informed consent monitored by the Partners Human Subject Protection Committee (Protocol #2007P00165). The patient had an 8 × 8 64-contact peri-sylvian grid and four strips implanted subdurally over the left hemisphere, sampling temporal, Rolandic, parietal and frontal cortices at 1 cm spacing. After rejection of electrodes that either had poor contact or showed epileptic interictal activity, such as interictal spikes, 42 electrodes remained for analysis. We analyzed a recording period of 9 hr at night and sleep scored the patient based on 30 s epochs. Staging was scored using information from scalp, EOG and intracranial electrodes and standard criteria. Awake, REM, N2, and N3 periods were used for further analysis. A fast Fourier transform with a window size of 30 s was performed and power in the SO (0.01–2 Hz) and spindle (9–17 Hz) band frequencies was obtained, and their correlations calculated for 30 s epochs during N2 and N3 periods The average Pearson's correlation coefficient R value was −0.2592 and standard deviation of 0.2281. The average R-squared was 0.1181. 15 electrodes showed significant negative correlation compared to 0 electrodes showing positive correlation after Bonferroni Correction at the p=0.05 level (n = 42, 0.05/42 = 0.0012 = corrected p-value for significance). When focusing on electrodes over the postcentral gyrus (GR7, GR8, GR12, GR14, GR15, GR16, GR21, GR22, GR23, GR24), the average R was -0.4657 with a standard deviation of 0.1721.

Cluster analysis

To compare results of network simulations based on different combinations of model parameters, we performed cluster analysis on the phase locking value (in frequency range 0.5–20Hz) for spindle (7–15Hz) and delta (0.5–4Hz) power. Each trial was 20 s duration; however, the first 10 s was discarded to avoid transient activity. We used Akaike information criteria (AIC) (Akaike, 1974) to measure the least number of components that are required to fit the data. This analysis resulted in a minimum of 10 components beyond which the AIC information measure did not increase. Next a Gaussian mixture model fitting with 10 components using EM method with 500 iterations was applied to obtain the clusters (corresponds to different colored points in Figure 4A). We then determined the projection of the clusters to the space of neuromodulators and finally fitted an ellipsoid to provide an approximate location of the clusters in the neuromodulator space. The ellipsoid fit was performed using geom3d library (https://github.com/dlegland/matGeom/wiki/geom3d). The power, cluster analyses and 3D plotting were performed using MATLAB (MathWorks Inc.). Simulations for various values of neuromodulators were supported by computational resources from the neuroscience gateway (Sivagnanam et al., 2013).

References

  1. 1
  2. 2
  3. 3
  4. 4
  5. 5
  6. 6
    Activity of norepinephrine-containing locus coeruleus neurons in behaving rats anticipates fluctuations in the sleep-waking cycle
    1. G Aston-Jones
    2. FE Bloom
    (1981)
    Journal of Neuroscience 1:876–886.
  7. 7
    Basic Neurochemistry
    1. HA Baghdoyan
    2. R Lydic
    (2012)
    982–999, Chapter 57 - The neurochemistry of sleep and wakefulness, Basic Neurochemistry, 8th edn, New York, Academic Press.
  8. 8
    Spiking-bursting activity in the thalamic reticular nucleus initiates sequences of spindle oscillations in thalamic networks
    1. M Bazhenov
    2. I Timofeev
    3. M Steriade
    4. T Sejnowski
    (2000)
    Journal of Neurophysiology 84:1076–1087.
  9. 9
    Computational models of thalamocortical augmenting responses
    1. M Bazhenov
    2. I Timofeev
    3. M Steriade
    4. TJ Sejnowski
    (1998)
    Journal of Neuroscience 18:6444–6465.
  10. 10
  11. 11
    Model of thalamocortical slow-wave sleep oscillations and transitions to activated States
    1. M Bazhenov
    2. I Timofeev
    3. M Steriade
    4. TJ Sejnowski
    (2002)
    Journal of Neuroscience 22:8691–8704.
  12. 12
  13. 13
  14. 14
  15. 15
  16. 16
  17. 17
  18. 18
  19. 19
  20. 20
  21. 21
  22. 22
  23. 23
    Cellular basis of EEG slow rhythms: a study of dynamic corticothalamic relationships
    1. D Contreras
    2. M Steriade
    (1995)
    Journal of Neuroscience 15:604–622.
  24. 24
  25. 25
  26. 26
    Ionic mechanisms underlying synchronized oscillations and propagating waves in a model of ferret thalamic slices
    1. A Destexhe
    2. T Bal
    3. DA McCormick
    4. TJ Sejnowski
    (1996)
    Journal of Neurophysiology 76:2049–2070.
  27. 27
  28. 28
  29. 29
  30. 30
  31. 31
  32. 32
  33. 33
  34. 34
  35. 35
  36. 36
    The AASM Manual for the Scoring of Sleep and Associated Events: Rules, Terminology and Technical Specifications: American Academy of Sleep Medicine
    1. C Iber
    2. M AAoS
    (2007)
    The AASM Manual for the Scoring of Sleep and Associated Events: Rules, Terminology and Technical Specifications: American Academy of Sleep Medicine.
  37. 37
  38. 38
  39. 39
  40. 40
  41. 41
  42. 42
  43. 43
  44. 44
  45. 45
  46. 46
  47. 47
  48. 48
    Modulation of neuronal firing mode in cat and guinea pig LGNd by histamine: possible cellular mechanisms of histaminergic control of arousal
    1. DA McCormick
    2. A Williamson
    (1991)
    Journal of Neuroscience 11:3188–3199.
  49. 49
  50. 50
  51. 51
  52. 52
    Muscarinic inhibition of persistent Na+ current in rat neocortical pyramidal neurons
    1. T Mittmann
    2. C Alzheimer
    (1998)
    Journal of Neurophysiology 79:1579–1582.
  53. 53
  54. 54
  55. 55
    Grouping of spindle activity during slow oscillations in human non-rapid eye movement sleep
    1. M Mölle
    2. L Marshall
    3. S Gais
    4. J Born
    (2002)
    Journal of Neuroscience 22:10941–10947.
  56. 56
    A taxonomic analysis of sleep stages
    1. B Müller
    2. WD Gäbelein
    3. H Schulz
    (2006)
    Sleep 29:967–974.
  57. 57
  58. 58
    Electroencephalography: Basic Principles, Clinical Applications, and Related Fields
    1. E Niedermeyer
    2. FL da Silva
    (2005)
    Electroencephalography: Basic Principles, Clinical Applications, and Related Fields, Lippincott Williams & Wilkins.
  59. 59
  60. 60
    Replay and time compression of recurring spike sequences in the hippocampus
    1. Z Nádasdy
    2. H Hirase
    3. A Czurkó
    4. J Csicsvari
    5. G Buzsáki
    (1999)
    Journal of Neuroscience 19:9497–9507.
  61. 61
  62. 62
  63. 63
  64. 64
    Elevated levels of histamine metabolites in cerebrospinal fluid of aging, healthy humans
    1. GD Prell
    2. JK Khandelwal
    3. RS Burns
    4. PA LeWitt
    5. JP Green
    (1988)
    Comprehensive Gerontology. Section A, Clinical and Laboratory Sciences 2:114–119.
  65. 65
  66. 66
    A Manual of Standardized Terminology, Techniques and Scoring System for Sleep Stages of Human Subjects
    1. A Rechtschaffen
    2. A Kales
    (1968)
    Bethesda, Md: U. S. National Institute of Neurological Diseases and Blindness, Neurological Information Network.
  67. 67
    Quantitative modelling of sleep dynamics
    1. PA Robinson
    2. AJK Phillips
    3. BD Fulcher
    4. M Puckeridge
    5. JA Roberts
    (2011)
    Philosophical Transactions of the Royal Society A: Mathematical, Physical and Engineering Sciences 369:3840–3854.
    https://doi.org/10.1098/rsta.2011.0120
  68. 68
  69. 69
  70. 70
  71. 71
  72. 72
  73. 73
    The visual scoring of sleep in adults
    1. MH Silber
    2. S Ancoli-Israel
    3. MH Bonnet
    4. S Chokroverty
    5. MM Grigg-Damberger
    6. M Hirshkowitz
    7. S Kapen
    8. SA Keenan
    9. MH Kryger
    10. T Penzel
    11. MR Pressman
    12. C Iber
    (2007)
    Journal of Clinical Sleep Medicine 3:121–131.
  74. 74
  75. 75
  76. 76
    Abolition of spindle oscillations in thalamic neurons disconnected from nucleus reticularis thalami
    1. M Steriade
    2. M Deschênes
    3. L Domich
    4. C Mulle
    (1985)
    Journal of Neurophysiology 54:1473–1497.
  77. 77
    Network modulation of a slow intrinsic oscillation of cat thalamocortical neurons implicated in sleep delta waves: cortically induced synchronization and brainstem cholinergic suppression
    1. M Steriade
    2. RC Dossi
    3. A Nuñez
    (1991)
    Journal of Neuroscience 11:3200–3217.
  78. 78
    Thalmic Oscillations and Signaling
    1. M Steriade
    2. EG Jones
    3. R Llinas
    (1990)
    New York: John Wiley & Sons.
  79. 79
    The functional states of the thalamus and the associated neuronal interplay
    1. M Steriade
    2. RR Llinás
    (1988)
    Physiological Reviews 68:649–742.
  80. 80
  81. 81
    Intracellular analysis of relations between the slow (< 1 Hz) neocortical oscillation and other sleep rhythms of the electroencephalogram
    1. M Steriade
    2. A Nuñez
    3. F Amzica
    (1993b)
    Journal of Neuroscience  13:3266–3283.
  82. 82
    Natural waking and sleep states: a view from inside neocortical neurons
    1. M Steriade
    2. I Timofeev
    3. F Grenier
    (2001)
    Journal of Neurophysiology 85:1969–1985.
  83. 83
  84. 84
  85. 85
  86. 86
    Low-frequency rhythms in the thalamus of intact-cortex and decorticated cats
    1. I Timofeev
    2. M Steriade
    (1996)
    Journal of Neurophysiology 76:4152–4168.
  87. 87
  88. 88
  89. 89
    Basal forebrain acetylcholine release during REM sleep is significantly greater than during waking
    1. J Vazquez
    2. HA Baghdoyan
    (2001)
    American Journal of Physiology. Regulatory, Integrative and Comparative Physiology 280:R598–601.
  90. 90
  91. 91
  92. 92
  93. 93
  94. 94

Decision letter

  1. Ronald L Calabrese
    Reviewing Editor; Emory University, United States

In the interests of transparency, eLife includes the editorial decision letter and accompanying author responses. A lightly edited version of the letter sent to the authors after peer review is shown, indicating the most substantive concerns; minor comments are not usually included.

Thank you for submitting your article "Cellular and Neurochemical Basis of Sleep Stages in the Thalamocortical network" for consideration by eLife. Your article has been reviewed by three peer reviewers, one of whom, Ronald L Calabrese (Reviewer #1), is a member of our Board of Reviewing Editors and the evaluation has been overseen by David Van Essen as the Senior Editor. The following individuals involved in review of your submission have agreed to reveal their identity: Paul Garcia (Reviewer #2); Miles Whittington (Reviewer #3).

The reviewers have discussed the reviews with one another and the Reviewing Editor has drafted this decision to help you prepare a revised submission.

Summary:

In this very interesting study that combines computational modeling and some in vivo recordings of brain activity during sleep and wakefulness and under ketamine anesthesia, the authors propose a comprehensive model of how major neuromodulators may control sleep states and their transitions in mammals. The model is grounded in previous experimental evidence and covers the major neuromodulators ACh, GABA, and HA. They show how these modulators might synchronize the spindle (N2) and slow oscillations, SO, (N3) states, as well as REM. Using ECoG recordings from humans and LFP recordings from cats and mice, during NREM sleep they show that the power of spindle and SO is negatively correlated in human and positively correlated in cats and mice, and explain this discrepancy by the differences in the relative levels of ACh. They also explore similarities and differences between SO during N3 and ketamine anesthesia and explain them in terms of influences of ketamine on neuromodulator levels and effects. The study identifies potential intrinsic and synaptic mechanisms through which neuromodulators acting in combination may mediate transitions between sleep stages. These findings should be of wide interest to the sleep, and anesthesia communities and to those interested in coordinated whole brain activity.

Essential revisions:

There are some major concerns about claims but overall this appears to be a strong study of general significance.

1) The authors have gone too far when suggesting that they are modelling the sleep-related effects of the neuromodulators they have chosen. The authors are really basing their conclusions on an oversimplified (and under-developed) link between ACh and effects on membrane currents. Their models do show remarkable state changes with manipulation of the 4 main parameters (K-leak, Ih, GABA, AMPA) and these finding are both interesting and challenging to the field – but the relation back to ACh and HA is not as strong as they imply. There are many precedents for ACh and HA having the opposite effect on the main model parameters on certain cell subtypes cf. the ones they use.

The authors should clarify that they are making a selection of potential modulator effects, and in Discussion suggest this MAY relate to known endogenous brain-state modulators like ACh, Histamine (and many others).

Alternatively the authors can revisit their model and represent more accurately the known effects of HA and ACh on each of the 'loose' cell types included in the model.

2) A similar concern applies to the handling of data gathered with ketamine anesthesia. The authors attempt to mimic ketamine in their model by decreasing NMDA's influence – but ketamine is a promiscuous molecule that has influence on a great number of channels. Additionally, the unconsciousness produced by the "dissociative" anesthetic, ketamine, has a very different phenotype than that of propofol or isoflurane which more or less produce a quiescence resembling sleep. As opposed to these GABAergic anesthetics, ketamine does not activate the sleep nuclei of the hypothalamus (VLPO), does not depress thalamic activity, and activates wake-promoting nuclei in the brain stem. The adjunct agent used in their in vivo experiments (xylazine) does on the other hand promote sleep, but does not do this through NMDA antagonism. Therefore the specificity of their model is really in question. George Mashour has done work on ketamine and his work might be better referenced. http://www.ncbi.nlm.nih.gov/pmc/articles/PMC4076669/

Further, we would expect that an increase in inhibition in the model would produce slow oscillations similar to sleep. This model state could be compared to in vivo data using propofol or isoflurane. Perhaps, this model is "tuned" to produce these stable oscillations in a variety of different parameter regimes – but as mentioned in the expert reviews (Reviewer #2), without an evaluation of the stimulation protocols it is difficult to determine what is a network effect vs. a reflection of the periodicity already present in the stimulation.

The expert reviews of reviewers #2 and #3 are appended to aid in the revision.

Reviewer #2:

This manuscript reports the results of a computational model of sleep cycling based on variations in neurotransmitters. The model was also compared to neuronal oscillations during sleep using human, cat and mouse datasets. The main conclusion is that the initiation of the oscillatory pattern indicative of N2 sleep can be simulated by a decrease in histamine and acetylcholine. A further decrease in acetylcholine can also explain a transition to deeper sleep. In general, this model contributes much to our understanding of the neurochemical basis of sleep macroarchitecture. Additionally, the results challenge theories that increases in the inhibitory neurotransmitter GABA mediate transitions among sleep stages. Although some of the conclusions need to be amended and there are several places that need some grammatical corrections, there is much to be admired in this paper. However some of the grammatical errors and non-sequitors (non-explored ideas) make me concerned that it was hastily prepared.

Summary:

– Some parts of the summary are awkward in regard to sentence structure and punctuation. The second sentence does not sound like proper subject verb agreement, and the third sentence requires the word and between the two subjects to make sense. Α/theta bursts are mentioned in the Summary but not emphasized in the paper (previous version?) please remove.

Introduction:

– First sentence of second paragraph should be reworded so characterize is not repeated. The first sentence of the third paragraph should start with: "Neuromodulators".

Results

– The TC stimulation is not represented in Figure 2. This input is also not described in the methods, given the strong excitatory connections between thalamus and cortex the nature of the stimulation ("simulating sensory input") could be driving the oscillations and must be included to evaluate the main conclusions of the paper.

– Why did they examine such a large frequency range to evaluate phase-locking? When the connections are not completely known (in vivo recordings) phase-locking among narrow frequency bands (especially among specific regions, e.g., thalamus and cortex) is much more important for synchronization. Very little of the quantitative results rely on phase-locking.

– Why is ketamine compared to SWS? Others report it as sharing similarities with REM. Propofol or isoflurane anesthesia would be a better comparison.

Methods

– Extrasynaptic and Intrasynaptic GABA concentrations are mentioned but not explored in the paper, please remove (unless added to discussion as to a potential limitation of the model).

– Network description: radius is not defined explicitly.

– The in vivo description mentions non-anesthetized cats and ketamine is given?

Discussion

–- Much of the Discussion can be considerably shortened.

– The Discussion does not adequately address other potential influences on spindle power (e.g., age and cortical volume).

Reviewer #3:

This is a very interesting paper that uses a fairly detailed computational model, constrained by experimental observations in a variety of species, which is used to demonstrate the dependence of sleeps stage on a variety of network and intrinsic cell properties (leak current, AMPA, GABA, I(h)). The paper extends the findings to link sleep stages to two main neuromodulators known to be involved in controlling sleep – Ach and Histamine.

In general the paper is well written and the data beautifully presented. I very much liked the demonstration of different sleep stages with different intrinsic properties but it would have been good to refer to the very detailed model of Hill & Tononi which, a while ago, came to very similar conclusions.

My only major problem with the manuscript comes from the attempt to relate the model parameter changes to histamine and Acetylcholine. I think a bit more detail and biological realism needs to be included here. For example, Ach is stated to affect AMPA-mediated transmission and leak current. Of the parameters used to construct the model cells it also affects Km and GABA release. Reading the methods it seems Ach was modelled as a scale factor for all intrinsic and synaptic properties – but this scale factor was increased for decreased Ach levels. This means, for example, that properties shown not to be affected by Ach were affected in the model, and properties of Ach well documented (such as depolarisation of main interneuron subtypes) were omitted or modelled conversely.

On a similar note, it is not clear how the 2 modulators chosen fit with the GABA levels used. There is a wealth of evidence linking Ach and Histamine to GABA release and postsynaptic effects – some which fit the basic hypothesis the authors propose (decreased Ach, increased GABA release), and some antagonistic to it (increased histamine, increased tonic GABA). This could be clarified with a more rigorous modelling of the known effects of the two agents to the conductances used in the model. If no changes in sleep stages are seen then this would provide evidence for an involvement of the other major neuromodulators involved in sleep not considered here.

[Editors' note: further revisions were requested prior to acceptance, as described below.]

Thank you for resubmitting your work entitled "Cellular and Neurochemical Basis of Sleep Stages in the Thalamocortical network" for further consideration at eLife. Your revised article has been favorably evaluated by David Van Essen (Senior editor), a Reviewing editor, and two reviewers.

The reviewers consider the manuscript to be much improved. There is only one minor issue that you may want to address before final acceptance, as outlined below:

Figure 1C and D Vertical axes might be clearer if they read "Normalized Conductance (%)" It should be obvious they are normalized to the Awake condition by the graph.

Reviewer #2:

This is a re-submission on a computational model of transitions among sleep stages that is firmly based in neurophysiology. The model is compared to neurophysiologic recordings from LFP records obtained in animal models and human ECoG data. Finally, the model's predictions are tested with simulations designed to mimic propofol anesthesia. This revised manuscript has much improved readability in both the precision of writing and clarity of thought. The main conclusions are that a reduction in acetylcholine and histamine levels in the thalamus corresponds to N2 sleep, characterized by spindle activity. And that an increase in GABA alone does not produce the oscillations characteristic of N2. The experiments with propofol anesthesia demonstrate this nicely. The model also predicts that reduced cortical acetylcholine influences the transition to N3 sleep. Overall, this is a nice contribution to the sleep and anesthesia literature.

Reviewer #3:

This one of the best, most comprehensive revisions of a paper I have seen. The authors have taken on every major comment and provided a great deal more clarity, focus and detail. I have no hesitation in recommending it for publication now.

https://doi.org/10.7554/eLife.18607.015

Author response

[…]

Essential revisions:

There are some major concerns about claims but overall this appears to be a strong study of general significance.

1) The authors have gone too far when suggesting that they are modelling the sleep-related effects of the neuromodulators they have chosen. The authors are really basing their conclusions on an oversimplified (and under-developed) link between ACh and effects on membrane currents. Their models do show remarkable state changes with manipulation of the 4 main parameters (K-leak, Ih, GABA, AMPA) and these finding are both interesting and challenging to the field – but the relation back to ACh and HA is not as strong as they imply. There are many precedents for ACh and HA having the opposite effect on the main model parameters on certain cell subtypes cf. the ones they use.

The authors should clarify that they are making a selection of potential modulator effects, and in Discussion suggest this MAY relate to known endogenous brain-state modulators like ACh, Histamine (and many others).

As requested, we made a number of changes to the text to clearly indicate that we only tested a subset of the effects of neuromodulators and many known effects of neuromodulators were not included in our model. As pointed out by the reviewer, the four main parameters (K-leak, AMPA, Ih, GABA) were sufficient to explain transitions between sleep stages in the model. These cellular and synaptic properties were chosen based on the past experiments where they have been linked to the effects of neuromodulators. This putative (and certainly incomplete) link allowed us to interpret our results and to make specific predictions for future studies on the variations in the neuromodulators during different sleep stages and on the differences in the neuromodulatory response between humans and cats/mice. To further address reviewer concerns, we also examined transition between sleep stages in a more detailed model (see our response below), which revealed similar qualitative results as the reduced “four-parameter” model. The way we want to interpret this is that a simple model identifies the core or minimally required changes, which are likely to be linked to the known effects of the neuromodulators in the thalamocortical network. We include the relevant discussion in the modified manuscript, please see below.

The text in the result section was modified as follows:

“While our model does not capture the entire spectrum of changes associated with known effects of neuromodulators, we identified a minimal set sufficient to account for characteristic changes of brain electrical activity across the sleep-wake cycle. […] We also tested in the model the effect of extracellular GABA concentration on tonic inhibition.”

The following paragraph was included to the Discussion section:

“In our study, the effects of neuromodulators were limited by changing the following four model properties – maximal conductances of the K+-leak, Ih, GABA and AMPA currents. […] In contrast, after removing the changes of the K+-leak, Ih, GABA and AMPA currents but leaving the action of neuromodulators on the other channels, the detailed model failed to show transitions between sleep stages (Figure 11B), suggesting that the cellular and synaptic properties identified in our study are critical in determining patterns of electrical activity in the thalamocortical network. “

Alternatively the authors can revisit their model and represent more accurately the known effects of HA and ACh on each of the 'loose' cell types included in the model.

As suggested by the reviewer, we examined a more detailed model to include additional effects of neuromodulators on the intrinsic currents that are specific to cell types. In these simulations, ACh influenced the K+-leak, M-current, H-current, persistent Na+ current and Ca2+-dependent K+ current in the cortical neurons. Histamine levels determined the K+-leak current and H-currents in all neuron types. The main qualitative predictions of the detailed model were similar to the reduced model. The four core model properties (K+-leak, Ih, GABA and AMPA) were critical to model transitions between sleep stages and selectively removing any of them, as shown in the Figure 5, led to the failure to explain sleep stage transitions as observed in vivo. Below we present the figure (now Figure 11) showing state transitions in a detailed model (Figure 11A) and the failure of modeling transitions between sleep stages when the core properties were prevented to change in the detailed model (Figure 11B).

We would like to note that the distributions of the intrinsic and synaptic currents in our model were varied across cell types. For instance, the cortical inhibitory neurons had lower conductance or absence of M-current, persistent Na+ and Ca2+-dependent K+ currents. In thalamic cells, intrinsic currents were different and were based on the data. Thus, the same neuromodulators influenced these different neuron types differently in our network model.

The following changes were implemented to the text:

We included Figure 11 which presents more detailed model of the effects of neuromodulators.

We also included the following text relevant to this figure in the Discussion section:

“We also tested a more detailed model, implementing the effects of ACh and HA on the other K+ currents as well as on the persistent Na+ and M-currents (Constanti and Galvan, 1983; McCormick and Williamson, 1989; Mittmann and Alzheimer, 1998). We observed that this detailed model qualitatively captured results similar to the minimal model (Figure 11A). In contrast, after removing the changes of the K+-leak, Ih, GABA and AMPA currents but leaving the action of neuromodulators on the other channels, the detailed model failed to show transitions between sleep stages (Figure 11B), suggesting that the cellular and synaptic properties identified in our study are critical in determining patterns of electrical activity in the thalamocortical network.”

2) A similar concern applies to the handling of data gathered with ketamine anesthesia. The authors attempt to mimic ketamine in their model by decreasing NMDA's influence – but ketamine is a promiscuous molecule that has influence on a great number of channels. Additionally, the unconsciousness produced by the "dissociative" anesthetic, ketamine, has a very different phenotype than that of propofol or isoflurane which more or less produce a quiescence resembling sleep. As opposed to these GABAergic anesthetics, ketamine does not activate the sleep nuclei of the hypothalamus (VLPO), does not depress thalamic activity, and activates wake-promoting nuclei in the brain stem. The adjunct agent used in their in vivo experiments (xylazine) does on the other hand promote sleep, but does not do this through NMDA antagonism. Therefore the specificity of their model is really in question. George Mashour has done work on ketamine and his work might be better referenced.http://www.ncbi.nlm.nih.gov/pmc/articles/PMC4076669/

Further, we would expect that an increase in inhibition in the model would produce slow oscillations similar to sleep. This model state could be compared to in vivo data using propofol or isoflurane. Perhaps, this model is "tuned" to produce these stable oscillations in a variety of different parameter regimes – but as mentioned in the expert reviews (Reviewer #2), without an evaluation of the stimulation protocols it is difficult to determine what is a network effect vs. a reflection of the periodicity already present in the stimulation.

We again thank the reviewer for this valuable comment. We choose to model Ketamine condition based on experimental evidences suggesting that ketamine inhibits NMDA receptors and ACh nicotinic and muscarinic receptors (Rudolph and Antkowiak, 2004; Alkire et al., 2008). However, we agree with the reviewer that ketamine could have many other effects, which are not fully understood. Thus, we replaced the section on the effect of ketamine with a simulation of the propofol anesthesia, as suggested by the reviewer.

Propofol increases decay time constant of the inhibitory post-synaptic potential (Kitamura et al., 2003). It is also known to reduce levels of ACh (Kikuchi et al., 1998) and inhibit TMN activity which is the source for HA (Nelson et al., 2002). In both animals and humans, propofol has shown to increases 8-15 Hz activity in the frontal regions and appearance of slow oscillations (Murphy et al., 2011; Purdon et al., 2013). In agreement with these previous studies, we present below an example of our in-vivo recordings from cat (this was a pilot experiment with N=1, so we decided not to include these data to the main text), that demonstrate increase in 8-15 Hz oscillations following the application of propofol (Author response image 1). Note that the amplitude of 8-15 Hz oscillations in propofol condition was higher compared to natural sleep.

Based on these in vivo studies, we simulated effect of propofol by (a) increasing decay time constant of synaptic GABAergic inhibition as well as (b) modeling inhibition of ACh and HA release. Following such manipulations, we observed an increase in the 8-15 Hz activity (compare to our baseline model of the natural sleep). We then compared effect of increasing GABA time constant alone versus combined action of increasing GABA time constant and effect of reduction of ACh and HA. In the first case, we found no slow oscillations or oscillatory activity in 8-15 Hz. These findings suggest that propofol’s action may not be fully explained by its effect on synaptic inhibition and considering its indirect effect on release of other neuromodulator is required to observe the changes seen in vivo in the thalamocortical network.

Author response image 1
Multichannel recording from a cat during natural sleep (control) and propofol induced anesthesia condition.

Top plot shows the recording electrode and picture of the implanted electrode. Middle plot shows the activity during natural sleep (raw data on top and filtered activity on bottom). Bottom plot show the activity during propofol anesthesia.

https://doi.org/10.7554/eLife.18607.014

We now included the section “Modeling the effects of propofol anesthesia requires the combined action of ACh, HA and GABA”.

We have also replaced Figure 9.

In response to the reviewer comment about external stimulation to the model (please also see 2d reviewer comment below), we need to state that we have not been clear about nature of the stimulation in the model. We never applied an external persistent periodic or DC stimulation in any of the model simulations presented in the paper. All the model activity reported in this study results from intrinsic dynamics of the thalamocortical network. The only stimulation we applied, and what was mentioned in the original text, was a brief DC pulse to probe stability of the network. We have clarified this issue in the revised version of the manuscript.

To clarify the stimulation protocol, the following text was updated:

“During awake period, activity in the cortical PY and IN and thalamic RE cells was sparse, i.e., mean firing was around 2-3 Hz and 4-5 Hz in PY and IN cells, respectively, […] then the network returned to the baseline, as observed in experiments with a sensory input during alert awake conditions (Bruno and Sakmann, 2006; Hengen et al., 2016).”

And we described the stimulation in the Methods:

“Stimulation Protocol: A brief DC input lasting for 100 ms and amplitude of 1.5 nA was applied to TC neurons during awake state in some simulations to test stability of the network dynamics. The DC pulse was modeled using two Heaviside functions. No periodic stimulation was applied in any of the modeling experiments.”

The expert reviews of reviewers #2 and #3 are appended to aid in the revision.

Reviewer #2:

[…]

Summary:

- Some parts of the summary are awkward in regard to sentence structure and punctuation. The second sentence does not sound like proper subject verb agreement, and the third sentence requires the word and between the two subjects to make sense. Α/theta bursts are mentioned in the Summary but not emphasized in the paper (previous version?) please remove.

Thank you for pointing us to those issues. The summary section was edited to address reviewer’s critics and to clarify presentation. We also removed the mentioning of the α bursts in the summary. In the original manuscript, we still discuss α bursts when discussing the activity during REM sleep, however, it is scattered across the text and, we agree, that it is not significant enough to be mentioned in the summary.

Introduction:

– First sentence of second paragraph should be reworded so characterize is not repeated. The first sentence of the third paragraph should start with: "Neuromodulators".

We thank the reviewer for pointing out these changes. We revised the manuscript as it was suggested.

Results

– The TC stimulation is not represented in Figure 2. This input is also not described in the methods, given the strong excitatory connections between thalamus and cortex the nature of the stimulation ("simulating sensory input") could be driving the oscillations and must be included to evaluate the main conclusions of the paper.

Unfortunately, we missed to describe stimulation in sufficient details. We never applied an external persistent periodic or DC stimulation in any of the model simulations presented in the paper. All the model activity reported in this study results from intrinsic dynamics of the thalamocortical network. The TC stimulation, mentioned in the text, was applied as a brief DC pulse during awake stage to demonstrate the stability of the awake state activity. There was no external stimulation applied at any other time during the simulations. Below is the figure showing the response to brief TC stimulation.

We modified the text in the modified manuscript to clarify this further as shown below.

“During awake period, activity in the cortical PY and IN and thalamic RE cells was sparse, i.e., mean firing was around 2-3 Hz and 4-5 Hz in PY and IN cells, respectively, while the TC neurons primarily remained silent and RE cells were busting intermittently (Figure 2A and Figure 3A). To test stability of the network dynamics, we applied thalamic stimulation in the form of a brief (100 msec) DC pulse once during awake state. It triggered transient increase of firing in cortical PY neurons and INs, and then the network returned to the baseline, as observed in experiments with a sensory input during alert awake conditions (Bruno and Sakmann, 2006; Hengen et al., 2016).”

We also added to the Discussion:

“All the network rhythms reported in this study result from intrinsic (autonomous) dynamics of the thalamocortical network; no external stimulation was applied to entrain any of the oscillatory activities in the model.”

And we described the stimulation in the Methods:

“Stimulation Protocol: A brief DC input lasting for 100 ms and amplitude of 1.5 nA was applied to TC neurons during awake state in some simulations to test stability of the network dynamics. The DC pulse was modeled using two Heaviside functions. No periodic stimulation was applied in any of the modeling experiments.”

– Why did they examine such a large frequency range to evaluate phase-locking? When the connections are not completely known (in vivo recordings) phase-locking among narrow frequency bands (especially among specific regions, e.g., thalamus and cortex) is much more important for synchronization. Very little of the quantitative results rely on phase-locking.

We thank the reviewer for point to this issue. It was a typo in our manuscript and the phase locking was computed for 0.5-20 Hz and not for 0.5-100 Hz. The only part of the result that utilized this phase locking analysis was for identifying the clusters of activity in Figure 4A and the 4B. Using a range from 0.5 to 20 Hz allowed for better detection of the clusters. We found that using very narrow power bands corresponding to slow oscillation (0.5-4Hz) or spindle activity (7-15 Hz) resulted in less clear distinction between the clusters and adding a dimension by using separate phase-locking value for δ and spindle bands did not improve the separation of clusters or a change in AIC curves. We now corrected the text to indicate the correct frequency range used for the phase locking.

– Why is ketamine compared to SWS? Others report it as sharing similarities with REM. Propofol or isoflurane anesthesia would be a better comparison.

Methods

We thank the reviewer for this important comment. Similar critique was also found in the reviewer 1 comments. We now removed the results on the ketamine and instead present results on the propofol anesthesia condition. We included new figure and a new section to discuss these new results. All the details about this section are presented in our response to question 2 of the “Essential revisions” (please see above).

– Extrasynaptic and Intrasynaptic GABA concentrations are mentioned but not explored in the paper, please remove (unless added to discussion as to a potential limitation of the model).

Indeed, in the Discussion section of the original manuscript we had a brief discussion of the Extrasynaptic vs Intrasynaptic GABA effects. Since the origin of change in GABA is not fully understood, we had to examine various possibilities whether the measured GABA levels reflect the extrasynaptic GABA concentration or the tonic GABA (Figure 10). As pointed by the reviewer, this is one of the limitations of the model due to the lack of clear experimental evidence on this topic. Therefore, we kept this text in the revised manuscript, however, we emphasized the reasoning for discussing this topic along with new supplementary figure. A detailed discussion on this topic was also provided to response to reviewer 3 (see below, last comment from reviewer 3). Particularly, the modified text was included in subsection “Role of GABA in sleep stage transitions.”

– Network description: radius is not defined explicitly.

Thank you. Similar comment was also made by Reviewer 1. We now modified the method section “Network structure” to include the details about network connectivity and radius.

– The in vivo description mentions non-anesthetized cats and ketamine is given?

We now removed the section on Ketamine (both animal data and model) and instead reported the modeling results on propofol induced anesthesia. Thus, in this revised version, all the cats recordings are from non-anesthetized animals.

Discussion

– Much of the Discussion can be considerably shortened.

We thank the reviewer for this suggestion. We removed several paragraphs from the Discussion section.

- The Discussion does not adequately address other potential influences on spindle power (e.g., age and cortical volume).

Indeed, many factors including age and cortical volume may impact the spindle and slow oscillation power. While the model is too simplified to address potential impact of the cortical volume, we included a brief discussion about possible changes with aging as shown below:

Furthermore, our results on the putative relationship between levels of neuromodulators and characteristic sleep electrical activities could be potentially relevant to understanding the changes observed with aging. Indeed, the level of the HA in cerebrospinal fluid is known to be elevated with aging (Prell et al., 1988), where spindle power is known to be reduced (Crowley et al., 2002). In our model, the HA level is one of the major predictors of spindle power. Thus, we can speculate that a reduction of HA with age could contribute to the decline of the spindle power.

Reviewer #3:

This is a very interesting paper that uses a fairly detailed computational model, constrained by experimental observations in a variety of species, which is used to demonstrate the dependence of sleeps stage on a variety of network and intrinsic cell properties (leak current, AMPA, GABA, I(h)). The paper extends the findings to link sleep stages to two main neuromodulators known to be involved in controlling sleep – Ach and Histamine.

In general the paper is well written and the data beautifully presented. I very much liked the demonstration of different sleep stages with different intrinsic properties but it would have been good to refer to the very detailed model of Hill & Tononi which, a while ago, came to very similar conclusions.

We greatly appreciate the reviewer’s positive comments about the manuscript. We thank the reviewer for recognizing the quality of our results and our approach. The reviewer’s critique has helped us make the manuscript more accessible to readers.

We thank the reviewer for pointing us about missing Hill et al. 2005 citation. Unfortunately, it disappeared by mistake after one of the internal revisions. We now cite this paper in the revised manuscript. In our new study we extend many of the early findings reported in Hill et al. 2005. Main advances of our new study include: (a) modeling all major sleep stages (including stage 2 and REM sleep) based on the biophysical mechanisms; (b) identifying role of the many other intrinsic and synaptic properties in sleep stage transitions (e.g., role of AMPA strength or Ih conductance strength); (c) reporting intermediate stages with spindle and slow oscillations that resulted from changing both ACh and HA; this was required to replicate the recordings from humans and animals.

We included this relevant citation in our discussion as follows:

“Previous computational studies proposed a minimal set of the mechanisms sufficient to model sleep spindles and slow oscillations (Destexhe et al., 1996; Bazhenov et al., 1998, 1999; Timofeev et al., 2000; Bazhenov et al., 2002; Compte et al., 2003b; Hill and Tononi, 2005; Bonjean et al., 2011). An increase in the K+ leak current was identified as a critical component for the transition between awake and slow wave sleep (Bazhenov et al., 2002; Hill and Tononi, 2005), and it was predicted that synchronization of the slow oscillation depends on cortico-cortical connections (Hill and Tononi, 2005). Several phenomenological and reduced mathematical models have been proposed to test the effects of neuromodulators on the sleep-wake cycle (McCarley and Hobson, 1975) (Robinson et al., 2011). Nevertheless, none of these previous models examined transitions among all major sleep stages based on biophysical mechanisms.”

My only major problem with the manuscript comes from the attempt to relate the model parameter changes to histamine and Acetylcholine. I think a bit more detail and biological realism needs to be included here. For example, Ach is stated to affect AMPA-mediated transmission and leak current. Of the parameters used to construct the model cells it also affects Km and GABA release. Reading the methods it seems Ach was modelled as a scale factor for all intrinsic and synaptic properties – but this scale factor was increased for decreased Ach levels. This means, for example, that properties shown not to be affected by Ach were affected in the model, and properties of Ach well documented (such as depolarisation of main interneuron subtypes) were omitted or modelled conversely.

We believe that we failed to present model equations properly in the Methods so it made an impression that”Ach was modelled as a scale factor for all intrinsic and synaptic properties”. Among all intrinsic currents ACh had only influenced the K+ leak current and had no effect on the other intrinsic currents. It also had effect on the synaptic AMPA-type transmission. We now modified our equations in Methods to clearly show the effect of ACh. Further, we would like to note that different cell types had different strength of intrinsic currents, including K+ leak, thus the same neuromodulator had varying effects on the different neuron types based on the cell specific combination of intrinsic currents.

The Methods section reads now:

“The cortical neuron model included dendritic and axo-somatic compartments, similar to the previous studies (Timofeev et al., 2000; Bazhenov et al., 2002); (Chen et al., 2012) and was described by the following equations,

CmdVddt=[AchPY.IdKleak]IdleakIdNaIdNapIdCaIdKCaIdKmIsym
gcs(VdVs)=IsNaIsKIsNap

where the subscripts s and d correspond to axo-somatic and dendritic compartments, is the variable which is inversely proportional to the level of ACh in the cortex (AchPY=1 for awake state, AChPY= 1.25 for N2, AChPY=1.8 for N3 and =0.85 during REM) and determines the strength of the potassium leak current (IdKleak)[…]”

We would like to add that, as also was suggested by the first reviewer, we now included the results from a detailed model (Figure 11), which implements the action of ACh and HA on the all K+ currents, Ih and persistent Na+ current. While these additional actions of ACh and HA were not sufficient by themselves to cause transitions between sleep stages, including these effects kept the qualitative features of the stage transitions that were reported in the minimal model (Figure 11).

On a similar note, it is not clear how the 2 modulators chosen fit with the GABA levels used. There is a wealth of evidence linking Ach and Histamine to GABA release and postsynaptic effects – some which fit the basic hypothesis the authors propose (decreased Ach, increased GABA release), and some antagonistic to it (increased histamine, increased tonic GABA). This could be clarified with a more rigorous modelling of the known effects of the two agents to the conductances used in the model. If no changes in sleep stages are seen then this would provide evidence for an involvement of the other major neuromodulators involved in sleep not considered here.

We again thank the reviewer for this important comment. Indeed, the origin of the changes in the GABA level during sleep is not well understood. In this study, we assumed that phasic GABA release is increased during stage 2 and 3 sleep, when ACh and HA are reduced. However, as suggested by the reviewer, there is a complex interaction between GABA and other neuromodulators. In our study, we examined in detail how varying levels of GABA influence the activity in the thalamocortical network for different combinations of ACh and HA (Figure 4). We identified that change in the level of GABA has relatively less effect compared to the changes in ACh and HA (note the clusters in Figure 4B do not change a lot in the z-axis as much as it changes in x and y-axis). However, for intermediate levels of ACh and HA, GABA determined the density of spindling activity (note the cluster 1,5 and 6 are located at the higher GABA levels). Further, we tested an alternative idea that experimental measurements of GABA reflect the extrasynaptic GABA (see figure below and Figure 10), which reflects amount of tonic inhibition. In these simulations, we found that changes in tonic inhibition did not alter the transition between stages but there were changes in amount of oscillatory activity in N2, N3 and REM stages. Due to the lack of clear experimental evidence on the origin of GABA changes across stages of vigilance, we were unable to make clear predictions about role of GABA in sleep stage transitions. Nevertheless, our results suggest that GABA through phasic and tonic inhibition alter the nature of activity within each sleep stage in the thalamocortical network. Furthermore, by examining the variations across wide range of GABA levels (eg in Figure 4) and in conditions of propofol anesthesia (Figure 9) our study could help guiding future experimental studies.

We now include text related to this issue in subsection “Role of GABA in sleep stage transitions.”

https://doi.org/10.7554/eLife.18607.016

Article and author information

Author details

  1. Giri P Krishnan

    Department of Medicine, University of California, San Diego, La Jolla, CA, United States
    Contribution
    GPK, Conception and design, Analysis and interpretation of data, Drafting or revising the article
    For correspondence
    gkrishnan@ucsd.edu
    Competing interests
    The authors declare that no competing interests exist.
    ORCID icon "This ORCID iD identifies the author of this article:" 0000-0002-3931-7633
  2. Sylvain Chauvette

    1. Department of Psychiatry and Neuroscience, Université Laval, Québec, Canada
    2. Centre de Recherche de l’Institut Universitaire en Santé Mentale de Québec, Université Laval, Québec, Canada
    Contribution
    SC, Acquisition of data, Analysis and interpretation of data, Drafting or revising the article
    Competing interests
    The authors declare that no competing interests exist.
  3. Isaac Shamie

    Departments of Radiology and Neurosciences, University of California, San Diego, La Jolla, CA, United States
    Contribution
    IS, Acquisition of data, Analysis and interpretation of data, Drafting or revising the article
    Competing interests
    The authors declare that no competing interests exist.
  4. Sara Soltani

    1. Department of Psychiatry and Neuroscience, Université Laval, Québec, Canada
    2. Centre de Recherche de l’Institut Universitaire en Santé Mentale de Québec, Université Laval, Québec, Canada
    Contribution
    SS, Acquisition of data, Analysis and interpretation of data, Drafting or revising the article
    Competing interests
    The authors declare that no competing interests exist.
  5. Igor Timofeev

    1. Department of Psychiatry and Neuroscience, Université Laval, Québec, Canada
    2. Centre de Recherche de l’Institut Universitaire en Santé Mentale de Québec, Université Laval, Québec, Canada
    Contribution
    IT, Acquisition of data, Analysis and interpretation of data, Drafting or revising the article
    Competing interests
    The authors declare that no competing interests exist.
  6. Sydney S Cash

    Department of Neurology, Massachusetts General Hospital and Harvard Medical School, Boston, United States
    Contribution
    SSC, Acquisition of data, Analysis and interpretation of data, Drafting or revising the article
    Competing interests
    The authors declare that no competing interests exist.
  7. Eric Halgren

    Departments of Radiology and Neurosciences, University of California, San Diego, La Jolla, CA, United States
    Contribution
    EH, Acquisition of data, Analysis and interpretation of data, Drafting or revising the article
    Competing interests
    The authors declare that no competing interests exist.
  8. Maxim Bazhenov

    Department of Medicine, University of California, San Diego, La Jolla, CA, United States
    Contribution
    MB, Conception and design, Analysis and interpretation of data, Drafting or revising the article
    Competing interests
    The authors declare that no competing interests exist.

Funding

Canadian Institutes of Health Research (MOP-136969)

  • Igor Timofeev

Canadian Institutes of Health Research (MOP-136967)

  • Igor Timofeev

Office of Naval Research (MURI: N000141310672)

  • Sydney S Cash
  • Eric Halgren
  • Maxim Bazhenov

National Institutes of Health (MH099645)

  • Eric Halgren
  • Maxim Bazhenov

National Institutes of Health (EB009282)

  • Eric Halgren
  • Maxim Bazhenov

National Sciences and Engineering Research Council of Canada (298475)

  • Igor Timofeev

The funders had no role in study design, data collection and interpretation, or the decision to submit the work for publication.

Acknowledgements

This work was supported by grants from ONR (MURI: N000141310672) and NIH (MH099645 and EB009282), Canadian Institutes of Health Research (grants MOP-136969, MOP-136967) and by National Sciences and Engineering Research Council of Canada (grant 298475).

Ethics

Human subjects: The data were collected part of other study where informed consents were obtained. This is primarily a modeling paper and the human sleep data are used to illustrate a point.

Animal experimentation: The data were collected part of other study where ethical guidelines were followed. This is primarily a modeling paper and the cat and mice sleep data are used to illustrate a point.

Reviewing Editor

  1. Ronald L Calabrese, Emory University, United States

Publication history

  1. Received: June 9, 2016
  2. Accepted: October 11, 2016
  3. Version of Record published: November 16, 2016 (version 1)

Copyright

© 2016, Krishnan et al.

This article is distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use and redistribution provided that the original author and source are credited.

Metrics

  • 1,745
    Page views
  • 429
    Downloads
  • 11
    Citations

Article citation count generated by polling the highest count across the following sources: Scopus, Crossref, PubMed Central.

Download links

A two-part list of links to download the article, or parts of the article, in various formats.

Downloads (link to download the article as PDF)

Download citations (links to download the citations from this article in formats compatible with various reference manager tools)

Open citations (links to open the citations from this article in various online reference manager services)

Further reading

    1. Neuroscience
    2. Structural Biology and Molecular Biophysics
    Dipak N Patil et al.
    Research Article
    1. Biochemistry and Chemical Biology
    2. Neuroscience
    Apurwa M Sharma et al.
    Research Advance