https://orcid.org/0000-0003-3194-3258
Figures
Abstract
Theta and gamma rhythms and their cross-frequency coupling play critical roles in perception, attention, learning, and memory. Available data suggest that forebrain acetylcholine (ACh) signaling promotes theta-gamma coupling, although the mechanism has not been identified. Recent evidence suggests that cholinergic signaling is both temporally and spatially constrained, in contrast to the traditional notion of slow, spatially homogeneous, and diffuse neuromodulation. Here, we find that spatially constrained cholinergic stimulation can generate theta-modulated gamma rhythms. Using biophysically-based excitatory-inhibitory (E-I) neural network models, we simulate the effects of ACh on neural excitability by varying the conductance of a muscarinic receptor-regulated K+ current. In E-I networks with local excitatory connectivity and global inhibitory connectivity, we demonstrate that theta-gamma-coupled firing patterns emerge in ACh modulated network regions. Stable gamma-modulated firing arises within regions with high ACh signaling, while theta or mixed theta-gamma activity occurs at the peripheries of these regions. High gamma activity also alternates between different high-ACh regions, at theta frequency. Our results are the first to indicate a causal role for spatially heterogenous ACh signaling in the emergence of localized theta-gamma rhythmicity. Our findings also provide novel insights into mechanisms by which ACh signaling supports the brain region-specific attentional processing of sensory information.
Author summary
Recent evidence indicating that ACh signaling is both transient and spatially circumscribed raises the question of how this feature impacts information processing in cortical networks. Here we demonstrate that spatially segregated ACh modulation of excitatory-inhibitory neural networks generates theta-modulated gamma rhythms, a hallmark of attention and information processing. Theta-gamma coupling arises naturally as neuronal activity traverses high-ACh regions of the network, and gamma activity alternates between distinct high-ACh sites at theta frequency. These findings provide novel insights into neurophysiological mechanisms for ACh-regulated theta-gamma coupling and uncoupling.
Citation: Yang Y, Gritton H, Sarter M, Aton SJ, Booth V, Zochowski M (2021) Theta-gamma coupling emerges from spatially heterogeneous cholinergic neuromodulation. PLoS Comput Biol 17(7): e1009235. https://doi.org/10.1371/journal.pcbi.1009235
Editor: Jonathan Rubin, University of Pittsburgh, UNITED STATES
Received: December 18, 2020; Accepted: July 1, 2021; Published: July 30, 2021
Copyright: © 2021 Yang et al. This is an open access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.
Data Availability: Data is held at https://github.com/YihaoYang/HeterACh-FreqCoupling.
Funding: This work was supported by National Institute on Drug Abuse R01DA045063 https://www.drugabuse.gov (MS), National Institute of Neurological Disorders and Stroke PO50NS091856 https://www.ninds.nih.gov (MS), National Institute of Neurological Disorders and Stroke R01 NS118440 https://www.ninds.nih.gov (SJA), National Institute of Mental Health DP2 MH 104119 https://www.nimh.nih.gov (SJA), National Science Foundation BCS-1749430 https://nsf.gov (VB) and by the National Institute of Biomedical Imaging and Bioengineering R01EB018297 https://www.nibib.nih.gov (MZ). The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.
Competing interests: The authors have declared that no competing interests exist.
Introduction
Acetylcholine (ACh) signaling in neocortex emanates from the basal forebrain (BF). Recent anatomical studies indicate that in contrast to more traditional views of the BF projection system as “diffusely” organized, afferent and efferent projections of the BF ACh system are highly topographically organized [1–4]. Prior notions of BF ACh activity as having relatively low temporal resolution and spatial selectivity, and characterized by predominantly extra-synaptic actions (i.e., “volume” transmission) [5–8], have also been refuted by more recent evidence indicating fast and spatially discrete ACh spread [8,9]. As a functional corollary of these developments, the prior conceptualization of ACh as acting diffusely and globally has been challenged by studies indicating event- or task trial-specific ACh signaling in specific neocortical regions [8,10]. Fig 1 shows data from multiple electrochemical recording sites in proximity, demonstrating asynchronous, spatially segregated neocortical ACh signaling.
Cholinergic signaling is largely asynchronous and can influence target circuitry in a temporally and spatially highly heterogenous manner. Data recorded in prelimbic cortex (A) are shown here. The four platinum (Pt) recording sites fabricated onto a ceramic backbone electrode are illustrated in B and the placement of these recording sites in prelimbic cortex are shown in A. C depicts the dimensions of an individual recording site. The upper and lower pairs of recording sites were separated by 100 μm. The data shown in D-G were recorded via an upper sensor (“sensor 1) and a lower sensor (“sensor 2”). The neurochemical recording scheme, shown in C, was previously described in detail. Amperometric measures were validated in terms of reflecting newly released acetylcholine (ACh) release [10,52–54]. D. Fixed-potential amperometry signals from a representative animal during a 23.5-hour recording. Scored sleep states are identified above the raw amperometry signals shown below and include rapid eye movement (REM), Slow-wave sleep (SWS) and waking (WAKE) periods. Transients are denoted by arrows with unique transients shown in the same color as the corresponding trace and common transients shown as double blue arrows. E. Cholinergic transient event rates as a function of sleep state. Absolute event rates across the two sites show a high degree of similarity between one another. REM sleep states show the highest event rates/minute while SWS states show the fewest. F, G. Event rasters contrasting the timing of transient events for the two sites shown in d and e. Data from the opposite channel is shown relative to the onset of transient events for sensor 2 (F) or sensor 1 (G). A total of 524 events were detected from sensor 2 and 382 events were detected from sensor 1 in the 23.5-hour recording. Note that while the highest concentration of activity is coincident across the two sensors, only a fraction of each channels events occurred in close proximity (±2.5 s) to one another (sensor 1: 40.3%; sensor 2: 29.9%).
Here, we investigate in silico the effects of temporally static but spatially circumscribed, heterogeneous activation of muscarinic ACh receptors on the activity patterns of excitatory-inhibitory (E-I) neural networks. Our simulation results indicate that localized theta (∼ 5 − 10Hz) and gamma (∼ 30 − 100Hz) band activity rhythms emerge in response to spatially segregated ACh modulation of neural excitability. Here, the modeled spatial cholinergic distributions are meant to represent a short snapshot from the evidence of spatially circumscribed ACh signaling in recording studies in rodents (Fig 1D), where discrete locations of high levels of cholinergic signaling were observed adjacent to locations with low levels of cholinergic activity. We analyzed the emerging neuronal activity patterns in the presence of stationary high levels of cholinergic signaling in a single versus in multiple locations of the network. Localized, gamma band activity rhythms emerged in cells undergoing high levels of cholinergic stimulation. Moreover, for multiple high-ACh ‘hotspots’, these gamma oscillations appeared only within the currently active network regions, resulting in their modulation at theta frequency.
Our results postulate that theta-gamma coupling is an emergent property of spatially segregated ACh modulation of neural response properties. We further identified mechanisms underlying the dependence of theta-gamma coupled activity on the spatial distribution of simulated ACh neuromodulation. In particular, gamma-band activity was supported in high-ACh regions via the pyramidal-interneuron gamma (PING) mechanism [11], where inhibitory interneurons strongly modulate and synchronize activity of pyramidal cells [11,12]. Theta band modulation of gamma activity within or between high-ACh regions was associated with spike frequency adaptation, linked to effects of muscarinic receptor activation on M-type K+ currents [13]. These mechanisms led to intrinsically tight coupling between gamma and theta band activity where the degree of theta-gamma coupling correlated with proximity to high-ACh regions. Additionally, we investigated the consequences of spatially heterogeneous ACh modulation on the attentional processing of external (sensory) stimuli.
Theta-gamma coupled activity in cortical and hippocampal areas is thought to be a hallmark of attentive cognitive processing [14] and multiple experimental studies have shown that ACh signaling promotes theta-gamma coupling in these circuits [15,16] (see Discussion). Our modeling results propose that this cognitively significant firing pattern is directly caused by spatially heterogeneous modulation of neural properties due to spatially circumscribed release of ACh.
Results
New experimental data indicates that cholinergic signaling is largely asynchronous and can influence target circuitry in a temporally and spatially highly heterogenous manner. An example of such evidence is presented on Fig 1. In this experiment four platinum recording sites were fabricated onto a ceramic backbone electrode where the upper and lower pairs of recording sites were separated by 100 μm (Fig 1A). Fig 1D depicts sample measurements (in terms of currents) depicting localized temporal changes of ACh concentrations. Further analysis of these transients indicates that, while their overall number follows standard notion that the highest ACh release is during REM and the lowest happens during the SWS (Fig 1E), this release is localized and highly asynchronous (Fig 1F and 1G).
Here we elucidated, using in silico modeling, the dynamical changes in neuronal activity patterns stemming from such a localized and asynchronous mode of ACh signaling. Namely, we investigated how spatially localized regions of ACh neuromodulation generated network-wide oscillatory activity in the gamma and theta frequency bands, in two-dimensional E-I networks. Using Hodgkin-Huxley type model neurons, ACh effects on the slow, adapting M-type K+ current were simulated by varying its maximal conductance, , in individual cells across the network. Through the action of muscarinic receptors, ACh blocks the M-type K+ current, thus low values of
corresponded to high ACh modulation. Spatially heterogeneous ACh modulation was achieved by constructing spatial mappings of
values for individual cells in the network. The
map mimicked the post-synaptic effects of spatially localized, asynchronous transients of ACh release in a short time window (as suggested by results shown in Fig 1). To this end, we postulate that relevant network dynamics occur on two separate temporal scales: 1) fast dynamics on a timescale of milliseconds associated with neuronal firing, and 2) slow timescales on the order of 5–10 s associated with localized ACh release and subsequent degradation or uptake. Here, for simplicity, we investigated fast scale neuronal dynamics in the presence of fixed spatial distributions of ACh. In the
maps, each unit square corresponded to a single cell in the model network (i.e. the unit length, being the side of the square, corresponded to the minimal distance between modeled cells), and the dimensions of
modulated regions we consider encompass, on average, tens of neurons.
With network connectivity fixed in a local excitation/global inhibition topology, we show that different spatial mappings, which affect the excitability of both E and I cells, result in gamma and theta band rhythmic activity. While there is evidence that the M-current is present not only in excitatory pyramidal cells, but also in various classes of inhibitory interneurons, namely somatostatin positive (SST+) neurons [17], we obtained qualitatively similar results when neuromodulatory effects of the M-current were included only in excitatory neurons (S1 Fig). Consistent with previous work [11,12], continuously active E and I cells in regions of low
values, corresponding to high ACh modulation, show gamma band oscillatory activity. Here, theta band modulation of gamma activity is an emergent property of the network, generated as firing activity traverses within or between spatially segregated low
regions. Below, we first illustrate this novel mechanism for theta-gamma coupling with a randomly generated, spatial
mapping and then analyze the mechanism in more detail for simple spatial
mappings.
Coexistence of theta/gamma rhythms is caused by spatially segregated 
distributions
To illustrate the emergence of theta-gamma coupled firing activity due to a spatially complex ACh landscape, we generated a spatially heterogeneous mapping by randomly assigning n = 9 center positions of low
regions or ’hotspots’ of radius r = 4.2 (units in minimal distance between model cells, Fig 2A).
A, Randomly generated spatial map of values on the 20×20 E cell lattice. B, Spike raster plot illustrating portion of E cell (cells 1−400) and I cell (cells 401−500) firing patterns. The pixel color indicates cell
value (same color scale as in A). The shaded area indicates the time range of snapshots in F. C, Average E cell firing frequencies plotted as cell position on the E cell lattice showing highest firing within the regions or ’hotspots’ of low
. D, Dominant rhythmic activity of individual E cells (dark blue = none, light blue = theta band, green = gamma band and yellow = mixed, both gamma and theta) displayed on E cell lattice illustrating that E cells with lower
values showed stronger and more stable gamma activity while cells with moderate
values exhibited stronger theta rhythm. E, Power spectrum analysis of E cell network firing. High power exhibited in both theta (5-12Hz) and gamma (40-60Hz) rhythm bands. F, Distribution of
values for cells showing specific rhythmic activity illustrating correlation of activity type with
value (same color code as in D). G, Snapshots of E cell firing rates from 4500–4850 ms during the simulation shown in B (shaded area). Cells inside the orange contour lines have
values less than 0.6 mS/cm2.
Cells within the ’hotspots’ exhibited gamma band firing activity that was intermittently interrupted as activity moved to other regions with low
values, with network activity cycling periodically through the hotspots with theta band frequency. As illustrated in the raster plot in Fig 2B, this activity pattern resulted in subpopulations of E cells showing theta-band modulated gamma activity, and high power in both the theta and gamma frequency ranges in overall network activity (Fig 2E). Individual E cell’s firing frequencies were directly correlated with their
values (Fig 2C, 2D and 2F). Namely, E cells near the centers of low
regions (i.e. corresponding to highest concentrations of ACh) had the highest firing frequencies (Fig 2C) with higher power gamma band activity (Fig 2D and 2F), while E cells with moderate
values fired in theta band ranges. Cells outside the
hotspots showed little activity.
For comparison, in S2 Fig we show sample raster plots and analysis of network wide activity patterns and cell spiking frequencies for homogeneous spatial maps. In these simulations, we have kept
at the same value for all the neurons and varied the value between 0 and 1.5 mS/cm2. We observed emergence of networkwide gamma band oscillations, and cellular spiking in that frequency, for low values of
, while systematic theta band oscillations were not present across the entire
range.
Characterization of gamma and theta band activity caused by a single peak gKs distribution
To understand how spatially heterogeneous distributions generate coupled theta-gamma activity, we next analyzed simple spatial gKs distributions. For a single low gKs region (Fig 3), firing activity was restricted to within the gKs hotspot by the higher neuronal excitability elicited by low gKs values and by the global inhibition provided to the rest of the network by the firing cells within the hotspot. Gamma range firing frequency of E cells within the gKs hotspot was generated by gating of E cell firing by local I cells through the PING mechanism [18–20]. As the radius r of the low gKs hotspot increased (without changing its minimum value; S3 Fig), the number of neurons exhibiting gamma slightly increased and theta band activity remained low (Fig 3B). As the hotspot radius increased past r = 5.5, however, the number of neurons primarily exhibiting gamma band activity started to decrease while the number of neurons firing at mixed theta and gamma ranges started to increase. While gamma neurons were located within the center of the gKs hotspot, theta neurons or neurons showing mixed theta-gamma band activity were located towards the hotspot’s outer edges (Fig 3C and 3D). Thus, theta band activity was generated due to firing activity spatially drifting within the gKs hotspot (Fig 3E).
A, Example single peak gKs mapping projected on the 20×20 E cell lattice where r indicates the radius of the approximately circular region (hotspot) with low gKs values. B, The number of neurons primarily exhibiting gamma (green curve), theta (blue curve) and mixed rhythms (yellow curve) as a function of the radius of the gKs hotspot. C, D, Dominant rhythmic activity of individual E cells (dark blue = none, light blue = theta band, green = gamma band and yellow = mixed, both gamma and theta) plotted at cell position on the E cell lattice for gKs hotspot radius of r = 5.5 and r = 6.1, respectively. E, Snapshots of E cell firing activities from 4500–4750 ms during the simulation shown in D (gKs hotspot radius r = 6.1). Cells inside the orange contour lines have values less than 0.6 mS/cm2.
This result echoed our previous findings [13] in similar 2D E-I networks with spatially uniform gKs values. In those networks, when gKs values were low, firing activity was spatially localized in a subregion of the network with surrounding cells showing intermittent firing, leading to drifting of an activity bump. When gKs values were higher (gKs>0.2 mS/cm2), the bump exhibited translational motion across the whole network promoted by spike frequency adaptation mediated by the gKs regulated M-type K+ current. Here, similarly to those results, the movement of firing activity within gKs hotspots with larger radius resulted in theta band rhythmicity.
We further investigated how the emergence of theta/gamma oscillations due to a single gKs hotspot depended on the constant current applied to individual cells (which corresponds to external input) and the level of simulated ACh concentration (i.e. lower bound of gKs value within the hotspot, S4 Fig). The emergence of theta/gamma oscillations was very robust across these parameters. We generally observed that the oscillations emerged within the hotspot for lower minimal values of gKs, except when neuronal external current was high and the oscillations were not confined to a hotspot but spread throughout the network. Another exception was for the lowest value of gKs, when firing activity was stationary within the hotspot and only gamma oscillations were observed.
We additionally investigated how the theta and gamma oscillatory frequency changes with the size of the hot spot (S5 Fig). We observed that frequencies in both bands tend to decrease with the increase of hotspot radius. This is most likely due to the fact that as the hotspot increases in size the local feedback inhibition increases reducing the excitation individual cells receive.
Finally, we investigated how network topology (i.e. rewiring of the excitatory connectivity) affected the emergence of theta/gamma oscillations (S6 Fig). We observed that increasingly random E-E connectivity abolished theta band oscillations, while the presence of gamma band activity remained largely unchanged. This result suggests that the generation of theta-gamma coupling requires that synaptic excitation which is generally more local than synaptic inhibition in the network. We also observed a decrease in gamma frequency as a function of network rewiring. This is due to the fact that increasingly random connectivity promotes zero phase synchrony causing EPSPs mediated by spiking of presynaptic neurons to partially fall within the spike and/or refractory time of their postsynaptic targets, reducing their excitatory effect and the overall level of excitation in the network (S6B Fig).
Characterization of theta/gamma band oscillations and their coupling emerging from a double peak gKs distribution
To investigate the emergence of rhythmic network activity when two spatially adjacent locations experienced cholinergic modulation at the same time, we considered the presence of two gKs hotspots in the network (Fig 4). Fig 4A shows an example of this kind of gKs mapping, in which d represents the distance between the two gKs hotspot centers while r denotes their radius (units in minimal distance between model cells). Here, cells located in the center of the hotspots exhibited theta-modulated gamma rhythms (mixed) while those on the peripheries of the hotspots showed primarily theta activity (Fig 4B). This occurred because spiking activity alternated between the hotspots at theta frequency. When a given hotspot was active it predominantly exhibited gamma band oscillation. This resulted in strong theta-gamma coupling as the gamma band oscillations appeared at a given site with theta band frequency (Fig 4B).
A, Example of double peaked spatial mapping of gKs values for corresponding neurons on the 20×20 E cell lattice for hotspot radius r = 4.6 and distance between hotspot centers d = 8. B, Dominant rhythmic activity of individual E cells (dark blue = none, light blue = theta band, green = gamma band and yellow = mixed, both gamma and theta) plotted at cell position on the E cell lattice due to double peaked gKs spatial distribution shown in A. C, D and E, Numbers of neurons exhibiting dominant theta (C), dominant gamma (D), or theta-gamma coupled activity (E) as a function of radius of gKs hotspots r and the distance between hotspot centers d for a double peaked gKs spatial distribution. F, G and H, Spike raster plots (top panels) and histograms for the dominant rhythmic activity exhibited by neurons based on their individual gKs values (bottom panels) for networks with double peaked gKs distributions with r and d values indicated by labels in C, D, E. In the raster plots, E cells are numbered 1 to 400, and I cells are numbered 401 to 500. Color indicates gKs values of cells with the scale in A. I, Snapshots of E cell firing rates from 4500–4800 ms during the simulation in H (shaded area in H; r and d values indicated by label H in C). Cells inside the orange contour lines have values less than 0.6 mS/cm2.
To analyze the properties of theta-gamma coupling as a function of the size and relative position of the hotspots, we varied the parameters r and d of the gKs spatial mapping and identified the numbers of cells showing primarily gamma or theta band activity, or theta-gamma coupled firing (Fig 4C, 4D and 4E). We observed that the numbers of neurons with theta, gamma and mixed rhythms increased monotonically with larger r (Fig 4C and 4E). At the same time, the highest number of cells representing all three types of dynamics (i.e. theta, gamma and mixed) happened around d = 6 with gamma and theta cell numbers declining faster for larger d than the mixed rhythm population. (Fig 4D).
To better understand these effects, we closely analyzed network dynamics for chosen sets of parameters. The parametric locations of the raster plots and cell statistics for different rhythms depicted on Fig 4F, 4G and 4H are indicated in Fig 4C, 4D and 4E. When r = 4.2 and d = 4 (Fig 4F), two gKs hotspots were small enough and sufficiently close to each other that cells in both spots fired simultaneously. Active neurons mainly fired with gamma rhythmicity except a few neurons with intermediate gKs values carried a theta rhythm. No neurons exhibited both rhythms in this case.
In contrast, for larger r and d values, cells in the two spots started to show alternating firing patterns (Fig 4G and even more so Fig 4H and 4I) with mixed theta-gamma rhythmicity. As described above, firing activity alternated between the two spots such that cells in each spot fired intermittently at gamma frequency with the alternation in activity between the spots occurring at theta frequencies. This alternation in firing between the hotspots occurred due to competition between the local excitation among cells within the hotspot, global inhibition generated by cells in the other hotspot and the magnitude of spike frequency adaptation (SFA) mediated by the activation of the M-current. Firing within the active hotspot mediated inhibitory signaling received by E cells in the silent hotspot. Due to their high excitability with low gKs values, feedback excitation with neighboring cells and small effects of SFA due to previously low activation, the silent hotspot E cells could start firing and inhibit the active hotspot. Subsequently as the SFA accumulated in the activated region the other hotspot can get triggered and take over. For larger radius values, theta-gamma activity dominated the network as activity moved consistently between the hotspots and only a few cells at the centers of the hotspots fired primarily at gamma frequency (Fig 4H).
All the simulation results presented here were performed for an M-current activation time constant of τz = 75ms following [21]. To better understand how the M-current timescale affects frequency of theta and gamma band oscillations, we performed additional simulations scanning different values of τz∈[25ms, 125ms] (S7 Fig). As expected, we observe a significant slowdown of theta band frequency with increased τz, while no systematic changes were detected in gamma band oscillations.
We also analyzed how theta and gamma frequency changed as a function of gKs hotspot size and the distance between hotspot centers (S8 Fig) and as a function of connectivity of the excitatory subnetwork (S9 Fig). We observed that frequencies of both theta and gamma oscillations decreased as the hotspots’ radius increased. However, the frequencies of the two oscillations showed opposite trends with the increase of hotspot center distance–frequency of the theta band decreased while the gamma frequency generally increased. We attribute this effect to the fact that as the two hotspots are positioned farther apart it takes longer for them to switch their activations. When excitatory connectivity was made increasingly random (S9 Fig), activity switching between the two hotspots stopped and firing became synchronous with cells in both hotspots exhibiting gamma band oscillations.
Finally, we investigated whether similar dynamical switching between the hotspots would occur if, instead of gKs hotspots, increased excitability was driven with additional external current in two hotspots (; see Materials and Methods), with a homogeneous gKs map across the network (S10 Fig). We observed that only for very narrow range of gKs values is it possible to get localized, selective activity switching between the two stimulation hotspots with theta frequency (namely for gKs ~0.2). For other values of gKs, activity switching does not take place, or is not specific to the neurons within the hotspots (i.e other neurons outside the hotspots become activated; S10 Fig–black colored spikes in addition to red spikes). This non-specific activation obliterated theta band oscillations. Thus, we interpret these results that while emergence of theta/gamma oscillations may be possible in response to heterogeneous patterns of external drive, they occur for a narrow band of network gKs modulation and they are significantly less robust.
Variability of theta/gamma band oscillations with spatially random gKs distributions
The analysis in the preceding sections illustrates that, in networks with local excitation and global inhibition, theta-modulated gamma band firing occurred in cells on the outer edges of individual gKs hotspots or, if the gKs hotspot was large compared to the excitation range, within the hotspot itself. With multiple, spatially separated gKs hotspots, all cells within gKs hotspots can exhibit theta-modulated gamma band firing as competition between local excitation and global inhibition within and between hotspot cells causes alternation of firing episodes. These results indicate that both sizes as well as relative positions of the hotspots matter. In the more biologically realistic scenario of spatially random gKs distributions consisting of multiple hotspots, these same mechanisms contribute to generating theta-gamma coupled firing, however we found that resulting strengths of theta and gamma band activity were highly variable, depending on the specific realization of the gKs spatial map (i.e. positions of the hotspots).
As an example, Fig 5 shows different realizations of gKs spatial mappings with 6 hotspots of radius r = 5.4 (Fig 5A and 5B, star in S5 Fig) and r = 2.8 (Fig 5C and 5D, square in S8 Fig). Despite the same hotspot properties, the number and size of effective low gKs regions was varied as hotspots could overlap and coalesce if their centers were next to each other. Network dynamics were likewise highly variable with some mappings showing clear coexistence of theta/gamma rhythms (Fig 5A and 5C) while in others gamma (Fig 5B) or theta (Fig 5D) power dominated. This is due to the fact that for the gKs mappings with coexisting theta/gamma rhythms, individual hotspots coalesced into larger but spatially distinct low gKs regions (Fig 5A and 5C; note that networks have periodic boundary conditions so hotspots may wrap around the lattice). In these cases, theta-gamma coupling was generated similarly as in the double peaked gKs spatial distribution. On the other hand, gamma power dominated when the individual gKs hotspots coalesced into an effective single low gKs region (wrapped around the corners of the network, Fig 5B). The strongest theta power was generated when hotspot centers were more evenly dispersed across the network and activity traveled successively along the low gKs regions (Fig 5D). This observed variability in theta/gamma power provides a possible insight into experimental inter-animal variability, as the relative locations of ACh release sites could be highly individualized within the experimental subjects.
The gKs spatial mappings (color plots), network spike raster plots (right panels, top) and network frequency power spectrum (right panels, bottom) for two pairs (A, B and C, D) of gKs map realizations with the same sets of parameters. A, B The results of 2 random gKs mappings with 6 hotspots of radius r = 5.4 (marked by ’star’ in S5 Fig). C, D The results of 2 random gKs mappings with 6 hotspots of radius r = 2.8 (marked by ’square’ in S5 Fig).
To gain a better understanding of the variability in observed network rhythmicity, we simulated effects of random gKs mappings having different numbers of hotspots of different sizes (S11 Fig) and measured mean theta and gamma power (and their ratio) and their relative standard error (RSE) across simulation runs with different instantiations of the mappings. We additionally measured peak theta and gamma frequency for each of these hotspot distributions (S12 Fig). As shown in the example above, the relative power amplitudes varied significantly for maps consisting of the same number of hotspots having the same sizes, reporting relatively high RSE, except for the situation when only one hotspot was present.
Vicinity to high ACh region modulates strength of theta/gamma coupling
It has been shown experimentally that ACh release promotes theta-gamma coupling, and furthermore, this coupling is specifically mediated via M1 muscarinic receptors [8,15]. In our simulations, cells showing the strongest theta-gamma coupling were located within gKs hotspots, when a single large hotspot (Fig 3), or more than one hotspot were present (Fig 4). To investigate how the strength of theta-gamma coupling varied with location relative to the position of the hotspots, we constructed local field potential (LFP) signals at different distances from a gKs hotspot (Fig 6). Here, specifically, we used a double peaked gKs spatial mapping (Fig 6A) since this network showed the most robust theta/gamma coupling. The locations at which LFP was calculated is bounded by “C” and “D”, corresponding to the respective panels showing rasterplots from those locations (Fig 6C and 6D, red dots correspond to spikes of the neurons used for LFP calculation at these respective locations). The LFP signal was computed from the spike trains of the 12 E cells closest to the marked location (see Materials and Methods). For completeness, S13 Fig shows the results for an LFP signal generated directly from the voltage traces of the selected cells, with results remaining qualitatively the same. The locations were chosen to be progressively further away from the center of one of the gKs hotspots (x-axis in Fig 6B). Separately filtering the LFP signal in the theta and gamma bands (Fig 6E corresponding to location C and Fig 6F corresponding to location D), we calculated the modulation index (MI) to quantify the coupling strength between the two signals as described in the Materials and Methods section. MI values as a function of distance from the gKs hotspot center decreased (Fig 6B), showing that theta-gamma coupling strength decreased with distance from the gKs hotspot.
A local field potential (LFP) signal was constructed from spike trains of E cells at different distances from the center of a gKs hotspot. A, Double peaked gKs spatial mapping with locations marked to calculate LFPs. B, Modulation Index (MI) between gamma and theta filtered LFP traces as a function of the distance from the center of the gKs hotspot (as indicated in A.) C and D, Example raster plot of the network-wide spiking near the hotspot (C, as marked on A) and away from the hotspot (D, as marked on A). E cells are numbered 1 to 400, and I cells are numbered 401 to 500. Cells used to compute the LFP signal are marked in red. E and F, Examples of theta (blue curve) and gamma (orange curve) filtered LFP signals computed at locations near (C) and far (D) from the hotspot, marked on A, B.
Gamma frequency spiking is evident primarily at the locations which have low gKs (high ACh; In Fig 6C, blue dots correspond to spikes of neurons located within the gKs hotspot as color coded on Fig 6A). This is due to the fact that only these locations can generate reliable PING dynamics. Away from these regions, even if activity of the network briefly traverses a given location, the gamma oscillations will be greatly reduced or not present at all (due to lack of excitability mediated by low gKs). Thus, the diminished theta-gamma coupling away from the gKs hotspot is due to reduction of local gamma oscillations. Examples of the two filtered LFP signals at different network locations are shown in Fig 6D and 6E.
Effect of gKs modulation on network response to external stimuli
Behavioral attentional tasks, particularly visual attention tasks, report differences in responses to sensory-mediated stimuli depending on whether such stimuli undergo attentional processing [22,23]. Localized cholinergic signaling was shown to be necessary and sufficient for the attentional elevation of the processing of stimuli [24,25].
To investigate how spatially heterogeneous gKs modulation (simulating ACh mediated attentional drive) affects responses to external (sensory) stimuli, we measured relative changes in firing frequency of a subset of excitatory cells targeted by an external excitatory drive, , when the targeted E cells were inside or outside of a gKs hotspot (i.e. within or outside attentional drive, respectively). We compared three situations corresponding to three behavioral conditions: 1) targeted E cells are inside the gKs hotspot corresponding to presentation of a sensory stimulus which is attended to; 2) targeted E cells are outside the gKs hotspot corresponding to presentation of a sensory stimulus but attention is directed elsewhere; and finally 3) the subset of E cells is targeted by the external drive but there is no gKs modulation in the network (spatially homogeneous gKs at the default value) corresponding to presentation of a sensory stimulus in the absence of attention. These conditions were simulated for single (Fig 7; left column) and double (Fig 7; right column) peaked gKs spatial mappings.
A, B, C, D, Average E cell firing frequencies plotted as cell position on the E cell lattice showing varied firing responses to an external excitatory stimulus targeting a subset of E cells. Dashed circles indicate location and radius of gKs hotspots and arrows indicate location of targeted cells. Targeted cells were either inside (A, B) or outside (C, D) of the gKs hotspot. E, The change in average firing responses of the targeted E cells to the external stimulus relative to their responses in the absence of gKs modulation (homogeneous gKs at default value). Yellow (blue) bars correspond to when targeted E cells are inside (outside) the gKs hotspot. The four bars from left to right correspond to cases shown in A, C, B and D, respectively.
We observed that for the spatial mapping with single peak gKs, the firing response of the targeted neurons (Fig 7A; location marked by arrow) was significantly higher when they were located in the gKs hotspot (condition 1), relative to their response to the stimulus in the absence of gKs modulation (condition 3) (Fig 7A). When targeted cells were outside the gKs hotspot (Fig 7C location marked by arrow), their firing response to the stimulus (condition 2) was significantly suppressed compared to their response in the absence of gKs modulation (condition 3) (Fig 7C and 7E). This suppression was due to the global inhibition induced by the cells firing in the gKs hotspot that attenuated the response to the external stimulus.
For the spatial mapping having two gKs peaks, responses of targeted cells were generally the same as for the single hotspot case whether the targeted cells were located inside or outside of a gKs hotspot (Fig 7B and 7D). We additionally observed that when the targeted cells were in one of the gKs hotspots, their firing dominated the network, shutting down activation at the other hotspot. This led to abolition of theta-gamma coupling as the neurons in the targeted hotspot fired at gamma frequency.
Discussion
Cholinergic signaling is necessary and sufficient for the detection of cues in attentional contexts. Moreover, cholinergic signaling influences the degree of forebrain desynchronicity across circadian stages. Until recently it has been thought that cholinergic signaling occurs at a relatively low temporal resolution but also with highly limited spatial selectivity. In contrast, recent results indicated that ACh release is more localized and asynchronous within activated brain modalities (for example, Fig 1). Such evidence for spatially heterogenous ACh signaling has remained rare, in part because until the advent of biosensors and, more recently, the GRABACh sensor, prior methods available for monitoring ACh release did not allow for measurements at a relevant spatial resolution. However, using the ACh3.0 GRAB sensor, spatially heterogenous cholinergic signaling was recently shown in the somatosensory cortex of walking mice. Cholinergic hotspots with a diameter of about 40 μm were selectively activated by runs and appeared to be surrounded by inactive areas of similar diameters (see Fig 3O in [26]). Using a prior version of this sensor to monitor ACh release in hippocampal slices, the size of hotspots was concluded to be even smaller, about 16.5 μm in diameter [27]. The aim of this modeling study was to understand how a highly heterogenous distribution of ACh affects rhythmic network firing activity.
We have previously shown [12], using biophysical computational modeling, that a network of randomly coupled excitatory and inhibitory neurons can generate transient gamma oscillatory activity in response to simulated spatially global but temporally brief pulses of ACh. This effect was mediated by blockage of the M-type K+ current (i.e. gKs conductance). Depending on network connectivity, gamma activity decayed with the simulated gKs transient modulation or was sustained in the network after the gKs transient completely dissipated.
In contrast, the present study demonstrates that in a network with local excitation/global inhibition connectivity, spatially heterogenous ACh modulation, modeled by varying the M-current differently in different cells, leads to the emergence of spatially localized theta and gamma band activity rhythms.
The appearance of the individual rhythms and theta-modulated gamma oscillations strongly depended on the specific features of the spatial distribution of the maximal conductance of the M-current, gKs (i.e. the number of gKs hotspots and their radius). Furthermore, we identified two basic mechanisms mediating strong theta-gamma coupling. First, if a single gKs hotpot is present and its size is larger spatially than the scope of the local excitation (i.e. range of E-E connectivity), the neuronal activation moves throughout the population encompassed by the hotspot. Second, when two or more gKs hotspots are present, the neuronal activity alternates between the hotspots. In both cases, the activated population exhibits gamma band activity. The emerging gamma rhythm is mediated via the PING mechanism, whereas the emerging theta rhythm is caused by activity traversing gKs hotspots, mediated by SFA.
Together, these simulations revealed that spatially heterogenous release of ACh leads to the manifestation of theta and gamma oscillations and their coupling. For spatially homogeneous gKs distributions, when gKs is low there is effectively no spike frequency adaptation in the network and thus gamma oscillations may be supported indefinitely (S2 Fig). On the other hand, homogeneous gKs distributions at high values lead to random neuronal spike-frequency adaptation patterns, without consistent formation of theta band oscillations. However, spatially localized and constrained regions of low gKs allow for theta band modulation of the activity within these regions with gamma oscillations emerging during phases of firing activity. Furthermore, the specific pattern and frequencies of observed oscillatory activity is highly dependent on the specific distribution of ACh release sites and neuronal network structure. For example, the specific frequencies observed in the theta and gamma range depend on network wiring topology and the cellular signaling properties (i.e. excitatory cell input or inhibitory time constant; S14 Fig, supplemental material). These in silico network dynamics reproduce important in vivo results, as discussed below.
Network connectivity plays a role in supporting the emergence of theta-band activity. The local excitation/global inhibition connectivity structure allows for firing activity in discrete and localized network sites, with the gKs distribution controlling the location and sequence of activated sites. We have shown, that similar results are obtained when inhibitory neurons had sparser random connectivity to both excitatory and inhibitory targets (S15 Fig, supplemental material). In general, as long as inhibition has a larger footprint than excitation in the network, qualitatively similar dynamics will be obtained. Such connectivity is thought to be found on the meso-scale in the connectivity of functional modules [28]. In addition we believe that a larger inhibitory footprint could be potentially obtained via long range excitatory connections that specifically target locally connected inhibitory cells.
Brain rhythms in the theta (∼ 5 − 12Hz) and gamma (∼ 30 − 100Hz) frequency bands have been shown to critically contribute to essential cognitive processes in many brain regions, including the cortex and hippocampus [19,20,29–32]. Recently, the coupling of these rhythms, such that the amplitude of gamma band activity is modulated by the phase of theta band activity, has been identified as a key component of the neuronal local field potential (LFP) or electrocorticogram (ECoG) observed during perceptual, attentional and other cognitive processes [22,33–37]. Such cross-frequency coupling has been proposed to signify the cooperation of diverse circuits and neuronal populations in order to integrate multiple cognitive operations and force the execution of stimuli-bound cognitive or behavioral outputs [38,39]. For example, in rats performing a signal-detection based attention task, successful performance of the task was shown to be dependent on the appearance of theta-modulated gamma band activity in the LFP measured in the prefrontal cortex (PFC) [15]. Generation of coupled theta-gamma activity was associated with fast (or transient), cue-bound cholinergic signaling in the PFC. Importantly, disruption of post-synaptic ACh signaling, by blocking muscarinic M1 acetylcholine receptors (mAChRs), attenuated gamma synchronicity, disrupted theta-gamma coupling and caused detection failures [15]. Studies in hippocampal and entorhinal cortical circuits have similarly shown that ACh signaling promotes theta-gamma coupling [16]. Our simulation results suggest that such muscarinic-dependent theta and gamma activity may be generated by spatially heterogeneous modulation of neural properties due to spatially circumscribed release of ACh.
The above cited results concentrate on theta and gamma oscillatory components found in the LFP signal, rather than directly in neural spiking activity. There are a limited number of studies suggesting that local increases in ACh signaling can lead to increases in gamma periodic multiunit spiking in the neocortex [40] and theta periodic spiking in the hippocampus [41]. However, to our knowledge the relationship of ACh and theta-gamma LFP coupling to spiking of neurons has not been studied in vivo. Clearly, there is a substantial need for such research.
Previous computational models for theta-gamma coupling activity in E-I networks relied on the presence of two populations of inhibitory cells that gate the firing of E cells with different time scales [42–46]. Namely, a fast I population generates E cell firing at gamma frequencies through the PING mechanism, while a slower I population gates the gamma oscillatory firing at theta frequencies. In the novel mechanism we describe here, obtaining theta band rhythmicity is mediated by activation of muscarinic receptors that modulate activity on timescales corresponding to theta band oscillation via SFA. Here, both excitatory and inhibitory cell populations are endowed with the M-current, however similar results are observed when muscarinic receptors are expressed (within the model) on excitatory pyramidal cells only.
In our networks, the strength of theta-gamma coupling varied with distance from the center of a gKs hotspot. In LFP measurements within the hotspot, gamma oscillations were tightly aligned with the peaks of the theta rhythm. This functional coupling decreased for locations away from the hotspot, leaving the theta and gamma band activity largely uncoupled. Although the available data from in vivo recordings already support the view that theta-gamma coupling is caused by, and occurs in the region of, elevated cholinergic signaling and muscarinic M1 receptor stimulation [15], prior electrochemical recording techniques have not achieved levels of spatial resolution that would allow the characterization of hotspots and their "colder" boundaries. More recent G-protein coupled ACh sensors appear capable of differentiating levels of cholinergic signaling on a scale of tens of micrometers [27]. However, measuring LFPs simultaneously with fluorescence imaging remains challenging and, for the observation of theta-gamma coupling, may require simultaneous recordings during task performance [15]. Progress in the development of in vivo recording techniques may soon allow a direct test of the neurobiological validity of our findings.
We observed that spatially localized gKs hotspots acted to gate responses to external excitatory input to the network. External stimuli applied to E cells within a hotspot generated a dramatic increase in spiking frequency of stimulated neurons and, to a smaller degree, around the network, as compared to the response of the stimulus in the absence of gKs modulation. Conversely, when external stimuli were applied to cells located outside the gKs hotspot, stimulated cells showed a significantly reduced response compared to their response in the absence of gKs modulation. These results mirror electrophysiological findings from primate visual cortex, where attention localized to neurons’ retinotopic field augment firing rate responses [47] as well as coherence of firing with theta and gamma oscillations [22,37]. This effect provides a mechanism not only to aid in detection and discrimination of sensory stimuli in attended locations [39], but also may aid in selective communication of responses to attended stimuli between distant cortical areas [22,37].
Our present results suggest that cholinergic activity influences target circuitry in a highly spatially heterogenous manner which influences the locations of cellular rhythmic activity. Results presented here consider “time-frozen” gKs (ACh) spatial distributions while in reality, they will be continually changing. The time scale of ACh diffusion/uptake can be inferred from Fig 1C to be somewhere between 5-15s, which is clearly long enough to establish stable theta/gamma band rhythms by the presented mechanism. Further development of our model by integrating the impact of such temporal dynamics, may have interesting implications for understanding the neuronal mechanisms contributing to, for example, the ‘Biased Competition’ model of attention, that is, the mechanisms that allow behaviorally significant stimuli to undergo feature extraction, while the processing of competing stimuli, even if placed in the same receptive field, is suppressed [48,49].
With our prior neurobiological and computational findings about the functional significance of fast, transient cholinergic signals, conceptualizations about the regulation and function of cholinergic neurons have advanced from traditional views about monolithic neuromodulator actions to temporally and spatially differentiated signaling across, for example, cortical layers and microcolumns. Together with our present results, this new paradigm for cholinergic signaling suggests novel underlying mechanisms for how cholinergic activity can rapidly re-direct information flow in target circuitry and thus play an essential role for maintaining behavioral and cognitive flexibility [8,50].
Materials and methods
Experimental recordings
The evidence depicted in Fig 1 was adopted from experiments which measured fast, transient cholinergic signals (“transients”) across circadian cycles in the cortex and hippocampus [51]. Data recorded in prelimbic cortex are shown here. The four platinum (Pt) recording sites were fabricated onto a ceramic backbone electrode where the upper and lower pairs of recording sites were separated by 100 μm. The data shown were recorded via an upper sensor (“sensor 1) and a lower sensor (“sensor 2”). The neurochemical recording scheme, was previously described in detail, while amperometric measures were validated in terms of reflecting newly released acetylcholine (ACh) [10,52–54]. Briefly, newly released ACh is hydrolyzed by endogenous acetylcholinesterase (AChE), and the resulting choline is oxidized by choline oxidase immobilized onto the Pt electrodes. The resulting hydrogen peroxide is oxidized electrochemically and current yields are recorded amperometrically. m-phenylenediamine (mPD) was electropolymerized onto the electrodes to enhance the selectivity for detecting analyte relative to currents produced by potential electroactive interference. Electrodes were calibrated prior to implantation into the brain. Electrochemistry data were collected via the FAST-16 recording system at a sampling frequency of 20 Hz. Electrophysiological signals were acquired at 128 Hz for sleep scoring analysis as previously reported (Opp ICELUS Acquisition program: [55]. Scored sleep states include rapid eye movement (REM) sleep, slow-wave sleep (SWS) and waking (WAKE) periods.
Cortical neuron model
Neuron membrane potential dynamics were described by a Hodgkin-Huxley based model of cholinergic modulation in pyramidal cells [21,56]. The effects of ACh as mediated through muscarinic ACh receptors were shown to be well modeled by varying the maximum conductance gKs of a slow, low-threshold K+ mediated adaptation current. The model also featured a fast, inward Na+ current, a delayed rectifier K+ current and a leakage current. With C = 1μF /cm2, units of Vi being millivolts and units of t being milliseconds, the current balance equation for the ith cell was
where a constant current
was externally applied and
represented the synaptic current received by the ith neuron.
For some simulations we added noisy input current pulses dictated by a Poisson process (Poisson Rate
) with amplitude of 6μA/cm2 and duration 1ms.
Activation of the inward Na+ current was instantaneous and governed by the steady-state activation function mi,∞(Vi) = {1+exp[(−Vi−30.0)/9.5]}−1. The Na+ inactivation gating variable hi was governed by
where h∞(V) = {1+exp[(V+40.5)/6.0]}−1 and τh(V) = 0.37+2.78{1+exp[(V+40.5)/6.0]}−1.
The delayed rectifier K+ current was gated by ni, the dynamics of which was given by
where n∞(V) = {1+exp[(−V−30.0)/10.0]−1 and τn(V) = 0.37+1.85{1+exp[(V+27.0)/15.0]}−1.
To model ACh blockade of the muscarine-sensitive M-current observed in cortical neurons, the maximum conductance of the slow, low-threshold K+ current in the ith cell, , was varied between 1.5 mS/cm2 for no ACh modulation and 0 mS/cm2 for strong ACh modulation. In this model neuron, decreasing values of gKs increase membrane excitability as reflected in the frequency-current relation (Fig 8), as well as affect spike-frequency adaptation and the neural phase response curve [56,57]. The dynamics of the corresponding gating variable zi were governed by
where z∞(V) = {1+exp[(−V−39.0)/5.0]}−1. These M-current kinetics are similar to M-current models in [58] and [59]. Values of other parameters were: gNa = 24.0mS/cm2,
, gL = 0.02mS/cm2, VNa = 55.0mV, VK = −90.0mV and VL = −60.0mV, τz = 75ms.
A, Frequency-current (f-I) curve of the cortical neuron model with different gKs values simulating different levels of ACh neuromodulation (gKs = 0 mS/cm2 for high ACh modulation and gKs = 1.5 mS/cm2 for no ACh neuromodulation). B, Spike-Frequency-Adaptation (SFA) showed by the voltage traces of different gKs values with current of 1.5μA/cm2 being applied (Same color code with A).
Network model
We simulated two-dimensional networks with 400 excitatory (E) neurons and 100 inhibitory (I) neurons evenly distributed over separate square lattices (20×20 E cell lattice and 10×10 I cell lattice, Fig 9). The inhibitory cells accounted for 20% of cells similar to what has been reported experimentally in the cortex [60]. A local excitation-global inhibition network topology (similar to center-surround or lateral inhibition topologies) was used in which E cells sent outgoing connections to their 40 nearest neighbors on the E cell lattice and to their 10 nearest neighbors on the I cell lattice. Inhibitory cells sent outgoing connections to all E cells and all I cells. Periodic boundary conditions were imposed on cells near the lattice edges. This is an established network scheme for cortical connectivity with the ability to balance short-range excitation and global inhibition [61].
Network consisted of 100 inhibitory neurons (blue top layer) and 400 excitatory neurons (red middle layer) evenly distributed over square lattices periodic boundary conditions. E-E and E-I synapses were short range (red arrows) while I-I and I-E synapses (blue arrows) were global. Bottom layer shows an example of a heterogeneous spatial distribution mapping (yellow to blue indicates high to low
values).
To investigate how the observed results persist with random but sparse inhibitory connectivity, we performed additional simulations where we reduced the density of inhibitory connectivity to as low as 40% (S15 Fig). The observed network dynamics did not change qualitatively.
To illustrate network dynamics on a raster plot, we indexed neurons by lattice column such that a neuron’s index. IDi, was set to the sum of its lattice y-coordinate and the product of its lattice x-coordinate with the length of lattice network, IDi = yi+xi×L (L =20 for E-cells and L = 10 for I-cells). The first 400 indices were assigned to E-cells while I-cells’ indices ranged from 401 to 500.
The synaptic current represented the total synaptic current received by neuron i and was given by
where
at times t>tjk (spike time of jth neuron’s kth spike). wij is the ijth element in the adjacency matrix for the weighted directed graph for synaptic connections in our network model. We used 0.01 mS/cm2, 0.05 mS/cm2, 0.04 mS/cm2 and 0.04 mS/cm2 for E-E, E-I, I-I and I-E synaptic strengths, respectively. For all synapses we used the same decay time constant τ = 3.0 ms. The reversal potential of the synaptic current (
) was set to 0 mV for excitatory synapses and -75mV for inhibitory synapses.
Generation of heterogenous ACh spatial maps
To simulate spatially heterogeneous distribution of ACh levels, we constructed a mapping of maximal conductance values gKs across the E cell and I cell lattices (Fig 9, bottom layer). The values for E cells, based on their i, j position in the 20×20 lattice, were computed by an iterative process that approximates locally diffusive spread of
in response to a point source release with subsequent decay. The iterative process and resulting
mapping were computed prior to simulation of neural network activity and remained constant throughout the network simulation. In the iterative process, initially, all
values were set to 1.5 mS/cm2 and the sites of point source release were chosen.
The value of at iteration step n+1 was computed by:
where the first term with coefficient D as the diffusion constant represented discrete diffusion on the 2D lattice network, Rij simulated the effect of ACh release (modeled by a decreased
level) and B was the decay constant to represent the effect of ACh decay. Periodic boundary conditions were imposed for cells at the lattice edges. The iteration time-step was Δt. For simplification we set Δt, Δx, Δy = 1 for unit time step and unit distances on the lattice. The
levels for the I cells were assigned as the average value of the 4 nearest E cells in the 2×2 block centered on the I cell. Simulations shown in S1 Fig hold gKs fixed at 0 mS/cm2 across all I cells, without qualitatively changing results.
In all mappings, we set the lower bound of
to be 0.2 mS/cm2 and the upper bound to 1.5 mS/cm2, if not otherwise specified. The iterative process was frozen at different time steps to yield low
(high ACh) regions with different radii. The frozen
mappings dictated the
values used during the simulations of neural network activity.
Measurements of dynamics
All results are averages over four simulation replications from random initial conditions, if not specified otherwise. To detect the network-wide oscillatory rhythms, we performed Fourier transforms of the correlogram of the spiking raster plot. We used peak power in the range 2.5−20 Hz to detect theta band power and in the 25−100 Hz range to detect gamma band activity.
We note that gamma oscillations specifically constitute a wide range of frequencies. The specific realization of gamma frequency will depend on specific inhibitory synaptic time-constant used and the drive that the cell receives (internal, coming from the local network, or external from other brain modality).
To detect an individual cell’s rhythmic activity, we performed Fourier transforms of the autocorrelation function of neuronal spiking time-series data. Power in the theta and gamma frequency bands were compared to the average power over all frequencies to identify if a cell exhibited theta or gamma rhythms.
Local field potentials and coupling strength
To simulate the local field potential (LFP) measurements at different locations in the network, we first convolved the discrete neuronal spiking times with a Gaussian function to generate continuous spike traces for each neuron (Vi(t) for ith neuron). The Gaussian filter had a σ of 1.5 ms. To compute the LFP at a particular location in the network, spike traces for the E cell at the location and for its 12 nearest neighbor E cells were summed. For the results shown in S13 Fig, we constructed LFP traces by directly summing the voltage traces of these cells.
The LFP trace was filtered at the two frequency ranges determined by the peak frequencies within theta and gamma bands, respectively. We denoted the theta and gamma band filtered LFP signals by xθ(t) and xγ(t) respectively. The amplitude envelope of xγ(t), denoted as Aγ(t), was extracted using a standard Hilbert transform. To quantify the phase coupling of theta-gamma rhythms, the modulation index (MI) for phase-amplitude coupling between theta and gamma bands was computed as described in [62].
Supporting information
S1 Fig. Network dynamics with no M-current in the I-population.
For this simulation, the gKs values of I-population are 0 mS/cm2 and E-cells have the same gKs values as in Fig 4H (radius at 6.1 and distance between two spots is 8 units.) The I-I and I-E synaptic strengths were adjusted to 0.06 and 0.035 mS/cm2 respectively. A: Spike raster plot illustrating E cell (cells 1–400) and I cell (401–500) firing patterns. B: Dominant rhythmic activity of individual E-cells (dark blue = none, light blue = theta band, green = gamma band and yellow = mixed, both gamma and theta) plotted at cell position on the E-cell lattice.
https://doi.org/10.1371/journal.pcbi.1009235.s001
(TIF)
S2 Fig. Network dynamics with homogeneous gKs maps.
A: The number of neurons primarily exhibiting gamma (green curve) or theta (blue curve) rhythms as a function of the gKs value that is uniform in all cells in the network. B: The rhythm power of network dynamics in theta band (blue curve) and gamma band (green curve) as a function of the network gKs value. C, D: Spike raster plots illustrating E cell (cells 1–400) and I cell (401–500) firing patterns with network gKs values at 0.2 (C) and 1.2 (D) mS/cm2.
https://doi.org/10.1371/journal.pcbi.1009235.s002
(TIF)
S3 Fig. Cross-sections of single peak gKs distribution.
A plot illustrating the cross-sections of gKs values with different radii in the simulations of single peak gKs distributions; (Fig 3.)
https://doi.org/10.1371/journal.pcbi.1009235.s003
(TIF)
S4 Fig. Gamma and theta band firing varies with gKs level within a single peak distribution and with neural external drive current.
For a single gKs hotspot spatial mapping (r = 5.6), the lower bound of gKs reached at the center of the hotspot (x-axis, in mS/cm2, gKs upper bound was set to 1.5 mS/cm2) and the external input current to all neurons, (y-axis) was varied. Panels show measures of network theta power (A), network gamma power (B), ratio of theta to gamma power (C), and numbers of cells primarily firing in the theta frequency band (D), in the gamma frequency band (E) and with power in both bands (mixed, F). Network theta power and the number of cells primarily exhibiting theta rhythmicity were sensitive to level of gKs in its spatial distribution. Specifically, for a single gKs hotspot with radius r = 5.6 gamma/theta band activity depended on the minimum gKs value within the hotspot and the external input current
applied to the neurons. Smaller values for the lower bound of gKs increased the difference in neuron modulation within the hotspot compared to outside the hotspot, and larger values of the input current promoted network-wide excitability (i.e. not limited to gKs hot spots), leading to global strengthening in theta/gamma power (top/right rows/columns). For dynamics localized to discrete spots of activity (bottom/left-center rows/columns), increased network power in the gamma band, and higher numbers of cells primarily exhibiting gamma activity, occurred for the lower minimum gKs values.
https://doi.org/10.1371/journal.pcbi.1009235.s004
(TIF)
S5 Fig. Theta and gamma band frequencies vary with the size of gKs hotspot within a single peak distribution.
Here the gKs lower bound (i.e. within the hotspot) is gKs = 0.2 mS/cm2, while its upper bound (i.e background value) is gKs = 1.5mS/cm2. The external input current to all neurons, .
https://doi.org/10.1371/journal.pcbi.1009235.s005
(TIF)
S6 Fig. Effect of E cell connectivity topology on network dynamics with single peak gKs distribution.
A, Power of network theta and gamma rhythms as a function of E cell synaptic connection rewiring probability with a single gKs hotspot (r = 5.7). It shows the drastic decrease of theta rhythm across the network as network E-cells’ connections became less local. B, The prominent frequency in gamma band as a function of E cell rewiring probability. Inset shows gKs spatial mapping on the E cell 2D lattice for the single hotspot (r = 5.7). C and D, Examples of spiking raster plots with E cell rewiring probability at p = 0 and p = 0.875, respectively, demonstrating the shift from theta-gamma coupled activity to a synchronized gamma rhythm. E cells are numbered 1 to 400, and I cells are numbered 401 to 500. Color indicates gKs values of cells with the scale in the inset in B. We investigated how network topology affected the observed oscillatory rhythms with a single peak gKs spatial distribution. To this effect, we progressively rewired initially local E-E connections to random E cells across the network. The rewiring probability, p, (x-axis on S6A Fig and B) denotes the fraction of E-E connections rewired: when p = 0 the network has the original local excitation/global inhibition connectivity, whereas for p = 1 the network has random excitation spanning the whole network. We observed two major effects as a function of the increased rewiring. First, theta power was significantly diminished for p > 0.1 while gamma power remained relatively high across all rewired connectivities (S6A Fig). This observation underscores the importance of local excitatory connectivity in supporting localized firing within the gKs hotspot. Secondly, the frequency of gamma oscillations almost linearly decreased with increasing rewiring (S6B Fig). This is due to the fact, that the random connectivity mediates emergence of zero-phase synchrony between the excitatory neurons. This in turn makes the presynaptic spike arrive on the target excitatory cells at the time when they are (partially) in their refractory times, reducing the impact the EPSP has on these cells. This finally decreases the amount of overall excitation in the network reducing the frequency of gamma. The two rasterplots (S6C and S6D Fig) exemplify these observations.
https://doi.org/10.1371/journal.pcbi.1009235.s006
(TIF)
S7 Fig. Changes in network dynamics as a function of M-current time constant, τz.
For the simulations, the double peaked gKs spatial mapping is the same as in Fig 4H (the spot radius is r = 6.1 and distance between two spots is d = 8 units.) (A) The frequency in theta band decreased as the M-current time constant, τz, increased. (B) The frequency in gamma band largely didn’t change as the M-current time constant was increased. (C), (D) Spike raster plots illustrating E cell (cells 1–400) and I cell (401–500) firing patterns when M-current time constant was 25 ms and 125 ms, respectively.
https://doi.org/10.1371/journal.pcbi.1009235.s007
(TIF)
S8 Fig. Theta and gamma band frequencies vary with features of a double peaked gKs spatial distribution.
Heatmaps showing the most prominent frequency in the power spectrum of network firing in the theta band (A) and in the gamma band (B) as the radius r of the gKs hotspots and distance d between hotspot centers is varied. Labels F, G and H correspond to same labels in Fig 4(main). White squares indicate networks without significant power in theta or gamma frequency bands. The frequency of network theta and gamma band activity changed as the radius r and distance d between gKs hotspot centers were varied in a double peaked gKs spatial distribution (S8 Fig). Theta band activity decreased in frequency as the distance between the hotspot centers increased. This was due to fact that the amount of bleed-over excitation from the active hotspot to the silent hotspot decreased as a function of distance between the hotspots. This bleed-over excitation increased excitatory input to the inactive hotspot, subsequently allowing for faster switching between the hot spots when they were close. Gamma band frequency was approximately inversely correlated to the number of cells in the network predominantly exhibiting gamma firing frequency (compare S8B Fig with Fig 4 in main part). In particular, when many cells fired predominantly at gamma frequency, more E cells recruited their local I cells into the PING inhibitory gating that is signaled globally in the network. This stronger inhibitory gating slowed the release of E cells from inhibition and thus network gamma band activity.
https://doi.org/10.1371/journal.pcbi.1009235.s008
(TIF)
S9 Fig. Effect of network topology on theta-gamma coupling in a double peak gKs distribution.
A, Power of network theta band activity decreased with increased probability of rewiring synapses between E cells for a double peak gKs mapping with hotspot radius r = 5.4 and distance between hotspot centers d = 6. The introduction of E-E synaptic connections between the different hotspots allowed cells to fire at the same time and synchronously. B and C, Spike raster plots for E cell rewiring probability at 0 and 0.875, respectively. E cells are numbered 1 to 400, and I cells are 401 to 500. Color indicates gKs values of cells with the scale in the inset in A. As in the single peak gKs spatial distribution, theta-gamma coupled firing activity was sensitive to changes in the local excitation, global inhibition connectivity structure of the network. In this case, we considered a double peak gKs spatial mapping which exhibited strong theta-gamma coupling: (r = 5.4, d = 6) and randomly rewired E-E synapses with varying probability. We observed drastic decreases in network theta power with increased E cell rewiring probability (S9 Fig). As rewiring probability increased, network dynamics changed from theta-modulated gamma band activity occurring in each gKs hotspot (rewiring probability at 0; S9B Fig) to synchronized gamma band activity in both gKs hotspots (rewiring probablity at 0.875; S9C Fig). This was caused by an increase in E-E connectivity between cells in different hotspots. The resulting additional excitation from the active hotspot enabled cells in the silent hotspot to overcome the global inhibition and fire at the same time as the active hotspot and in synchrony with those cells. Theta-gamma coupled activity was well maintained in the ’small-world’ network regime (when the rewiring probability was small, ~0.2). Generally, localized spots of spiking activity (and thus their switching) was obtained when the spatial extent of excitation was smaller than that of inhibition.
https://doi.org/10.1371/journal.pcbi.1009235.s009
(TIF)
S10 Fig. Network dynamics on homogeneous gKs spatial map with two hotspots created by increased constant current drive.
A) Illustration of applied current distribution for corresponding neurons on the 20 × 20 E cell lattice. Red color indicates neuros receiving and black color indicates cells receiving
. B-I, Spiking raster plots (top panels) and network frequency power spectrums (bottom panels) for homogeneous gKs values (as denoted above the panels) between 0 and 1.4 mS/cm2.
https://doi.org/10.1371/journal.pcbi.1009235.s010
(TIF)
S11 Fig. Variability in network rhythmic activity in randomly generated gKs maps.
A, B and C: Power of network theta (A) and gamma (B) rhythms and their ratio (C) computed from networks with randomly generated gKs spatial mappings when the number of hotspot centers was varied from 1 to 20 (y-axis) and hotspot radius r was varied from 2.8 to 5.4 (x-axis). Positions of hotspot centers were randomly chosen on the excitatory cell lattice. Results were averaged over 4 realizations of the gKs mapping with different positions of hotspot centers. The ’star’ marker corresponds to parameters for Fig 5A and 5B in main part and the ’square’ marker corresponds to parameters for Fig 5C and 5D in main part. D, E and F, The relative standard error (RSE) of the power of network theta (D) and gamma (E) rhythms and their ratio (F) computed from the 4 realizations of the randomly generated gKs spatial mapping with varying hotspot number (y-axis) and radius (x-axis). To systematically consider spatially random gKs distributions, we varied the number of gKs hotspots and their radius r, and then generated multiple gKs spatial mappings with different locations of hotspot centers. Network power in the theta and gamma bands, as well as theta-gamma power ratio, averaged over simulations from 4 realizations of the gKs mapping, varied widely. This was due to high variation in network rhythmic activity generated across the 4 gKs mappings realizations. Computation of the relative standard error (RSE) in power of network activity in the theta and gamma frequency bands across the gKs mapping realizations showed that gamma band power was generally similar, but theta band power and, thus, theta-gamma power ratio, showed higher variability across mapping realizations.
https://doi.org/10.1371/journal.pcbi.1009235.s011
(TIF)
S12 Fig. Theta and gamma band frequencies vary with features of randomly generated gKs maps.
Heatmaps show the most prominent frequency in the power spectrum of network firing in the theta band (A) and in the gamma band (B) as the number of hotspot centers was varied from 1 to 20 (y-axis) and hotspot radius r was varied from 2.8 to 5.4 (x-axis).
https://doi.org/10.1371/journal.pcbi.1009235.s012
(TIF)
S13 Fig. Strength of theta-gamma coupling as a function of the distance from the low gKs region.
Modulation Index (MI) between gamma and theta filtered LFP traces as a function of the distance from the center of the gKs hotspot as in Fig 6. The only difference is here we use the sum of cell voltage traces for the LFP calculation.
https://doi.org/10.1371/journal.pcbi.1009235.s013
(TIF)
S14 Fig. Changes in theta and gamma frequencies as a function of cellular signaling properties.
For these simulations, the double peaked gKs spatial mapping is the same as in Fig 4H (radii r = 6.1 and distance between two spots d = 8 units). A) The theta band frequency increased as the E-cell external current input was increased. B) The gamma band frequency increased as the E-cell external current
input was increased. C) The theta band frequency increased as the decay time constant for inhibitory synapses τI was increased. D) The gamma band frequency decreased as the decay time constant for inhibitory synapses τI was increased. For C and D we kept the product of inhibitory synaptic strength and τI constant in order to achieve similar inhibition level across the set of simulations.
https://doi.org/10.1371/journal.pcbi.1009235.s014
(TIF)
S15 Fig. Network dynamics with sparse inhibitory connectivity.
For the simulations, the double peaked gKs spatial mapping is the same as in Fig 4H (radii is 6.1 and distance between two spots is 8 units.) Spike raster plots illustrating E cell (cells 1–400) and I cell (401–500) firing patterns. A, Default inhibitory synaptic connectivity (all-to-all, 0.04 mS/cm2 as inhibitory synaptic strength, detailed in Materials and Methods section). B, Random connectivity with 80 percent of default inhibitory synaptic density, inhibitory synaptic strengths are adjusted to 0.045 mS/cm2. C, Random connectivity with 60 percent of default inhibitory synaptic density, inhibitory synaptic strengths are adjusted to 0.048 mS/cm2. D, Random connectivity with 40 percent of default inhibitory synaptic density, inhibitory synaptic strengths are adjusted to 0.075 mS/cm2.
https://doi.org/10.1371/journal.pcbi.1009235.s015
(TIF)
References
- 1. Zaborszky L, Buhl DL, Pobalashingham S, Bjaalie JG, Nadasdy Z. Three-dimensional chemoarchitecture of the basal forebrain: spatially specific association of cholinergic and calcium binding protein-containing neurons. Neuroscience. 2005;136: 697–713. pmid:16344145
- 2. Zaborszky L, Duque A, Gielow M, Gombkoto P, Nadasdy Z, Somogyi J. Organization of the basal forebrain cholinergic projection system: specific or diffuse? The rat nervous system. Elsevier; 2015. pp. 491–507.
- 3. Gielow MR, Zaborszky L. The input-output relationship of the cholinergic basal forebrain. Cell Rep. 2017;18: 1817–1830. pmid:28199851
- 4. Yuan R, Biswal BB, Zaborszky L. Functional subdivisions of magnocellular cell groups in human basal forebrain: test—retest resting-state study at ultra-high field, and meta-analysis. Cereb Cortex. 2019;29: 2844–2858. pmid:30137295
- 5. Hasselmo ME. Neuromodulation: Acetylcholine and memory consolidation. Trends in Cognitive Sciences. Elsevier Ltd; 1999. pp. 351–359. pmid:10461198
- 6. Hasselmo ME, Giocomo LM. Cholinergic modulation of cortical function. Journal of Molecular Neuroscience. J Mol Neurosci; 2006. pp. 133–135. pmid:17192659
- 7. Dayan P, Yu AJ. Expected and unexpected uncertainty: ACh and NE in the neocortex. Adv Neural Inf Process Syst. 2002;15: 173–180.
- 8. Sarter M, Lustig C. Forebrain Cholinergic Signaling: Wired and Phasic, Not Tonic, and Causing Behavior. J Neurosci. 2020;40: 712–719. pmid:31969489
- 9. Disney AA, Higley MJ. Diverse Spatiotemporal Scales of Cholinergic Signaling in the Neocortex. J Neurosci. 2020;40: 720–725. pmid:31969490
- 10. Parikh V, Kozak R, Martinez V, Sarter M. Prefrontal acetylcholine release controls cue detection on multiple timescales. Neuron. 2007;56: 141–154. pmid:17920021
- 11. Börgers C, Epstein S, Kopell NJ. Background gamma rhythmicity and attention in cortical local circuits: A computational study. Proc Natl Acad Sci U S A. 2005;102: 7002–7007. pmid:15870189
- 12. Lu Y, Sarter M, Zochowski M, Booth V. Phasic cholinergic signaling promotes emergence of local gamma rhythms in excitatory-inhibitory networks. Eur J Neurosci. 2020. pmid:32293081
- 13. Roach JP, Ben-Jacob E, Sander LM, Zochowski MR. Formation and dynamics of waves in a cortical model of cholinergic modulation. PLoS Comput Biol. 2015;11: e1004449. pmid:26295587
- 14. Colgin LL. Theta-gamma coupling in the entorhinal-hippocampal system. Current Opinion in Neurobiology. Elsevier Ltd; 2015. pp. 45–50. pmid:25168855
- 15. Howe WM, Gritton HJ, Lusk NA, Roberts EA, Hetrick VL, Berke JD, et al. Acetylcholine release in prefrontal cortex promotes gamma oscillations and theta–gamma coupling during cue detection. J Neurosci. 2017;37: 3215–3230. pmid:28213446
- 16. Newman EL, Gillet SN, Climer JR, Hasselmo ME. Cholinergic blockade reduces theta-gamma phase amplitude coupling and speed modulation of theta frequency consistent with behavioral effects on encoding. J Neurosci. 2013;33: 19635–19646. pmid:24336727
- 17. Obermayer J, Heistek TS, Kerkhofs A, Goriounova NA, Kroon T, Baayen JC, et al. Lateral inhibition by Martinotti interneurons is facilitated by cholinergic inputs in human and mouse neocortex. Nat Commun. 2018;9: 1–14. pmid:29317637
- 18. Whittington MA, Traub RD, Kopell N, Ermentrout B, Buhl EH. Inhibition-based rhythms: experimental and mathematical observations on network dynamics. Int J Psychophysiol. 2000;38: 315–336. pmid:11102670
- 19. Börgers C, Kopell N. Effects of noisy drive on rhythms in networks of excitatory and inhibitory neurons. Neural Comput. 2005;17: 557–608. pmid:15802007
- 20. Börgers C, Kopell NJ. Gamma oscillations and stimulus selection. Neural Comput. 2008;20: 383–414. pmid:18047409
- 21. Stiefel KM, Gutkin BS, Sejnowski TJ. The effects of cholinergic neuromodulation on neuronal phase-response curves of modeled cortical neurons. J Comput Neurosci. 2009;26: 289–301. pmid:18784991
- 22. Bosman CA, Schoffelen JM, Brunet N, Oostenveld R, Bastos AM, Womelsdorf T, et al. Attentional Stimulus Selection through Selective Synchronization between Monkey Visual Areas. Neuron. 2012;75: 875–888. pmid:22958827
- 23. Grothe I, Neitzel SD, Mandon S, Kreiter AK. Switching neuronal inputs by differential modulations of gamma-band phase-coherence. J Neurosci. 2012;32: 16172–16180. pmid:23152601
- 24. Gritton HJ, Howe WM, Mallory CS, Hetrick VL, Berke JD, Sarter M. Cortical cholinergic signaling controls the detection of cues. Proc Natl Acad Sci U S A. 2016;113: E1089–E1097. pmid:26787867
- 25. Parikh V, Sarter M. Cholinergic mediation of attention: Contributions of phasic and tonic increases in prefrontal cholinergic activity. Ann N Y Acad Sci. 2008;1129: 225–235. pmid:18591483
- 26. Jing M, Li Y, Zeng J, Huang P, Skirzewski M, Kljakic O, et al. An optimized acetylcholine sensor for monitoring in vivo cholinergic activity. Nat Methods. 2020;17: 1139–1146. pmid:32989318
- 27. Jing M, Zhang P, Wang G, Feng J, Mesik L, Zeng J, et al. A genetically encoded fluorescent acetylcholine indicator for in vitro and in vivo studies. Nat Biotechnol. 2018;36: 726–737. pmid:29985477
- 28. Fino E, Packer AM, Yuste R. The logic of inhibitory connectivity in the neocortex. Neuroscientist. SAGE PublicationsSage CA: Los Angeles, CA; 2013. pp. 228–237. pmid:22922685
- 29. Hasselmo ME, Anderson BP, Bower JM. Cholinergic modulation of cortical associative memory function. J Neurophysiol. 1992;67: 1230–1246. pmid:1597709
- 30. Tiesinga PHE, Fellous J-M, José J V, Sejnowski TJ. Computational model of carbachol-induced delta, theta, and gamma oscillations in the hippocampus. Hippocampus. 2001;11: 251–274. pmid:11769308
- 31. Buzsáki G. Theta oscillations in the hippocampus. Neuron. Cell Press; 2002. pp. 325–340. pmid:11832222
- 32. Deco G, Thiele A. Cholinergic control of cortical network interactions enables feedback-mediated attentional modulation. Eur J Neurosci. 2011;34: 146–157. pmid:21692884
- 33. Colgin LL. Theta—gamma coupling in the entorhinal—hippocampal system. Curr Opin Neurobiol. 2015;31: 45–50. pmid:25168855
- 34. Jensen O, Colgin LL. Cross-frequency coupling between neuronal oscillations. Trends Cogn Sci. 2007;11: 267–269. pmid:17548233
- 35. Colgin LL, Denninger T, Fyhn M, Hafting T, Bonnevie T, Jensen O, et al. Frequency of gamma oscillations routes flow of information in the hippocampus. Nature. 2009;462: 353. pmid:19924214
- 36. Tort ABL, Komorowski RW, Manns JR, Kopell NJ, Eichenbaum H. Theta—gamma coupling increases during the learning of item—context associations. Proc Natl Acad Sci. 2009;106: 20942–20947. pmid:19934062
- 37. Fries P, Reynolds JH, Rorie AE, Desimone R. Modulation of oscillatory neuronal synchronization by selective visual attention. Science (80-). 2001;291: 1560–1563. pmid:11222864
- 38. Canolty RT, Knight RT. The functional role of cross-frequency coupling. Trends Cogn Sci. 2010;14: 506–515. pmid:20932795
- 39. Aton SJ. Set and setting: How behavioral state regulates sensory function and plasticity. Neurobiology of Learning and Memory. Neurobiol Learn Mem; 2013. pp. 1–10. pmid:23792020
- 40. Rodriguez R, Kallenbach U, Singer W, Munk MHJ. Short- and long-term effects of cholinergic modulation on gamma oscillations and response synchronization in the visual cortex. J Neurosci. 2004;24: 10369–10378. pmid:15548651
- 41. Gulyás AI, Szabó GG, Ulbert I, Holderith N, Monyer H, Erdélyi F, et al. Parvalbumin-containing fast-spiking basket cells generate the field potential oscillations induced by cholinergic receptor activation in the hippocampus. J Neurosci. 2010;30: 15134–15145. pmid:21068319
- 42. White JA, Banks MI, Pearce RA, Kopell NJ. Networks of interneurons with fast and slow γ-aminobutyric acid type A (GABA(A)) kinetics provide substrate for mixed gamma-theta rhythm. Proc Natl Acad Sci U S A. 2000;97: 8128–8133. pmid:10869419
- 43. Wulff P, Ponomarenko AA, Bartos M, Korotkova TM, Fuchs EC, Bähner F, et al. Hippocampal theta rhythm and its coupling with gamma oscillations require fast inhibition onto parvalbumin-positive interneurons. Proc Natl Acad Sci U S A. 2009;106: 3561–3566. pmid:19204281
- 44. Gloveli T, Dugladze T, Rotstein HG, Traub RD, Monyer H, Heinemann U, et al. Orthogonal arrangement of rhythm-generating microcircuits in the hippocampus. Proc Natl Acad Sci U S A. 2005;102: 13295–13300. pmid:16141320
- 45. Tort ABL, Rotstein HG, Dugladze T, Gloveli T, Kopell NJ. On the formation of gamma-coherent cell assemblies by oriens lacunosum-moleculare interneurons in the hippocampus. Proc Natl Acad Sci U S A. 2007;104: 13490–13495. pmid:17679692
- 46.
Kopell N, Börgers C, Pervouchine D, Malerba P, Tort A. Gamma and Theta Rhythms in Biophysical Models of Hippocampal Circuits. In: Cutsuridis V, Graham B, Cobb S, Vida I, editors. Hippocampal Microcircuits: A Computational Modeler’s Resource Book. New York, NY: Springer New York; 2010. pp. 423–457. https://doi.org/10.1007/978-1-4419-0996-1_15
- 47. Chen Y, Martinez-Conde S, Macknik SL, Bereshpolova Y, Swadlow HA, Alonso JM. Task difficulty modulates the activity of specific neuronal populations in primary visual cortex. Nat Neurosci. 2008;11: 974–982. pmid:18604204
- 48. Chelazzi L, Miller EK, Duncan J, Desimone R. Responses of neurons in macaque area V4 during memory-guided visual search. Cereb Cortex. 2001;11: 761–772. pmid:11459766
- 49. Desimone R, Duncan J. Neural Mechanisms of Selective Visual Attention. Annu Rev Neurosci. 1995;18: 193–222. pmid:7605061
- 50. Xiang Z, Huguenard JR, Prince DA. Cholinergic switching within neocortical inhibitory networks. Science (80-). 1998;281: 985–988. pmid:9703513
- 51.
Gritton HJ, Lee TM, Opp M, M. S. Cholinergic neurotransmission during REM sleep: combined electrochemical recordings of cholinergic transients and field potentials challenge an old dogma. Society for Neuroscience Annual Meeting. 2009.
- 52. Parikh V, Pomerleau F, Huettl P, Gerhardt GA, Sarter M, Bruno JP. Rapid assessment of in vivo cholinergic transmission by amperometric detection of changes in extracellular choline levels. Eur J Neurosci. 2004;20: 1545–1554. pmid:15355321
- 53. Giuliano C, Parikh V, Ward JR, Chiamulera C, Sarter M. Increases in cholinergic neurotransmission measured by using choline-sensitive microelectrodes: Enhanced detection by hydrolysis of acetylcholine on recording sites? Neurochem Int. 2008;52: 1343–1350. pmid:18346819
- 54. Howe WM, Berry AS, Francois J, Gilmour G, Carp JM, Tricklebank M, et al. Prefrontal cholinergic mechanisms Instigating shifts from monitoring for cues to Cue-Guided performance: Converging electrochemical and fMRI evidence from rats and humans. J Neurosci. 2013;33: 8742–8752. pmid:23678117
- 55. Baracchi F, Opp MR. Sleep-wake behavior and responses to sleep deprivation of mice lacking both interleukin-1β receptor 1 and tumor necrosis factor-α receptor 1. Brain Behav Immun. 2008;22: 982–993. pmid:18329246
- 56. Fink CG, Murphy GG, Zochowski M, Booth V. A dynamical role for acetylcholine in synaptic renormalization. PLoS Comput Biol. 2013;9. pmid:23516342
- 57. Roach JP, Eniwaye B, Booth V, Sander LM, Zochowski MR. Acetylcholine Mediates Dynamic Switching Between Information Coding Schemes in Neuronal Networks. Front Syst Neurosci. 2019;13: 64. pmid:31780905
- 58. Crook SM, Ermentrout GB, Bower JM. Spike Frequency Adaptation Affects the Synchronization Properties of Networks of Cortical Oscillators. Neural Comput. 1998;10: 837–854. pmid:9573408
- 59. Olufsen MS, Whittington MA, Camperi M, Kopell N. New roles for the gamma rhythm: Population tuning and preprocessing for the beta rhythm. J Comput Neurosci. 2003;14: 33–54. pmid:12435923
- 60. Sahara S, Yanagawa Y, O’Leary DDM, Stevens CF. The fraction of cortical GABAergic neurons is constant from near the start of cortical neurogenesis to adulthood. J Neurosci. 2012;32: 4755–4761. pmid:22492031
- 61. Amari S. Dynamics of pattern formation in lateral-inhibition type neural fields. Biol Cybern. 1977;27: 77–87. pmid:911931
- 62. Tort ABL, Komorowski R, Eichenbaum H, Kopell N. Measuring phase-amplitude coupling between neuronal oscillations of different frequencies. J Neurophysiol. 2010;104: 1195–1210. pmid:20463205