Introduction

Cortical circuits are composed of sparsely connected networks of excitatory pyramidal neurons. To maintain a balanced level of excitability in the network, they form dense reciprocal connections with local inhibitory interneurons. The two largest subclasses of interneurons in the cortex are parvalbumin positive (PV) and somatostatin positive (SOM) cells [1]. A major distinction between PV and SOM interneurons is that they synapse onto different parts of the pyramidal neuron. PV tends to contact the soma and proximal dendrites (perisomatic), while SOM targets distal dendrites [2–5]. Many studies have identified differences in how they regulate the responsiveness of local excitatory principal neurons [6–8]. These differences could arise from their connectivity, intrinsic properties, or how they synapse onto pyramidal neurons. Indeed, the location of an inhibitory synapse qualitatively changes its effect on synaptic integration [5; 9-11]. Synaptic integration in pyramidal neurons arises from the interplay between passive and active ion channels and dendritic morphology [9; 12-14]. Dendrites produce regenerative spiking events, specifically Na⁺, NMDA, and Ca²⁺ spikes. Their initiation depends upon the spatiotemporal coordination of excitatory synapses, and their interaction with the depolarization of dendritic branches. Incorporating just NMDA spikes into a model neuron radically increases the complexity of its dendritic integration [15]. Na⁺ and NMDA spikes can drive somatic spiking in an in vivo-like model neuron [13]. A Ca²⁺ spike in the apical trunk converts a single somatic spike to a burst [9; 16], and may allow pyramidal neurons to act as multi-stage integrators [17; 18]. This complexity is only increased by including inhibition. The emission of bursts of action potentials is controlled by the timing and dendritic location of inhibition [9].

The rhythmicity of inhibition at distal dendrites can also affect somatic integration [19]. Reciprocal interactions between excitatory principal neurons and inhibitory interneurons gives rise to rhythmic activities (known as the PING mechanism). To initiate the rhythm, excitatory principal neurons activate inhibitory interneurons that deliver feedback inhibition to principal neurons. This transiently suppresses principal neuron firing. As the inhibition wanes, principal neurons resume their activity and reengage the inhibitory population, starting a new oscillatory cycle. Some evidence suggests that different interneuron subtypes pace rhythms at different frequencies, with gamma oscillations depending on PV [20–22], and beta oscillations on SOM [22].

Numerous cognitive processes are associated with these inhibitory rhythms in the cortex. Gamma oscillations are especially pronounced during stimulus presentation or behavioral initiation [23–25]. Beta rhythms occur during preparatory states or working memory [26; 27]. But fundamentally, these rhythms probably reflect different modes of local cortical activity and interregional communication [28–30]. For instance, gamma oscillations are associated with feedforward transmission of information in cortical circuits and beta oscillations may mediate feedback [28; 31]; but see [32].

Given all this, it is important to better understand how the complex integrative events in the dendritic tree are affected by the spatial targeting and rhythmicity of inhibition. Perhaps some rhythmic timescales are more effective in the dendrites, and others at the soma. Investigating this is beyond the reach of conventional experimental techniques. Moreover, computational models exploring how beta and gamma influence information processing have mostly relied upon simplified model neurons (e.g., [33; 34]). Thus, using a morphologically and biophysically detailed model layer 5 pyramidal neuron with Na⁺, NMDA, and Ca²⁺ spikes, we investigated how perisomatic and distal dendritic inhibition impacted dendritic integration, and the modulation of these effects by beta and gamma rhythmicity.

Results

Construction of a model cortical layer 5 pyramidal neuron

To study the effect of inhibitory rhythms on synaptic integration and somatic spiking, we adapted a previously published morphologically and biophysically detailed model of a cortical layer 5 pyramidal neuron (Fig. 1A; See methods for details). Modifications to the model were done in accordance with the published literature. In brief, this model featured a multicompartmental dendritic tree that produced dendritic Na⁺, NMDA, and Ca²⁺ spikes, along with somatic action potentials that could backpropagate (Fig. 1B). We distributed conductance-based synapses across the dendritic and somatic compartments with an average density of 2.16 excitatory and 0.22 inhibitory contacts per µm. Presynaptic drivers of excitatory synapses were drawn from a pool of 5,200 point process sources that emulated correlated afferent drive with 2-8 synaptic contacts from the same presynaptic neuron. Inhibitory synapses were divided into two populations, those targeting the soma and proximal 100 µm of the dendrites (referred to as perisomatic), and those synapsing outside that area (referred to as distal). To capture excitatory/inhibitory (E/I) balance, a hallmark of cortical activity, the rate of inhibitory synaptic drive was a rescaled version of the rate of excitatory drive, lagged by 4 ms to emulate feedforward inhibition. Naturalistic presynaptic drive elicited a median firing rate of 5.3 Hz, in agreement with in vivo rates in cortex (Fig. 1A inset, [35]).

Figure 1:

Dendrites were endowed with the following voltage-dependent conductances: a fast-inactivating Na⁺ current (I_NaT), muscarinic K⁺ current (I_m), fast non-inactivating K⁺ current (I_Kv3.1), high voltage-activated Ca²⁺ current (I_{Ca_HVA}), low voltage-activated Ca²⁺ current (I_{Ca_LVA}), and Ca²⁺ activated K⁺ current (I_SK). As a result, the basal and apical dendrites could generate Na⁺ and NMDA spikes (Fig. 1B; [13]). Dendritic Na⁺ spikes were regenerative events lasting less than 1 ms that were not preceded by somatic action potentials (Fig. 1C1; [36]). A Na⁺ spike was detected when the dendritic Na⁺ channel conductance (gNa) was higher than 0.3 mS/cm², except when this threshold was reached within 5 ms after a somatic action potential, to distinguish from backpropagating action potentials.

NMDA spikes occur when adjacent NMDA-bearing synapses were synergistically recruited by a combination of glutamatergic activation and local depolarization (Fig. 1C2; [18; 37]). They typically lasted between 20 and 80 ms. They were detected when a compartment’s membrane voltage exceeded -40 mV for at least 26 ms and NMDA current exceeded -100 pA ms of charge.

Ca²⁺ spikes are depolarizations generated at the nexus of the apical trunk upon activation of voltage-gated Ca²⁺ channels (Fig. 1C3; [38; 39]). They lasted between 20 and 50 ms. To detect them, the membrane potential must exceed -40 mV for at least 26 ms, and the combined Ca²⁺ currents (LVA and HVA) had to be 1.3 times higher than when the voltage criterion was reached (t_v). The Ca²⁺ spike ended when this value fell to 1.15 times its value at t_v.

Altogether, under conditions that mirror in vivo afferent drive, our model reproduces the dendritic spikes of a layer 5 pyramidal neuron.

Relationship between dendritic and somatic spikes

Dendritic spikes induced by synaptic activity are the principal drivers of somatic action potentials. Prior experimental and modeling work has found this to be the case for layer 2/3 and layer 5 pyramidal neurons [13; 18; 40-42]. Layer 5 pyramidal neurons have a substantially longer apical trunk, which increases the electrotonic distance of their apical tuft from the soma (Fig. 2A), and diminishes the ability of tuft synapses to elicit action potentials. Voltage-gated Ca²⁺ channels at the apical nexus compensate for this by producing a robust Ca²⁺ spike that drives a burst of action potentials at the soma [17; 18]. Underscoring that Ca²⁺ spikes compensate for morphology, pyramidal neurons with shorter apical dendrites have weaker Ca²⁺ spikes that only elicit a single spike [43–45]. We thus assessed how dendritic spikes at different electronic distances from the soma were related to action potential generation.

Figure 2:

Dendritic compartments differed in their degree of passive electrical coupling to the soma (i.e., electrotonic distance; Fig. 2A). We measured this by injecting a 20 Hz sinusoidal current in each dendritic compartment, and then calculating the ratio of the membrane voltage response at the soma over that at the dendrite. The apical trunk exhibited a relatively small attenuation ratio of ∼10%. Progressing distally into the apical tuft, attenuation reached 0.1%. Basal dendrites showed greater attenuation than the apical trunk due to their smaller diameter. However, with an attenuation ratio reaching ∼1%, their distal tips were still electrotonically closer than the apical tuft. Thus, large, and long-lasting changes in membrane potential, like those produced by dendritic spikes, would be required to have any effect on the soma.

We examined this by measuring the spike-triggered average between somatic action potentials and the presence of dendritic spikes across the dendritic tree. To simplify the complex geometry of our model neuron, dendritic compartments were grouped into deciles by their electrotonic distance, and whether they were on apical or basal branches. For each time lag from the somatic action potential, we measured the percent change in dendritic spike incidence from the mean rate across the entire simulation.

Dendritic Na⁺ spikes increased 2-3 ms prior to somatic action potentials in both basal and apical dendrites (Fig. 2B). This relationship was strongest for the compartments nearest the soma, with the rate of Na⁺ spikes in the apical trunk increasing 300% over baseline prior to somatic action potentials, and 100% in basal compartments. This relationship fell off as the dendritic spikes moved farther away from the soma, indicating that Na⁺ spikes in distal branches had little direct influence on somatic spiking.

The incidence of NMDA spikes increased ∼25 ms prior to somatic action potentials, much earlier than seen with dendritic Na⁺ spikes (Fig. 2C). But, like dendritic Na⁺ spikes, NMDA spikes in the apical branch had a stronger coupling with somatic spiking than those in basal branches, and this effect dropped with distance from the soma. NMDA spikes in the apical branch persisted after an action potential, while those in basal dendrites did not, potentially because of the action potential after-hyperpolarization.

Ca²⁺ spikes originate at the nexus of the apical dendrite. This region is electrotonically close to the entire apical trunk, facilitating the propagation of Ca²⁺ spikes (Fig. 3A). In our model, Ca²⁺ spike occurrence increased within 20 ms of somatic action potentials (Fig. 3B). Furthermore, we found that NMDA spikes in the apical dendrites tended to precede Ca²⁺ spikes (Fig. 3C). Our ionic current-based detection criteria allowed us to readily distinguish between these two phenomena, despite their similar membrane voltage profiles. Since NMDA spikes in the apical tuft normally have a weak relationship to somatic spiking (Fig. 2D), they may elicit somatic spiking indirectly by driving Ca²⁺ spikes. Put another way, the apical nexus may serve as a thresholded nonlinearity for NMDA spikes in the apical tuft to drive action potentials [18]. To test this, we measured how a Ca²⁺ spike changed the spike-triggered average between apical tuft NMDA spikes and action potentials (Fig. 3D, top). This revealed that action potentials preceded by a Ca²⁺ spike (by up to 20 ms) had increased coupling with apical NMDA spikes. No such change was seen in basal dendrites (Fig. 3D, bottom).

Figure 3:

Effect of perisomatic and distal dendritic inhibition on controlling response gain

Subtypes of inhibitory interneurons synapse on distinct dendritic zones. PV interneurons mainly synapse perisomatically, while SOM interneurons target distal dendrites [2–5]. These differences may affect their modulation of synaptic integration [46], either by shifting the threshold for evoking action potentials (a subtractive effect) or altering the slope of the relationship between excitation and firing rate (a divisive effect). There is some disagreement about the degree to which PV and SOM interneurons produce either of these effects [6–8], depending on the source of excitation and local circuitry [47]. So, before we applied beta and gamma rhythmic inhibition in our model, we studied the effect of tonically activating PV- and SOM-like inhibitory synapses.

To examine this, we doubled the rate of the inhibitory presynaptic drive onto either the PV-targeted perisomatic compartments (somatic and dendritic compartments within 100 µm), or the SOM-targeted distal (>100 µm) dendritic branches. Both decreased the firing rate of the pyramidal cell from 5.5 Hz to less than 1 Hz (Fig. 4A; control: 5.5±0.85 Hz; distal: 0.20±0.15 Hz; perisomatic: 0.70±0.31 Hz; mean±SD). These comparable changes may reflect either subtractive or divisive effects and could derive from different mechanisms.

Figure 4:

To isolate these factors, we first assessed how doubling inhibition affected action potential initiation at the soma. A series of current pulses were injected into the soma to measure the relationship between firing rate and injected current (f-I curve), which captures the gain function of the neuron (Fig. 4B). This revealed that both perisomatic and dendritic inhibition shifted the current threshold for action potential initiation. Such an effect is subtractive. In addition, perisomatic inhibition decreased the slope of the f-I relationship compared with the control, which is consistent with a divisive effect.

Although perisomatic inhibition produced the strongest subtractive effect, distal dendritic inhibition reduced firing rate the most (Fig. 4A). How can we reconcile these discordant findings? One possibility is that distal inhibition substantially reduces Na⁺, NMDA, or Ca²⁺ spikes. To determine this, we examined the overall rate of dendritic spike events. Perisomatic inhibition did not affect dendritic events compared to the control condition (Fig. 4C-E). By contrast, dendritic inhibition decreased NMDA and Ca²⁺ spikes (Fig. 4D,E). Na⁺ spikes were relatively unaffected (Fig. 4C). This lack of effect may arise from a shortening of the inactivation for voltage-gated Na⁺ channels balancing out the loss of excitatory drive.

Distinct excitation/inhibition balance effects of perisomatic and distal dendritic inhibition

The previous analysis found distinct effects on neuronal gain and dendritic spiking from tonic changes in perisomatic and distal dendritic inhibition. But in vivo inhibition is dynamic and time varying. PV and SOM interneurons form dense interconnections with local pyramidal neurons, supplying time-lagged feedback inhibition. This helps maintain the balance of excitation and inhibition (E/I) in the network by ensuring that an overall increase in excitatory activity is rapidly counterbalanced by proportionate inhibition.

To emulate situations where E/I balance is important, we increased the dynamic variation in excitatory drive (see Methods). As with the previous simulations, the rate of inhibitory synaptic drive was a lagged and rescaled version of the overall excitation rate. Normally that lag is 4 ms, in line with experimental estimates [48]. To probe whether perisomatic or distal dendritic inhibition have distinct effects on E/I balance, we independently varied their lags (Fig. 4F-I). One was kept at the nominal 4 ms lag, while the other was extended. Increasing either of their lags by 500 ms produced obvious differences in the emission of dendritic spikes and their coordination with action potentials (Fig. 4F). Normally dendritic spikes are evenly distributed in time and drive somatic action potentials. Lagging distal dendritic inhibition clustered dendritic spikes together in time, during which action potentials were emitted. Lagging perisomatic inhibition did not affect the spacing of dendritic spikes, but decreased their coupling with somatic spiking.

We systematically characterized these lag effects for the following spiking events modulated by tonic changes in inhibition: action potentials, Ca²⁺, and NMDA spikes (Fig. 4G). Increasing the lag of perisomatic inhibition lowered action potential firing, while for distal dendritic inhibition, the firing rate decreased out to a lag of 125 ms, and then returned to normal at 500 ms. To better understand these effects, we calculated the cross-correlation between dendritic spikes and action potentials. Increasing the lag decreased the coordination between Ca²⁺ and somatic spikes (Fig. 4H). For distal dendritic inhibition, this decrease was accompanied by a broadening of their temporal relationship, while a sharp temporal relationship was maintained in the perisomatic case. A similar pattern was observed with NMDA spikes on apical branches (Fig. 4I). Basal branches, on the other hand, were relatively unaffected.

In summary, distal and perisomatic inhibition modulate the coupling between synaptic drive and spiking through distinct processes. And while for both cases, the normal 4 ms E/I lag produced maximal conversion of dendritic spiking events into action potentials, extending this lag exerted distinct effects.

Effect of beta and gamma rhythmic inhibition on dendritic integration

The mechanisms modulating neuronal responsiveness during tonic inhibition of somatic and dendritic compartments may extend to rhythmic inhibition. Thus, we emulated beta and gamma rhythmic input (Fig. 5A,F). Depths of modulation were set to similarly entrain action potentials (Fig. 5B,G), and were comparable to spontaneous and optogenetically-induced gamma and beta bursts seen in vivo [49–52]. Beta rhythmic inhibition was modeled as a 16 Hz sinusoidally modulated rate (20% depth) of the Poisson processes driving inhibitory synapses. Gamma rhythmic inhibition was a 64 Hz sinusoidal modulation (40%) of inhibitory synapses. For both cases, excitatory synaptic drive was a stable Poisson process.

Figure 5:

Initially, we delivered beta rhythmic inhibition to the distal dendritic compartments, and gamma to the perisomatic. Even though both rhythms produced similar depths of modulation of somatic action potentials, the underlying causes were distinct. The phase of beta modulated the occurrence of Ca²⁺, NMDA, and Na⁺ spikes, with each showing an ∼75% depth of modulation with respect to their mean level (Fig. 5C-E). These changes were seen across the entire dendritic tree, spanning sites electrotonically close to and far from the soma. In addition, the impact on Na⁺ spikes was unexpected (Fig. 5E), since delivery of the same inhibition tonically had little effect. By contrast, gamma had virtually no effect on dendritic spikes (Fig. 5F-H). Its strongest impact was on Na⁺ spikes in basal dendritic compartments that were electronically close to the soma.

To uncover the modulation of action potentials by gamma, we turned our attention to somatic action potential initiation (Fig. 6). For these analyses we divided the inhibitory rhythm into two phases. During the peak phase, inhibition was greater than its mean rate, while during the trough phase inhibition was lower (see Fig. 5A,F). We found that the somatic action potential voltage threshold shifted lower during the ‘trough’ phase of gamma, when inhibition was at its weakest (Fig. 6A1), and without any change in the mean membrane voltage (Fig. 6A2). This is consistent with gamma phase modulating the shunting of voltage-gated Na⁺ channel currents, which occurs when GABAergic synapses co-locate with the channels mediating action potential initiation [53]. We also observed changes in action potential initiation to beta rhythmic inhibition, but through a different mechanism. During the ‘peak’ phase of beta, when inhibition was maximal, the threshold for evoking an action potential increased, which may reflect an ‘off-path’ shunting of excitatory current away from the soma and towards the dendrites (Fig. 6B1; [54]). Additionally, there was a decrease in membrane voltage during the peak phase, which may correspond to decreased excitation arising from the suppression of dendritic spikes (Fig. 6B2).

Figure 6:

Accompanying these effects were changes in the rate of dendritic spikes compared with the Poisson inhibition case. Beta increased the rate of Na⁺ (+16.6%) and Ca²⁺ spikes (+15.1%) but decreased the rate of NMDA spikes (-10.7%). Gamma caused no change in NMDA (+0.2%) and a weak increase in Na⁺ spikes (+6.4%), but a robust gain in Ca²⁺ spikes (+37.9%).

We next switched the portions of the pyramidal neuron targeted by the beta and gamma rhythms. Beta rhythmic inhibition was delivered perisomatically and gamma rhythmic inhibition to distal dendrites. While phase modulation of firing rate was maintained with both rhythms, the overall level of spiking was dramatically reduced (Fig. 5S1A,E). Neither rhythm provided effective modulation of Ca²⁺ or NMDA spikes (Fig. 5S1B,C,F,G). It is likely that the slow time scale of Ca²⁺ and NMDA spikes, ∼50 ms, is not optimal for the fast periodicity of the gamma rhythm, which cycles every ∼15 ms. In agreement with this, Na⁺ spikes, which last less than 1 ms, did show modulation by gamma rhythms delivered to the distal dendrites (Fig. 5S1D,H).

Swapping the location of beta and gamma synapses also altered their effects on somatic excitability. Gamma rhythmic inhibition on the dendrites had minimal or no impact on action potential threshold, but did shift the somatic membrane potential more negative (Fig. 6S1A). This hyperpolarization was not dependent on gamma phase, and likely reflected an overall decrease in the rate of dendritic spikes that could supply excitatory drive to the soma (NMDA: -26.1%, Na⁺:

-12.2%, and Ca²⁺: -39.9% compared with the Poisson inhibition case). By contrast, delivering beta rhythmic inhibition to the soma raised the action potential threshold and hyperpolarized the membrane potential during the peak phase (Fig. 6S1B). It also reduced the incidence of dendritic spikes (NMDA: -17.1%, Na⁺: -11.5%, and Ca²⁺: -7.6%).

Putting all this together, the effectiveness of rhythmic inhibition depends on where it impinges upon the neuron. Beta rhythms targeting the distal dendrites modulate the incidence of dendritic spikes in a phase-dependent manner, while gamma rhythms delivered perisomatically phase modulate somatic excitability. Swapping the locations of these rhythms diminishes these effects and lowers overall excitability. Because of this, the subsequent simulations will be restricted to the case of beta on distal dendrites and gamma on perisomatic compartments.

Modulation of dendritic spikes during oscillatory bursts

In vivo, beta and gamma rhythms occur as spontaneous bursts lasting less than a few hundred milliseconds. Since the previous analyses relied on tonically delivered rhythms, we verified that similar effects were observed with oscillatory bursts. Gamma and beta bursts were delivered to the same model with mean depth of modulation similar to the tonic case (Fig. 7A,F; see Methods for details).

Figure 7:

Under the burst regime, both rhythms mirrored their behavioral effects. Gamma bursts entrained spiking, with entrainment strongest during the middle of the burst (Fig. 7B). As with the tonically imposed rhythm, there was none or minimal modulation of Ca²⁺ (Fig. 7C), NMDA (Fig. 7D), and Na⁺ spikes (Fig. 7E). Beta rhythms entrained somatic action potentials (Fig. 7G), Ca²⁺ spikes (Fig. 7H), NMDA (Fig. 7I), and Na⁺ spikes (Fig. 7J). These modulations were evident within the first few cycles of a burst, suggesting that they did not require a buildup or evolving entrainment of an underlying process.

Effect of beta and gamma rhythms on responding to clustered synaptic drive

The results so far suggest that beta and gamma rhythms modulate synaptic integration through different mechanisms depending on their phase and the location of the synapses on the dendritic tree. To examine this further, we added patches of concentrated excitatory synaptic inputs onto either the apical or basal dendrites (Fig. 8A), with densities similar to functional clusters in vivo [55; 56]. Coactivated inputs were simulated by driving each synapse with a jittered (2 ms) Poisson process. We ran six separate simulations, either with the apical or basal clusters engaged, and under conditions of Poisson, beta, or gamma rhythmic inhibition.

Figure 8:

To capture how beta and gamma influence synaptic integration we measured the cross-correlation between presynaptic activations at the clustered input and somatic spiking. Separate cross-correlograms were calculated depending on whether the presynaptic spikes occurred during the peak or trough phase of the inhibitory rhythm. This necessarily introduced spurious periodicities into the cross-correlogram that were compensated (see Methods for details). Relative to the arhythmic Poisson inhibition case, beta rhythms enhanced the transmission of apical inputs when inhibition was low (trough phase) and suppressed them when inhibition was high (peak phase, Fig. 8B, top left). Basal inputs were either unaffected or moderately suppressed during the trough and suppressed during the peak (Fig. 8B, bottom left). The opposite was the case for gamma. It barely affected or moderately suppressed apical inputs (Fig. 8B top right), while basal inputs were enhanced during the trough and suppressed during the peak (Fig. 8B bottom right).

Summarizing these results (Fig. 8C), we found that somatic spiking driven by clustered basal synapses was bidirectionally modulated by gamma rhythms and suppressed by beta. On the other hand, spiking driven by apical clusters was bidirectionally modulated by the beta rhythm and suppressed by gamma. Thus, both rhythms regulate the sensitivity of pyramidal neurons to afferents throughout the dendritic tree, but in a counterposed location-dependent manner.

Discussion

Arising from multiple interneuron subtypes, inhibition sculpts pyramidal neuron activity by acting at different membrane regions and distinct rhythmic frequencies. Little was known about how these factors interact with the complexity of dendritic integration. To address this, we characterized the interaction between the location and rhythmicity of inhibition on integration in a morphologically and biophysically detailed layer 5 pyramidal neuron model with active dendrites. We found that distal dendritic inhibition modulated the occurrence of dendritic spikes, while perisomatic inhibition altered action potential generation. This translated into location-specific differences in the effectiveness of inhibitory rhythms. Beta rhythmic inhibition entrained dendritic spikes, focusing them into the phase when inhibition was at a minimum, but only when delivered to the distal dendrites. In contrast, gamma modulated the threshold for action potential initiation, but only when delivered perisomatically. As a likely result, beta oscillations bidirectionally controlled transmission in apical dendrites and suppressed those onto basal dendrites, while gamma oscillations did the opposite.

This indicates that the spatial targeting and timing of inhibition go hand in hand. Their alignment is thus fundamental to their distinct effects on neuronal integration. And in turn, this affects function at multiple levels, from synaptic responsiveness to neural coding, and microcircuit operation.

Synaptic responsiveness

Oscillatory rhythms potentially affect synaptic integration through several mechanisms. Early theorizing held that oscillations synchronize populations of excitatory neurons such that their efferent synapses benefit from spatiotemporal summation in eliciting action potentials [57]. Later it was proposed that if neurons in two regions exhibited coherent variations in their somatic membrane potentials, bringing them closer to or farther from action potential threshold, they could preferentially exchange activities [58]. Modeling and in vitro studies have shown that the phase of a sinusoidal current injected at the soma modulates the timing of evoked action potentials [59] or their probability of occurrence [60]. Similar gamma phase dependence has been found in vivo [61].

Rhythmic inhibition may also periodically modulate the sensitivity of the neuron to synaptic input. Inhibitory synapses achieve this by lowering the membrane resistance i.e., shunting [62]. This makes it harder for an excitatory synaptic current to drive the membrane voltage toward the action potential threshold. Indeed, optogenetically stimulating PV interneurons at gamma frequencies modulates the responsiveness of cortical neurons to sensory inputs in a phase-dependent manner [63]. However, normally PV interneurons are driven by local principal cells, so the inhibition they impose reflects the aggregate activity in the local network (leading to E/I balance, [64]). Model networks wired in this manner emit spontaneous gamma bursts that produce phase-dependent modulation of the relationship between EPSP amplitude and spiking probability, with a strong positive relationship during the trough of gamma, and a weak one during the peak [34]. Consistent with these results, putative monosynaptically connected pairs of single units recorded in vivo show an increase in spike transmission probability during the trough of gamma [65].

To summarize the present possibilities, oscillations could modulate a neuron’s responsiveness by (1) synchronizing the synaptic inputs impinging upon it, (2) modulate how close the membrane potential is to the action potential threshold, or (3) vary its sensitivity to synaptic inputs via shunting inhibition.

Our results fill out this picture in two ways. First, we add a fourth mechanism, which is (4) the modulation of dendritic spiking events. Second, we found that beta and gamma oscillations differentially engaged these mechanisms. Beta oscillations primarily caused 2 via 4, and affected 3. Gamma, on the other hand, operated primarily through 3. It should be noted that to isolate these effects, our model had Poisson excitatory synaptic inputs, which precluded mechanism 1, afferent synchrony. However, it is likely the degree of synchrony would impact both rhythms. Supporting this, in vivo whole-cell recordings from neurons in mouse visual cortex exhibit brief membrane depolarizations in synchrony with gamma [66].

Relevance to coding

Given that beta and gamma rhythms influence spiking via distinct mechanisms, they may also differentially impact coding. During gamma bursts in V1, neurons show a stronger modulation of their firing rate by the contrast of grating stimuli [67]. Gamma phase also modulates orientation selectivity and correlated variability between neurons [68]. This enhancement is modulated by attention [69]. Our results indicate that these effects reflect shifts in somatic sensitivity to excitatory synaptic inputs.

Much less is known about the influence of beta oscillations on coding. In the visual cortex of mice, activation of SOM interneurons increased with the visual stimulus size and homogeneity ([21], note that this study referred to the rhythm as low gamma). In a different study, using a sequence working memory task in primates, the order of presented items corresponded to their phase in the beta rhythm [70]. However, this does not indicate how coding varies with beta phase. Considering our finding that beta synchronized dendritic spikes, it may operate by enabling the summation of normally asynchronous afferents on distal dendrites, facilitating their ability to drive somatic spiking.

Interaction with microcircuitry

Pyramidal neurons are embedded in a cortical column where they are sparsely connected amongst themselves and densely connected with local interneurons. The dense connectivity with interneurons is crucial for the generation of local beta and gamma rhythms [21; 22]. How might this circuitry interact with the differential regulation of dendritic integration by beta and gamma?

We found that beta rhythms preferentially modulated a pyramidal neuron’s response to inputs on its apical dendrites. The apical tuft in the superficial layers of cortex primarily receives long-range inputs from regions higher up in the cortical hierarchy [71; 72]. This makes beta rhythms ideally positioned to modulate these top-down/feedback signals. Since the rhythm likely derives from local activation of SOM interneurons, it is worth considering what situations would activate those and how they might relate to the functioning of apical dendrites. SOM interneurons receive facilitating synapses from local pyramidal neurons [73; 74], which would make them especially sensitive to bursts of action potentials arising from Ca²⁺ spikes at the apical nexus [18]. Since Ca²⁺ spikes are driven by synaptic activation in the apical tuft, it is likely that beta rhythms regulate the generation of action potential bursts arising from long-range inputs.

Our results indicated that gamma rhythms mainly affect processing in basal dendrites. These are preferentially targeted by layer 5 pyramidal neurons, which sparsely interconnect with projection-type specificity [75–77]. A confluence of evidence implicates PV interneurons in the production of high-frequency inhibitory rhythms such as gamma [78; 79]. Given that PV interneurons diffusely interconnect with local pyramidal neurons [80; 81], and the synapses they receive from them are depressing [73; 74], they likely regulate abrupt increases in local ensemble activity. Interactions within ensembles that precede activating the PV population will be boosted, while those following will be attenuated.

Interneuron specializations and rhythmic timescales

Lastly, it is worth considering why beta and gamma rhythms are mediated by different types of interneurons. One view is that the relative timing of their rhythms arises from the electrotonic timescales of their feedback inhibition; SOM cells target distal dendrites and thus produce a slower inhibitory feedback signal, versus the rapid perisomatic inhibition of PV cells. And when it comes to explaining the emergence of PING rhythms at different timescales, this story is sufficient, but incomplete.

We found that swapping the cellular compartments receiving beta and gamma rhythmic inhibition impaired their effectiveness and lowered overall excitability. Thus, moving beyond mere timing, our results suggest that the timescales of these rhythms align with the specialized impacts of SOM and PV interneurons on neuronal integration.

Conclusion

For the most part, the study of beta and gamma rhythms has focused on their correlation with task-related events and entrainment of neuronal firing. Fundamentally, these effects derive from their influence on synaptic integration, which has multiple stages, from modulating the post-synaptic response, inducing a dendritic spike, and triggering an action potential. Since beta and gamma rhythms differentially impacted these processes, this invites a reappraisal of their role in coding and communication, along with unveiling a new hypothesis space for understanding their function.

Funding

This work was supported by the following grants: NIH NINDS 5R01NS123396 to DBH, and NIH MH122023 and NSF OAC-1730655 to SSN.

Note

This reviewed preprint has been updated to correct the figure images, which were switched around in the original publication, and to add funding information.

Methods

We adapted a previously published model of a L5 pyramidal neuron [1; 2] to the Brain Modeling Toolkit [3] format to facilitate reproducibility. This open-source software utilized NEURON 7.7 [4] with an integration time step of 0.1 ms for simulating the membrane potential. We added synapses and designed inputs in accordance with the literature to study synaptic integration in vivo. Details of the passive membrane properties and active conductances are only described briefly (full details can be found in [2]). Our focus here is on describing the enhancements we incorporated to mimic in vivo features in the model neuron.

L5 pyramidal cell compartmental model

Briefly, the cell was built on a morphological reconstruction of a L5 pyramidal tract neuron and included 10 active Hodgkin-Huxley conductances with kinetics taken from the literature [2]. These channels were a fast-inactivating Na⁺ current (I_NaT), persistent Na⁺ current (I_NaP), non-specific cation current (I_h), muscarinic K⁺ current (I_m), slow inactivating K⁺ current (I_Kp), fast inactivating K⁺ current (I_Kt), fast non-inactivating K⁺ current (I_Kv3.1), high voltage-activated Ca²⁺ current (I_{Ca_HVA}), low voltage-activated Ca²⁺ current (I_{Ca_LVA}), and Ca²⁺ activated K⁺ current (I_SK). The somatic compartments contained all of these channels except I_m. Both the apical and basal dendrite compartments contained I_{Ca_LVA,} I_{Ca_HVA,} I_SK, I_Kv3.1, I_NaT, I_m, and I_h. The membrane capacitance was set to 1 μF/cm² for the soma and axon and 2 μF/cm² for the basal and apical dendrites to correct for dendritic spine area. Axial resistance was set to 100 W-cm and leak reversal potential was -90 mV.

Divergence, probability of transmission, and functional clustering

Experimental studies indicate that each presynaptic excitatory neuron forms 2-8 synapses with its pyramidal partners, with a release probability between 0.16 and 0.9 [5]. Further, the 2-8 synaptic contacts occur over a span of 20-100 μm [6]). We constrained our model to have these properties. In addition, since synapses may cluster on dendritic branches such that functional presynaptic cell assemblies activate the same branch within 5-10 μm [7; 8], we forced synapses coming from the same presynaptic cell assembly to cluster in this way. In summary, a single presynaptic input fiber could have 2-8 synapses, drawn from a uniform distribution with a release probability between 0.16 and 0.9, also drawn from a uniform distribution [5]. And the synapses from a functional group formed clusters within 5-10μm on the same dendritic branch over a span of 100 μm.

Synaptic conductances

We randomly distributed excitatory synaptic conductances following previous modeling work [9]. Maximum conductance values followed a log-normal distribution with a mean of 0.2 nS and standard deviation of 0.345 nS for excitatory synaptic conductances. Inhibitory synaptic conductances were fixed at 1 nS [1]. We validated these conductances against excitatory and inhibitory post-synaptic currents measured in vitro. The model’s values are compared with experimental values in Table 1.

Table 1.

Table 2.

Synapse Density

Synapse density was based on a comprehensive report [10] that not only quantified the number of spines but also counted the synapses by identifying vesicle clouds and postsynaptic density using serial block-face electron microscopy. The authors found the apical tuft synapse density per micron as 2.16 ± 0.16 for excitatory and 0.22 ± 0.02 for inhibitory synapses (mean ± standard error of the mean). Other studies that only quantified spine density may undercount the true number of synapses by about 3 times because one can only count lateral spines, not those in the z-axis [11]. With that adjustment, the other values reported – 0.47 ± 0.01 spines/µm [12] and 0.53 ± 0.05 spines/µm [13]– are roughly in alignment with that in [10]. We used the values from [10] for excitatory synapse density because we found the higher value for excitatory synapse density was needed to achieve a reasonable baseline firing rate.

Morphological measurements (length of dendritic cable in µm) combined with these synapse density estimates yield the following estimates for total number of synapses: Excitatory: 16070 (2.16 syns/µm * 7440 µm) on apical dendrites and 10042 (2.16 syns/µm * 4649 µm) on basal dendrites. Inhibitory: 150 on soma [1; 14; 15] + (0.22 syns/µm *483 mm) = 256 perisomatic, 1637 (0.22 syns/µm * 7440 µm) on apical dendrites, 1023 (0.22 syns/µm * 4649 µm) on basal dendrites.

Synaptic kinetics and short-term plasticity

All excitatory transmission was mediated by AMPA/NMDA receptors and inhibitory transmission by GABA_A receptors as in our prior models. The rise/decay constants for the excitatory connections were as follows: Rise/decay time constants of 0.6/6.9 ms for AMPA, and 3.7/125 ms for NMDA synapses, and a reversal potential of 0 for both. And for the inhibitory GABA-ergic connections, the rise/decay time constants were 0.5 and 6.8 ms, respectively, with a reversal potential of -75 mV. We implemented short-term depression in both the excitatory and inhibitory synapses [16] via an exponential decay of synaptic weight over short time scales [17]. Following [16] we calculated a short-term plasticity (STP) train-induction value by applying 8 pulses at 50 Hz to the synapse while in voltage clamp and measured the resulting postsynaptic current (PSC). The 6^th, 7^th, and 8^th PSCs were averaged, and the 1^st PSC was subtracted from the result. This value was divided by the max PSC current to yield the STP train-induction value. Excitatory synapses to L5 pyramidal neurons had an STP train-induction value of -0.34. Inhibitory synapses arising from PV+ basket interneurons making perisomatic contact had a value of -0.60 and inhibitory synapses arising from SOM+ interneurons making dendritic contact had a value of -0.10. The equations to model the synaptic currents and short-term depression are provided in the next section.

Modeling of synaptic currents

The corresponding synaptic currents were modeled by dual exponential functions [18; 19]. All excitatory transmission was mediated by AMPA/NMDA receptors and inhibitory transmission by GABA_A receptors. The corresponding synaptic currents were modeled by dual exponential functions [18; 19], as shown in eqns. 1-3,

where V is the membrane potential (mV) of the compartment (dendrite or soma) where the synapse is located, I is the current injected into the compartment (nA), G is the synaptic conductance (µS), 𝑤𝑤 is the synaptic weight (unitless), and E is the reversal potential of the synapse (mV). g_x,max is the maximal conductance (µS), D implements short-term plasticity as defined below, and r_x determines the synaptic current rise and decay time constants based on the terms αTmax and β [18]. The voltage-dependent variable s(V) which implements the Mg²⁺ block was defined as: s(V) = [1 + 0.33 exp(-0.06 V)]^-1 [20]. The terms ON_NMDA and ON_AMPA are set to 1 if the corresponding receptor is open, else to 0. To calculate the factor D for depression: τ_D*dD/dt=1-D and D constrained to be ≤ 1. After each stimulus, D was multiplied by a constant d (≤ 1) representing the amount of depression per presynaptic action potential and updated as D→D*d. Between stimuli, D recovered exponentially back toward 1. We modelled depression using two factors d1 and d2 with d1 being fast and d2 being slow subtypes, and d=d1*d2 was constrained to be ≥ 1.

Synapse Distribution and Extrinsic inputs

We distributed 26,112 excitatory synapses over the cell. With an average of 5 synapses per pre- synaptic neuron, this yielded 5,222 presynaptic nodes. To represent pre-synaptic cell assemblies, we defined “functional groups” of presynaptic nodes. We put 100 nodes in each functional group, yielding 52 functional groups. We found that varying the number of neurons in the cell assembly did not qualitatively affect our results. To have exactly 100 nodes per functional group we rounded the total number of nodes down to 5,200. Inhibitory inputs were not defined by functional groups, but followed the literature cited in Table 1, with more than 1 synapse per pre-synaptic node. For example, with 1637 inhibitory synapses on the apical dendrites and an average of 12 synapses per pre-synaptic node, there were about 136 inhibitory nodes that contacted the apical dendrites. The numbers and distributions of model synapses along the dendrites are shown in Table 1.

The extrinsic inputs to the model were represented by spike trains arriving at the corresponding synapses along the dendritic tree. The L5 pyramidal model is spatially divided into the following sub-divisions: somatic, perisomatic, and proximal and distal components of the apical and basal trees. We used standard definitions for these sub-divisions except for perisomatic which is typically defined as the dendritic region within 100 μm of the soma [21].

Tonic inputs (Figs. 1-3)

The firing rate of each excitatory spike train was randomly sampled from an experimentally recorded distribution (Headley, personal communication). To reproduce the 1/f spectrum that is widely found in neuronal spike trains, the 100 spike trains in each functional group were modulated by a unique normalized (from 0.5 to 1.5) firing rate time series that had a 1/f spectrum (more slow modulation of firing than fast) that was created by filtering white noise using standard techniques. The normalized firing rate time series was then multiplied with the mean firing rate of each of the 100 pre-synaptic nodes. Spike trains are then generated by sampling from a Poisson distribution using the firing rate time series of each presynaptic node of the functional group. In this way, spike trains that received the same modulation trace would be statistically correlated, representing functional presynaptic cell assemblies. Since feed-forward inhibition follows excitation in time, we averaged the firing rate time series used to design the excitatory spike trains, shifted it forward in time by 4 milliseconds, and used this to modulate the firing rate of the inhibitory spike trains.

Enhanced tonic inhibition (Fig. 4)

Inhibitory synapses were divided into two groups, those targeting perisomatic region that originate mainly from PV+ interneurons, and those on more distal dendritic regions arising mainly from SOM+ interneurons [22; 23]. Perisomatic synapses were defined as those within 100 µm of the soma with the properties specified in Table 1 under “Inhibitory (Perisomatic and somatic).” Distal inhibitory synapses were greater than 100 µm from the soma and had their properties specified in sections “Inhibitory (Basal)” and “Inhibitory (Apical)” of the table To investigate if these inhibitory sources differentially regulated action potential initiation and dendritic spiking, we doubled the rate of presynaptic events to one of these inhibitory populations.

Investigation of E/I lag (Fig. 4F-I)

We used the same scheme to design the extrinsic inputs as in Figs. 1-3. However, to investigate the effect of E/I lag, we increased the variance of the overall excitatory firing rate. So, while the mean rate stayed the same, we raised the modulatory trace to the fourth power so that large deviations from the mean were exaggerated. This excitatory firing rate trace was used to modulate the inhibitory firing rates (same as for Figs. 1-3). The perisomatic and distal traces were shifted independently. One was fixed at 4 ms while the other was shifted by 4, 125, 250, or 500 ms.

Rhythmic inhibition (Figs. 5 & 6)

Gamma rhythms depend on PV+ interneurons, while beta rhythms arise from SOM+ interneurons [24]. Thus, to emulate an inhibitory rhythm we rhythmically modulated inhibitory synapses in the respective dendritic region (see description above in “Enhanced tonic inhibition”) at the appropriate frequency. Beta rhythms were emulated by modulating distal inhibitory synapses at 16 Hz with a depth of 20%. Gamma rhythms were emulated by modulating perisomatic synapses at 64 Hz with a depth of 40%. The different depths of modulation ensured that these rhythms produced similar entrainment of action potentials.

To generate inhibitory spike trains that produce rhythmic inhibition during oscillations, we modulated their collective firing rate with a sine wave:

where , 𝑓 is the frequency, 𝑡 is a vector representing time, 𝜙 is the phase, 𝑜𝑓𝑓 is the offset firing rate of the spike train being modulated, and 0 < 𝑑 < 1 is the depth of modulation which represents the amplitude of the sine wave relative to 𝑜𝑓𝑓. To generate the spike train, a random vector x[t] was created with values uniformly distributed between 0 and 1. A spike was generated if x[𝑡] ≤ 𝑟[𝑡]𝑑𝑡 where dt in our case was 1/1000.

Bursty inhibition (Fig. 7)

We applied an envelope to the rhythmic inhibitory modulation trace to simulate the bursty temporal nature of rhythms in vivo. The underlying sine wave was produced in the same way as in the previous section. The envelope was a Gaussian kernel with a standard deviation of two oscillatory cycles.

Clustered synaptic input (Fig. 8): For this experiment we kept everything the same as in the rhythmic inhibition simulations, but introduced additional synapses driven by 40 pre-synaptic spike trains. These spike trains were created by choosing one Poisson process and jittering the spike times by ±2 ms. Each spike train drove an average of 4.6 synapses for a total of 187 synapses. These were placed either on the apical dendrites near the nexus or the basal dendrites near the soma. For both nexus and basal clustered synapses, we targeted 140 μm of dendritic cable length with 187 synapses yielding a synaptic density of 1.3 synapses per micron. For the nexus case we targeted the 140 µm region starting at the bifurcation going toward the apical dendrites. For the basal case, we targeted 3 basal dendrites that were directly connected to the soma each one had a length of 40-50 µm adding up to 140 μm total (Figure 7A). The synaptic weights were drawn from the distribution in Table 1.

Detection of action potentials

Action potentials were detected by the NEURON software as being any time point when the somatic membrane potential crossed -10 mV with a positive first derivative.

Detection of dendritic spikes

Na⁺ and NMDA dendritic spike events were detected following the definitions developed in [9]. We developed our own detection algorithm for Ca²⁺ spikes since they were not previously studied. Na⁺ spikes were detected if the voltage-gated sodium conductance density for a compartment exceeded 0.3 mS/cm². If sodium spikes were within 2 ms of a somatic action potential, they were assumed to be a consequence of a backpropagating action potential and were not counted. An example of a Na⁺ spike can be seen in Fig. 1C1.

NMDA spikes were detected by two criteria – voltage and current. The voltage criteria were met if the membrane potential rose above -40 mV and remained there for 26 ms or more. The current criteria was met if the NMDA current exceeded 130% of the current at the -40 mV crossing and continued to be met until the current fell to 115% of that value. Whenever both the voltage and current criteria were met, an NMDA spike was recorded. An example of an NMDA spike is shown in Fig. 1C2.

Ca²⁺ spikes were detected in the same way as NMDA spikes – using a voltage and current threshold. The only difference was instead of NMDA current for the current criteria we used the sum of HVA and LVA calcium currents. An example is shown in Fig. 1C3.

Determination of electrotonic distance

We calculated electrotonic distance using NEURON’s built-in Impedance method. The method calculates signal attenuation between two compartments at a given frequency. We used 20 Hz because the wavelength is on a time scale comparable to the membrane time constant (10s of milliseconds). However, changing this frequency 10 Hz in either direction did not qualitatively affect our results.

Spike-triggered average (STA) of dendritic spikes (Figs. 2, 3)

For each dendritic compartment, we created a binary time series with 2 ms bins that was true during periods where a dendritic spike was present (based on their detected start and stop times). This was done separately for each type of dendritic spike (Na⁺, NMDA, and Ca²⁺ spikes). Each then had its spike-triggered average taken around the times of somatic action potentials and their percent change from the mean across the entire simulation calculated. The STAs were then grouped together based on their electrotonic distance percentile (10% bins) from the soma and type of dendrite (apical or basal), and the median taken for each time lag.

Generation of F-I curves (Fig. 4)

The model with Poisson excitation and inhibition was used. In the perisomatic condition, we doubled the rate of inhibition to the somatic and proximal dendritic compartments (<100µm from soma). In the distal case, we doubled the rate of inhibition to the distal dendritic compartments (>100µm from soma). In each case, we injected 20 s current steps ranging from -1 nA to 2.4 nA in 200 pA intervals.

Phase histogram of dendritic spikes (Fig. 5)

Entrainment of dendritic spikes to inhibitory rhythms was calculated by first converting dendritic spikes to a binary time series as was done for the STA calculation. The inhibitory rhythm was then binned by phase, with bin widths of π/4, yielding 8 phase bins. The mean of the dendritic spike binary series was taken for each phase bin. These were converted to percent change from mean and their median taken across compartments by their electronic distance and dendritic type (the same as in the STA calculation).

Determination of AP threshold and somatic membrane voltage distribution (Fig. 6)

AP threshold was taken to be the voltage 1 ms prior to the peak of the spike. Voltage distributions were calculated using 1 ms bins.

Phase histogram of dendritic spikes (Fig. 7)

For calculating entrainment of dendritic spikes to the rhythm, we followed the same procedure as detailed above for Fig. 5.

Cross-correlation (CC) between presynaptic spikes and action potentials (Fig. 8)

We measured the effect of inhibitory rhythms on the coupling between concentrated synaptic inputs and action potentials using a phase-stratified CC. Presynaptic spike times that drove the concentrated synaptic inputs were grouped together and separated based on whether they occurred during the positive or negative phase of the imposed inhibitory rhythm. Their CC with somatic action potentials was then calculated using 1 ms bins. Since both the presynaptic and action potential times had inherent periodicities, these had to be corrected for [25]. This was performed by converting the CC and both autocorrelations to the frequency domain using the fast Fourier transform, setting the coefficients at the inhibitory modulation rhythm to 0, performing element-wise division of the CC by the square-root of the autocorrelations multiplied together, and then converting back to the time domain using the inverse fast Fourier transform on that result. To measure the strength of the cross-correlation peak, we identified the positive peak of the CC nearest to zero time lag, and then integrated over its area with the edges defined by when the CC crossed from positive to negative values.