Research ArticleResearch Article: New Research, Cognition and Behavior

γ And β Band Oscillation in Working Memory Given Sequential or Concurrent Multiple Items: A Spiking Network Model

Shukuo Zhao, Jinpu Zhou, Yongwen Zhang and Da-Hui Wang
eNeuro 30 October 2023, 10 (11) ENEURO.0373-22.2023; https://doi.org/10.1523/ENEURO.0373-22.2023

Abstract

Working memory (WM) can maintain sequential and concurrent information, and the load enhances the γ band oscillation during the delay period. To provide a unified account for these phenomena in working memory, we investigated a continuous network model consisting of pyramidal cells, high-threshold fast-spiking interneurons (FS), and low-threshold nonfast-spiking interneurons (nFS) for working memory of sequential and concurrent directional cues. Our model exhibits the γ (30–100 Hz) and β (10–30 Hz) band oscillation during the retention of both concurrent cues and sequential cues. We found that the β oscillation results from the interaction between pyramidal cells and nFS, whereas the γ oscillation emerges from the interaction between pyramidal cells and FS because of the strong excitation elicited by cue presentation, shedding light on the mechanism underlying the enhancement of γ power in many cognitive executions.

Significance Statement

We constructed a spiking network to perform working memory (WM) tasks with sequentially or concurrently presented items. The model exhibits the coexistence of β (10–30 Hz) and γ (30–100 Hz) band oscillations during the delay period. We found that γ and β band oscillations recruit separate neural circuits. The low-threshold nonfast-spiking (nFS) neurons are involved in the β band oscillation, whereas the high-threshold fast-spiking (FS) neurons are involved in the γ band oscillation. Our results shed light on the well-known phenomenon that cognitive tasks enhance γ band oscillations.

Introduction

Working memory (WM), the ability to actively maintain and manipulate information in the absence of stimuli, plays a crucial role in cognitive function and executive control of behavior (Baddeley, 2003). The information maintained in WM can enter the brain concurrently or sequentially, meaning that WM cannot only hold multiple pieces of information arriving simultaneously as in visuospatial WM tasks (Bays and Husain, 2008; Zhang and Luck, 2008; Gorgoraptis et al., 2011), but also information presented sequentially, as in visual (Gorgoraptis et al., 2011) or speech processing (Cowan, 2010). The neural implementation of information maintenance is not well understood. The oscillatory model proposes that the information of one item is represented in WM by the reactivation of neurons in the γ cycle within nested γ-theta oscillations, mediated by a slow after-depolarization (ADP) with a time constant that should match the theta oscillation (Lisman and Idiart, 1995; Jensen and Lisman, 1998). The dynamic model posits that memory is maintained by item-specific patterns of synaptic plasticity and that neurons exhibit a nonstationary and short-lived attractor activity, in which only one memory representation can be active at a time, but successive reactivations of neuronal pools memorize different items (Mongillo et al., 2008; Lundqvist et al., 2016, 2018; Mi et al., 2017). The persistent activity model asserts that persistent activity of neurons induced by N-methyl-D-aspartate receptor (NMDAR)-mediated recurrent synaptic current encodes the corresponding items during the delay period (Compte et al., 2000; Edin et al., 2009; Wei et al., 2012). The ability of WM to maintain concurrent and sequential information challenges these models. On the one hand, the oscillatory model and dynamic models can use oscillatory activity to maintain sequential information, but it is difficult for them to represent concurrent information (Constantinidis et al., 2018); on the other hand, the persistent models can manipulate concurrent information but cannot address sequential stimuli and often focus on the persistent activity without oscillation. Experiments have shown that γ band oscillations are involved in WM maintenance, in particular, that WM load enhances the γ band oscillations but suppresses the alpha band oscillations. Here, we aim to construct a biophysically plausible network model to implement the storage of concurrent and sequential information and to investigate how WM load alters oscillations in the network.

The biological brain consists of microcircuits with three types of neurons that play key roles in generating oscillations: pyramidal cells, low-threshold nonfast-spiking (nFS), and high-threshold neuronal FS (Chen et al., 2017; Veit et al., 2017; Domhof and Tiesinga, 2021; Hahn et al., 2022). The pyramidal cells are excitatory, while the FS and nFS are GABAergic inhibitory. Pyramidal cells synapse on nFS and FS. In turn, nFS neurons synapse back onto the dendrites of the pyramidal cell, while FS neurons synapse back onto the soma of the pyramidal cell, creating a dynamic feedback loop that regulates excitatory and inhibitory activities. FS and nFS also inhibit each other, resulting in a competitive relationship. At rest, the nFS is more active, inhibiting the FS and forming a β oscillation with the pyramidal cells. When a perceptual stimulus is presented, the FS cells are activated, inhibiting the nFS cells and creating a γ oscillation with the pyramidal cells (Chen et al., 2017; Hahn et al., 2022). Given the importance of oscillations in higher cognitive functions, we aim to explore the role of oscillations and microcircuitry in working memory.

We proposed a spiking neural network to implement WM with concurrent and sequential directional information (Gorgoraptis et al., 2011). The network consists of two-compartment pyramidal cells, nFS, and FS. These cells and compartments are interconnected as in the biological brain. The cells are uniformly arranged in a ring according to their preferred direction (Compte et al., 2000). When activated, pyramidal cells activate pyramidal cells with similar preferences as well as nearby FS cells and nFS cells on the ring. FS and nFS inhibit neighboring neurons on the ring. We found that, regardless of whether the directions are presented concurrently or sequentially, the network elicits corresponding localized activities that persist throughout the delay period. The persistent activity suggests that the model can successfully maintain multiple stimuli presented simultaneously or sequentially in working memory. Furthermore, before cue presentation, the interaction between pyramidal cells and nFS dominates the activity of the network and causes a low-band oscillation (10–30 Hz); the cue presentation induces strong excitation and recruits FS into the network, and the interaction between pyramidal cells and FS dominates the activity of the network, enhancing the γ oscillation (35–100 Hz), which persists throughout the delay period.

Materials and Methods

Model architecture

Our model has 4096 excitatory pyramidal cells, 512 FS cells, and 512 nFS cells. We have two reasons to choose the number of neurons. First, the ratio of the number of excitatory neurons over that of inhibitory neurons is ∼4:1 (Braitenberg and Schutz, 1991). Thus, we chose 4906 pyramidal neurons and 1024 inhibitory interneurons (FS + nFS, 1024). Second, a network size of 2n neurons is desired for the fast Fourier transform (FFT) which recursively divides the input data into smaller subsets, and conducts Fourier transform computations on subsets. These neurons are evenly distributed in a ring configuration and are connected by AMPA, NMDA, and GABA synapses, forming an interconnected spiking network. According to our model, the strength of connections between neurons and the angular differences in their distribution across the ring configuration follows a Gaussian distribution. Neurons that are closer in angle on the ring have stronger connections. The footprint of the connections can be described as Equation 1: {W(θi−θj)=J− + (J+−J−)exp(−(θi−θj)22σ2)θi∈EW(θi−θj)=exp(−(θi−θj)24π2σ2)θi∈l. (1)

Js→s+=6.5 ,Js→fs+=35 ,Js→nfs+=70 . By normalizing the footprint 1360∫0360W(θi−θj)dθj=1 , we obtain J− in Equation 1. The standard deviations of the Gaussian distributions for connectivity footprint in our model were as follows:

- From pyramidal cells soma to pyramidal cells soma (σs→s ): 12.76°

- From pyramidal cells soma to FS cell (σs→fs ): 7.05°

- From pyramidal cells soma to nFS cell (σs→nfs ): 1.41°

- From FS cell to pyramidal cells soma (σfs→s ): 8.46°

- From FS cell to nFS cell (σfs→nfs ): 1.41°

- From nFS cell to pyramidal cells dendrite (σnfs→d ): 5.02°

The value of parameters Js→s+ and σs→s used in our study were adapted from Compte et al. (2000) and Wei et al. (2012). We further introduced structured excitatory projections from pyramidal cells to FS and nFS and inhibitory projections from FS to pyramidal soma and nFS cells or from nFS to pyramidal dendrite, referencing local circuit properties described by Wang et al. (2004) and Wang and Yang (2018). FS belongs to the large basket cell that projects to broader region, while nFS has narrow dendritic and axonal arbors and serves local inhibition. Generally, excitatory neurons project to a larger radius (around 200 μmphysiologically; as shown by Fitzpatrick et al., 1985) compared with inhibitory interneurons (around 100 μm; as demonstrated by Spiro et al., 1999) based on biological experiments. Nevertheless, the lateral innervation width of interneurons may be similar to that of excitatory neurons. It is because excitatory neurons project to a larger radius, activating inhibitory neurons within that radius to locally innervate inhibition (as assumed by Compte et al). The synaptic weights were calculated by the production of footprint and the maximum conductance as follows: gij=W(θi−θj)Gsyn. (2)

TheGsyn parameters for these connections are specified in Table 1 (Compte et al., 2000; Wei et al., 2012). To simplify the calculations, we have assumed that the remaining inhibitory connections in our model (e.g., inhibitory connection from nFS to FS or recurrent connection from FS to FS) follow a uniform distribution (a simple diagram is shown in Fig. 1).

View this table:
Table 1

Gsyn Parameters for soma, dendrite, FS, and nFS synaptic connections

Figure 1.
Figure 1.

Model structure and simulation protocol. A, The network consists of two-compartment pyramidal cells, fast-spiking (FS), and nonfast-spiking (nFS) interneurons. Pyramidal and inhibitory cells are uniformly placed on a ring and labeled by their preferred directions. B, Simulation protocol. In concurrent-cue tasks, an array of cues is presented to the network from 500 to 700 ms, followed by a delay period of 1300 ms (upper panel). In sequential-cue tasks, cues were presented one after another with an interval of 100 ms, and each cue lasted for 200 ms.

Neuron models

Pyramidal cells and inhibitory cells follow the leaky integrate-and-fire model (Tuckwell, 1988). The neurons generate spikes when their membrane potentials reach threshold (Vth=−50mV .) And their membrane potential is reset to Vreset=−60mV immediately after a spike for τs=2ms ; τd=2ms ; τfs= 1ms; τnfs= 1ms as refractory periods. The Pyramidal cells have two compartments (Wang, 1998): soma(s) and dendrite(d). The membrane potential of the somatic compartment (Vs) and dendrite compartment (Vd) of pyramidal cells obey Equations 3 and 4: CmdVsdt=−ICa−gsL(Vs−V0)−gc1(Vs−Vd)p1−Isyn, (3) CmdVddt=−gdL(Vd−V0)−gc1(Vd−Vs)p2−Isyn. (4)

V0=−70mV , Cm=0.5nF , gsL=0.025nS , gdL=0.025nS , gc1=0.25nS (Troyer and Miller, 1997; Wei et al., 2012). The parameters p1 and p2 characterize the difference of current intensity resulting from the identical current in soma and dendrite because of the difference in membrane surface area. Specifically, p1 represents the ratio of somatic area to the total neuronal area, and p2 represents the ratio of dendritic area to the total area (Wang, 1998). In the current model, p1 is set to 0.6, and p2 is set to 1-p1, which is 0.4. FS and nFS cells contain single compartments, and their membrane potentials obey Equations 5 and 6. Please note that when the neuron (soma) fires, the membrane potential of the dendrite component is not reset to −60 mV. The dendrite component generates its own spikes and resets the membrane potential after each spike. Cm_fsdVfsdt=−gfsL(Vfs−Vfs0)−Isyn (5) Cm_pvdVnfsdt=−ICa−gnfsL(Vnfs−Vnfs0)−Isyn. (6)

Cm_fs=0.2nF,Cm_nfs=0.8nF,gfsL=0.020nS,gnfsL=0.016nS,Vfs0=−86mV , Vnfs0=−76mV (Troyer and Miller, 1997; Wei, 2012). The term “Isyn” represents the total synaptic current. The ion current ICa in Equations 3 and 6 represents a high-threshold calcium current following Equations 7 and 8 (Wang, 1998): ICa=gCam∞2(V)(V−VCa) (7) m∞2(V)=11 + e−(V+20)/9. (8)

gCa−s=0.0015 nS, gCa−nfs=0.001 nS. V stands for membrane potential (equals Vs or Vnfs ).

Synapse models

Synaptic currents are mediated by AMPA, NMDA, and GABA transmissions. Three types of synaptic currents follow the equations: Ii,AMPA=(Vi−VAMPA)∑jgji,AMPAsj,AMPA (9) Ii,NMDA=(Vi−VNMDA)∑jgji,NMDAsj,NMDA1 + [Mg2+]exp(−0.026Vi/3.57) (10) Ii,GABA=(Vi−VGABA)∑jgji,GABAsj,GABA. (11)

The conductance “g” in our calculations primarily follows the Equation 2, VAMPA =VNMDA=0mV , VGABA=−70mV , [Mg2+] in Equation 10 equals 1 mm (Jahr and Stevens, 1990). The gating variable s for AMPA or GABA (sj,AMPA and sj,GABA in Eqs. 9, 11) follows the dynamic in Equation 12, and the dynamic of s (sj,NMDA ) for NMDA was modelled as Equations 13 and 14: ds(t)dt=−s(t)τs + ∑kδ(t−tk) (12) ds(t)dt=−s(t)τs,NMDA + αsx(t)(1−s(t)) (13) dx(t)dt=−x(t)τx,NMDA + ∑kδ(t−tk). (14)

tk is the spike sequence of the neuron. The rising time constant for NMDA (τx,NMDA ) is 2 ms, and the decay time constant τs is 2 ms, 10 ms, and 100 ms for AMPA, GABA, and NMDA, respectively (Wei, 2012).

Synapses from Pyrs to nFS mediated by NMDAR were facilitated in short term and follow the equation: dudt=U−uτf,NMDA + k(1−u)δ(t−tk). (15)

U is set to 0.4 in our simulation (Mongillo et al., 2008).

Simulation protocol

We simulated two types of WM tasks: concurrent-cue tasks (CCTs) and sequential-cue tasks (SCTs), similar to the experiments in (Gorgoraptis et al., 2011). In CCTs, cues are presented to neurons concurrently; each pyramidal cell indexed by θ received a current Iext(θ) during the concurrent cue presentation period from 500 to 700 ms (Fig. 1B). Iext(θ)=∑k=1nI02πσsexp(−(θ−θin,k)22σs2), (16)where θin,k is the kth direction and n is the number of directions in the cue array. I0=0.5nA,σs=2 degree.

In SCTs, pyramidal cells received multiple currents sequentially. Each current can be formulated as: Iext(θ)=I02πσsexp(−(θ−θin,k)22σs2). (17)

The kth current is applied at 500+300(k−1) to 700+300(k−1) ms (Fig. 1B) during the presentation of the kth cue (Wei et al., 2012).

Noise implementation

Noise is applied to pyramidal cells and inhibitory cells via AMPA synapses: Ii,noise=(Vi−VAMPA)gnoisesnoise−i (18)

gnoise−pyr=0.0039 nS,gnoise−pv=0.0019 nS,gnoise−cb=0.0022 nS, and snoise−i follows: ds(t)dt=−s(t)τs+η(t) (19)

τs = 2 ms, and η(t) is a Poisson process. p(η(t)=k)=e−λ*dt(λ * dt)kk!, k=0, 1... (20)where λ is the arrival rate of noise. We set λ = 1 kHz, which is equivalent to 1000 synaptic input at 1 Hz (Wei et al., 2012).

Starting parameters

To establish a more realistic and stable initial condition, we conducted a recording procedure to obtain a consistent value of multiple variables for each type of neurons. The following variables were used: the membrane potential (“Vm ”), the gating variables (“s” and “x”) for synapses, the time of the last firing activity (“LastTimeEachFired” minus 250 ms), and the AMPA synaptic current (“I_AMPA”) for the pyramidal soma (over a period of 10 ms), and these recorded values were then used as initial conditions for the algorithms. Initially, we performed a model run with constant initialization values of −51 mV for the membrane potential and 1e-34 for the “s” parameter without these recordings.

Numeric integration and code accessibility

The model was programmed in MATLAB code and the gating variables were integrated using a second-order Runge–Kutta (RK2) algorithm. xn+1=xn + 12dt(K1 + K2) (21) K1=f(xn,tn) (22) K2=f(xn + dtK1,tn + dt). (23)

The membrane potential is integrated with Euler scheme as Equations 3–6. This different choice in the integration method is because changes in gating variables during neuronal firing are more complex and usually nonlinear, whereas changes in the subthreshold membrane potential are approximately linear. We conducted a convergence test on two neurons with AMPA connections and provided a 2-mV constant stimulus to the presynaptic neuron for 2000 ms to examine the difference of spike time of postsynaptic neuron between Euler-RK2 and RK2 methods. We found that the spike time errors generated by the two methods are almost negligible for sufficiently small dt (such as 0.02, as in the experiments). By variation of the timesteps, we found the convergence order of integration method R(=log2||enew||−log2||eold||log2||Δtnew||−log2||Δtold||) equals 2.95, indicating that the error decays to 0.129 of its original error if the timestep is decreased by half.

This approach reduces the computational cost without compromising accuracy. The code/software described in the paper is freely available online at https://github.com/scientific-lab/Gamma-and-Beta-Band-Oscillation-in-Working-Memory. We performed our research using GPU(K80) with Ubuntu 18.04.6 LTS, Microsoft Surface Pro 5 with Windows 10 operating system, and MacBook Pro 11,1 with Windows 10 operating system.

Data analysis

Firing rate calculation

We recorded the firing activities of neurons (see Fig. 2) and counted the spikes of a certain type of neuron along the time in a single trial and approximated the firing rate by applying a one-dimension filter with a window size of 4 ms to the spike train.

Figure 2.
Figure 2.

Raster plots of spiking of neurons in the network. A, Raster plots of spiking of pyramidal cells given one cue. B, Raster plots of spiking of low-threshold nFS (blue) and high-threshold FS (red) given one cue. C, Two direction cues S1+S2 were concurrently presented to the network during the cue period (500–700 ms). D, Three directional cues S1+S2+S3 were concurrently presented to the network during the cue period. E, Two directional cues (S1 and S2) were sequentially presented to the network. F, Three directional cues (S1, S2, and S3) were sequentially presented to the network. The horizontal axis is the time in milliseconds, and the vertical axis is the preferred direction of the pyramidal cell. s1-s3 represent the indices of the cue, and the dotted lines denote the cue presented interval(s). Please note that these symbols have the same interpretations in the following figures.

Local field potential and spectrograms

We approximated the local field potential (LFP) using the total synaptic currents onto pyramidal cells following the method proposed by Mazzoni et al. (2015): LFP(t)=∑pyr[IAMPA(t−τAMPA) + INMDA(t)−1.65×IGABA(t)], (24)where τAMPA=6 ms. We applied a one-dimension filter with a window size of 4 ms to the LFP and calculated the power spectrogram of LFP using continuous wavelet transform with Morlet wavelet. The background energy of our model demonstrates the 1/f characteristic for high frequency oscillation (refer to Miller and Wang, 2006; Wang, 2010). We first computed the power-frequency correlation of our model without external input shown in Figure 3A. We observed that the power for low-frequency oscillation remains almost constant, but decays in accordance with a power-law trend. The power-frequency data were smoothed using an average of 250 points and represented by the red line to signify the background power of oscillation. We subtracted the value of the red line from the simulated power. Figure 3B illustrates the power of LFP, which has been normalized by the red curve in Figure 3A. This normalization is applied in subsequent spectral analysis to obtain the normalized power.

Figure 3.
Figure 3.

The oscillation in the network. A, The power of raw simulated LFP of the network. The black dots plot the power against frequency on a logarithmic scale for LFP and has a power-law distribution. The red curve represents the 250-point average smooth curve of the log-scaled power and the log-scaled frequency. B, The power of LFP normalized by the red curve in panel A. C, Example of the population firing rate of pyramidal cells in a single trial. D, Example of the spectrogram of normalized LFP shown in panel C based on Morlet wavelet. E, Approximated LFP in single trial as in panel C. F, Zoom-in on the approximated local field potential around the time of γ and β burst as seen in the spectrogram in panel D. Blue curve shows the local field potential filtered at 42 Hz, and red curve shows local field potential filtered at 17 Hz.

Burst rate calculation

A burst in a particular frequency is defined as an interval when spectral power exceeds twice the standard deviation above the average value of that frequency and lasts for at least three cycles (Lundqvist et al., 2016). After we located the bursts in a single trail, we applied the same process to multiple trails and calculated the frequency of the bursts in the β (10–30 Hz) or γ (35–100 Hz) band over time. Similar to how we extracted the firing rate, we filtered the calculated burst frequency with a one-dimensional Gaussian filter, with a window size of 60 ms, to approximate the burst rate for a certain frequency range.

Statistical comparison

In our study, we employed independent sample t tests and ANOVA to compare between groups after conducting tests for homogeneity of variance and making necessary corrections. We also provided the effect size and confidence interval for each significant test result. When necessary, we conducted multiple comparison analyses to offer comprehensive results. All statistical analyses were performed using MATLAB. Our model generated substantial data from multiple simulations with high power and narrow confidence intervals. Based on the results from the online power and sample size calculator (https://www.gigacalculator.com/calculators/), the power for all statistical tests were close to 100%.

Results

The continuous spiking network can hold concurrent and sequential cues

We considered a spiking neural network consisting of one population of 4098 pyramidal neurons, one population of 512 FS, and one population of 512 nFS. Pyramidal cells are uniformly distributed on a ring according to their preferred cue angles (Fig. 1A; Compte et al., 2000), mimicking the columnar organization of the monkey PFC (Goldman-Rakic, 1995; Rao et al., 1999; Constantinidis et al., 2001). Using this continuous spiking network model, we investigated how WM manipulates sequential and concurrent cues. We first presented one direction to the network for a brief period and then withdraw the cue to explore the spatiotemporal pattern of neural activity elicited by the cue presentation. As shown in Figure 2A, pyramidal cells showed relatively sparse and regular discharges in the absence of stimulus from 0 to 500 ms. At 500–700 ms, the stimulation elicited an intense firing in the cued direction. When the cue was withdrawn, this intense discharge in the preferred direction persisted and showed oscillatory activity throughout the 700- to 2000-ms delay period. We call this sustained oscillatory activity in the cued direction “oscillatory activity bump” and believe it maintains information about the cue. Pyramidal cells in other directions were unaffected during and after the cue presentation and continued to discharge sparsely and regularly. The inhibitory interneurons nFS (Fig. 2B, blue) and FS (Fig. 2B, red) showed opposite trends. In the prestimulus period (0–500 ms), the nFS discharged as sparsely and regularly as the pyramidal cell, whereas the FS barely discharged. Within 500–700 ms of stimulus presentation, the nFS stopped firing in the preferred direction, whereas the FS began to fire strongly and intensely in the preferred direction. This firing pattern continued throughout the 700- to 2000-ms delay period after stimulus withdrawal. The stimulus had little effect on the nonpreferred nFS and FS, and their spikes remained sparse or absent. Figure 4B provides additional information about the average membrane potential of the cells. It also shows examples of pyramidal cells, FS, and nFS in the preferred and nonpreferred directions. We then presented two or three directional cues to the network sequentially or concurrently. Given concurrent cues, the network simultaneously elicits two or three distinct oscillatory activity bumps that persist throughout the delay period (Fig. 2C,D), and each distinct activity bump maintains the information of one cued direction. For sequential cues, the network elicits distinct oscillatory activity bumps one after another (Fig. 2E,F). The later presented cue evokes a new persistent oscillatory activity bump and does not disrupt the persistent oscillatory activity bump elicited by the previous cue. These results show that our continuous network can store not only concurrent information but also sequential information of WM.

Figure 4.
Figure 4.

Voltage traces from neurons. A, Average membrane potential for pyramidal dendrites (top left), somas (top right), FS cells (bottom left), and nFS cells (bottom right) in response to one cue. Dotted lines indicate the presence of the stimulus. B, Example of single-cell membrane potential for pyramidal dendrites (top left), somas (top right), FS cells (bottom left), and nFS cells (bottom right) in the preferred and nonpreferred direction in response to a cue. The voltage trace of neurons in the preferred direction is represented by the black line and those in the nonpreferred direction by the blue line.

Labor division of interneurons in γ and β oscillation

To investigate the oscillatory behavior in WM, we calculated the firing rate of pyramidal cells and the LFP given one directional cue in a single trial (Fig. 3). The population firing rate and LFP oscillate before and after the cue presentation and throughout the delay period (Fig. 3C,E). The spectrogram of the LFP shows that the network oscillates in the β band before the cue presentation but oscillates in the β and γ bands during the delay period (Fig. 3D,F). Figures 4 and 5 provide an overview of the average membrane potential (Fig. 4A) and gating variables (Fig. 5A) during a single trial for different neurons or compartments. In addition, these figures provide examples of individual neurons in both preferred and nonpreferred directions (Figs. 4B, 5B). Overall, the mean membrane potential of pyramidal cells [Mpre = −55.054, Maft = −54.918, t = −43.399, p < 0.001, d = 0.430, cl = [−0.142,−0.130] (Fig. 4A), where Mpre is the average membrane potential before cue presentation and Maft represets that after the cue] and FS (Mpre = −67.212, Maft = −67.117, t = −24.981, p < 0.001, d = 0.224, cl = [−0.103,−0.088]; Fig. 4A) increased, and that of nFS cells decreased after stimulus presentation (Mpre = −53.999, Maft = −54.598, t = 542.997, p < 0.001, d = 3.464, cl = [0.596,0.601]; Fig. 4A). We calculated the total synaptic currents received by the interneurons (Fig. 6). We found that the cue presentation increases the synaptic input to FS close to the presented direction (Mpre = 0.057, Maft = 0.330, t = −271.328, p < 0.001, d = 1.360, cl = [−0.274,−0.270]; Fig. 6A) and activates these FS (Fig. 2B) during the cue presentation and throughout the delay period. Cue presentation has little effect on the synaptic input to FS far away from the cued direction (Mpre = 0.049, Maft = 0.041, t = 13.920, p < 0.001, d = 0.092, cl = [0.007,0.009]; Fig. 6B) and cannot activate them (Fig. 2B). We found that synaptic currents input to nFS near the cued direction decrease because of the inhibition from active FS during cue presentation and the delay period (Mpre = 0.081, Maft = −0.215, t = 188.834, p < 0.001, d = 0.820, cl = [0.293,0.299]; Fig. 6C), and these nFS are suppressed (Fig. 2B). However, the input to nFS far away from the cued direction does not change much (Mpre = 0.079, Maft = 0.080, t = −3.175, p = 0.002, d = 0.021, cl = [−0.002,−0.000]; Fig. 6D). In summary, before the cue presentation, weak excitation is insufficient to activate FS, and the interaction between pyramidal cells and nFS dominates the network activity, resulting in β band oscillation. The cue presentation induces strong excitation and activates FS, and the interaction between pyramidal cells and FS leads to γ band oscillation in the network. Therefore, these results show a division of labor between FS and nFS in γ band and the β band oscillations, as shown in Figure 6E.

Figure 5.
Figure 5.

Gating variables form neurons. A, Average AMPA, NMDA, and GABA gating variables associated with presynaptic pyramidal dendrites, somas, FS cells and nFS cells in response to one cue. Dashed lines indicate the presence of the stimulus. B, Example of single cell gating variable in response to one cue. The gating variables of individual neurons are shown at the same positions as the corresponding average gating variables in A. The gating variables of neurons in the preferred direction are represented by the black line and those in the nonpreferred direction by the blue line.

Figure 6.
Figure 6.

Labor division of FS and nFS interneuron. A, The total synaptic input to FS close to the cued direction. B, The total synaptic input to FS far away from the cued direction. C, The total synaptic input to nFS close to the cued direction. D, The total synaptic input to nFS far away from the cued direction. E, The schematic diagram of labor division of interneurons in the γ band and β band oscillation. Interaction between pyramidal cell and nFS results in β band oscillation, while interaction between pyramidal cell and FS leads to γ oscillation.

γ Power increases in the cued direction

Experiments (Pesaran et al., 2002) showed that the increase in γ band power was sustained during the delay period in the preferred direction, while γ band power did not change in the anti-preferred direction. We calculated the LFP based on the pyramidal cells close to the cued direction, which is the preferred direction, and the pyramidal cells far away from the cued direction, which is the anti-preferred direction. We found that in the preferred direction, the LFP mainly consists of β band oscillation before the cue presentation, and a strong γ band oscillation appears after the cue presentation in a single trial (Fig. 7A) or averaged over trials (Fig. 7C). The power of γ band increases significantly (Mpre = 23.942, Maft = 26.616, t = −80.835, p < 0.001, d = 0.784, cl = [−2.739,−2.609]), but the power of β band decreases (Mpre = 37.157, Maft = 35.381, t = 76.696, p < 0.001, d = 0.737, cl = [1.731,1.822]; Fig. 7E). There is little increase in γ band oscillation (Mpre = 26.370, Maft = 26.741, t = −54.777, p < 0.001, d = 0.479, cl = [−0.385,−0.358]; Fig. 7B) for the anti-preferred direction, and the β band oscillation is not decreased (Mpre = 37.752, Maft = 37.450, t = 47.799, p < 0.001, d = 0.348, cl = [0.289,0.314]; Fig. 7F). Note that there was an early increase in γ and β band power in Figure 7E,F. We think this was an artifact of the cwtft algorithm. CWTFT calculates the power of a signal by locally integrating the multiplication of signals with shiftable symmetric wavelets of different window sizes and center frequencies. In the same wavelet window, wavelet analysis is isotropic, thus unable to distinguish components before and after stimulus presentation. When the signal frequency undergoes significant changes within a certain window, the window that contains the changing signal cannot differentiate between the prechange and postchange components, showing an earlier increase in power (similar to Kaiju et al., 2017). In short, cue presentation increases the power of the γ oscillation and decreases β oscillation in the preferred direction (Fig. 7E,F), which are consistent with the experimental observations (Pesaran et al., 2002).

Figure 7.
Figure 7.

Oscillation power in the preferred direction and anti-preferred direction. A, The spectrogram for LFP approximated from the pyramidal cells preferring the cued direction in a single trial. B, The spectrogram for LFP approximated from the pyramidal cells far away from the cued direction in a single trial. C, The average of spectrogram for preferred direction over 900 trials. D, The average of the spectrogram for anti-preferred direction over 900 trials. E, The average normalized power (β band:10–30 Hz, γ band: 35–100 Hz) for the preferred direction. The red line denotes the γ power, and the blue line denotes the β power. F, The average normalized power for nonpreferred direction. The red line denotes the γ power, and the blue line denotes the β power.

WM load enhances the γ power during the delay period

We further analyzed the spectrogram of LFP given concurrent or sequential cues. The average spectrogram of LFP over 500 trials exhibits several characteristics. First, as previously shown, the network exhibits β band oscillation before the cue presentation, and the β power tend to decreases after cue presentation (two concurrent cues: Mpre = 107.469, Maft = 102.083, t = 39.865, p < 0.001, d = 0.412, cl = [5.121,5.651]; three concurrent cues: Mpre = 117.394, Maft = 108.605, t = 25.929, p < 0.001, d = 0.272, cl = [8.124,9.453]; two sequential cues: Mpre = 100.995, Maft = 97.639, t = 83.571, p < 0.001, d = 0.549, cl = [3.277,3.435]; three sequential cues: Mpre = 98.238, Maft = 98.990, t = −15.451, p < 0.001, d = 0.085, cl = [−0.847,−0.657]). Second, the network elicits higher γ band oscillation (around 60–70 Hz) on cue presentation (two concurrent cues: Mpre = 62.763, Maft = 68.745, t = −43.986, p < 0.001, d = 0.450, cl = [−6.249,−5.715]; three concurrent cues: Mpre = 68.370, Maft = 75.840, t = −21.950, p < 0.001, d = 0.233, cl = [−8.137,−6.803]; two sequential cues: Mpre = 58.800, Maft = 66.297, t = −133.050, p < 0.001, d = 1.081, cl = [−7.607,−7.386]; three sequential cues: Mpre = 58.202, Maft = 70.287, t = −193.409, p < 0.001, d = 1.299, cl = [−12.207,−11.963]) and γ band oscillation maintains during the delay period. Third, the γ power remains almost constant during the delay period for concurrent cues (Fig. 8A,B). In contrast, the γ power increases along with the arrival of new cues in sequential situations and the power is maintained throughout the delay period (two sequential cues: F = 13,057.29, p < 0.001; η2= 0.225, three sequential cues: F = 24,175.79, p < 0.001, η2= 0.446, with multiple comparisons p < 0.001; Fig. 8C,D). We show the average γ power in sequential trials in Figure 8E,F. The trajectory of γ power is similar to the results shown by Lundqvist et al. (2016; their Figure 6), albeit with a slight distinction: the γ band exhibited an increase following the presentation of the second stimulus before the presentation of the third stimulus.

Figure 8.
Figure 8.

The average spectrogram of LFP given concurrent or sequential cues over 500 trials (A–D) and the average normalized γ power for sequential cue conditions (E, F). A, The average spectrogram for two concurrent cues. B, The average spectrogram for three concurrent cues. C, The average spectrogram for two sequential cues. D, The average spectrogram for three sequential cues. E, The average normalized γ power for two sequential cues. F, The average normalized γ power for three sequential cues.

We performed a Fourier analysis on the LFP given concurrent cues, and the results are consistent with the spectrogram analysis. γ power increases with the WM load (F = 71.06, p < 0.001; Fig. 9). Three cues lead to the highest γ power and one cue results in the lowest γ power. The mean of the γ power are M1cue = 0.493, M2cue = 1.060, and M3cue = 2.357 for 1, 2, and 3 cues, respectively (multiple comparisons p12 = 0.002, p13 < 0.001, p23 < 0.001; Fig. 9). The dependence of γ power on WM load can be explained by network activity. Strong excitation induced by cue presentation activates high-threshold FS and the interaction between pyramidal cells and FS results in γ oscillation. One cue elicits one distinct activity bump with γ oscillation and more cues elicit more distinct activity bumps with γ oscillation (Fig. 2), suggesting that more neurons participate in γ oscillation and result in higher γ power.

Figure 9.
Figure 9.

The normalized power of local field potential. γ power induced by two concurrent cues (red line) is higher than that induced by one cue (blue line) but lower than that resulting from three concurrent cues (black line).

Brief bursts of narrow-band γ and β oscillation in a single trial

Although the raster plots of spiking in Figure 2 show oscillatory activity, the single spectrogram of LFP (Figs. 3B, 7A,B) shows a different scenario: γ and β oscillations occurred irregularly in the form of short bursts that are narrow and variable in frequency. This observation is similar to a previous monkey experiment (Lundqvist et al., 2016). Here, we identified the γ and β bursts using the algorithm proposed in the reference (Lundqvist et al., 2016). We identified each burst in 500 trials with sequential cues and plotted each burst as a point in Figure 10A. We also calculated the β and γ burst rates (Fig. 10B). The β burst rate decreases slightly after cue presentation (Mpre = 0.0207, Maft = 0.0196, t = 21.08, p < 0.001; Fig. 7F) and remains almost constant during the delay period (Fig. 10B). The γ burst rate increases when the cue is presented (Mpre = 0.0570, Maft = 0.0685, t = −48.93, p < 0.001; Fig. 7F) and increases with increasing WM load (F = 229.69, p < 0.001, multiple comparisons p < 0.001; Fig. 10B). During the delay period, the γ burst shows a slightly increasing trend, which is consistent with the experimental observation by Lundqvist et al. (2016; their Figure 8B).

Figure 10.
Figure 10.

The oscillatory bursts. A, The raster plot shows the occurrence of oscillatory bursts in 500 trials of simulation. B, The rate of β and γ bursts.

Discussion

We have developed a working memory model that utilizes a continuous attractor structure based on the bistable excitatory-inhibitory microcircuit consisting of FS, nFS, and pyramidal cells in PFC. The main contributions of our work are threefold. First, the proposed model can maintain the information of sequentially presented cues and concurrently presented cues. Previous oscillatory and dynamic models collapse in the presence of concurrent information because these models rely on the oscillatory phase or timing to encode the information (Constantinidis et al., 2018). The persistent activity model has not been used to manipulate sequential information (Edin et al., 2009; Wei et al., 2012). Second, we identified the mechanism underlying the enhancement of γ power in the WM task. Cue presentation elicits strong excitatory recurrent synaptic currents that activate high-threshold FS. The interaction between pyramidal cells and FS leads to γ oscillations and enhances γ power during the delay period. Third, we demonstrated the division of labor between high-threshold FS and nFS for oscillation in the WM. The interaction between nFS and pyramidal cells leads to β oscillation, while the interaction between FS and pyramidal cells leads to γ oscillation.

We showed that the transient and sparkling γ and β burst in a single spectrogram of LFP can be approximated from the persistent oscillatory activity of the network (Fig. 3B), and the γ burst rate increases with the WM load (Fig. 8B). The discrete and sparkling γ and β burst in a single LFP spectrogram has been considered as evidence that WM is manifested through transient or discrete oscillatory dynamics, rather than sustained firing (Lundqvist et al., 2016). However, our results indicate that transient and sparkling γ and β bursts in a single LFP spectrogram may indeed result from persistent oscillatory firing. Therefore, our results favor the theory of persistent activity underlying WM.

In this work, we have demonstrated the division of labor between FS and nFS in generating γ and β oscillation. The interplay between high-threshold FS and pyramidal cells results in the γ oscillation, whereas the interaction between nFS and pyramidal cells leads to the β oscillation. Before the cue presentation in WM tasks, the weak recurrent excitation from spontaneous discharges of pyramidal cells can activate low-threshold nFS but not the high-threshold FS, resulting in β oscillation. Cue presentation evokes strong activity of pyramidal cells and strong recurrent excitation to activate FS. As a result, the network exhibits γ oscillation. This observation is consistent with experimental results that optogenetic activation of FS induces γ rhythm (Cardin et al.,2009) and the spiking of somatostatin (SOM) and parvalbumin (PV) cells differentially correlates with β and γ oscillations and activation of PV cells enhances γ oscillation in V1 area (Chen et al., 2017).

We computationally identified that γ power can be enhanced by recruiting more high-threshold FS because of the stronger excitation evoked by the cue presentation in WM. This mechanism may be common to other cognitive processes, such as attention and sensory information processing in perceptual decision-making. It has been extensively observed in cognitive processes such as WM (Jensen et al., 2007; Yamamoto et al., 2014), attention (Engel et al., 2001; Fries et al., 2001; Jensen et al., 2007; Fries, 2009; Gregoriou et al., 2009; Kim et al., 2016), and sensory information processing in perceptual decision-making (Frien et al., 2000; Siegel and König, 2003; Kayser and König, 2004; Henrie and Shapley, 2005; Liu and Newsome, 2006; Berens et al., 2008). Studies had observed that the selected neuronal population increased their firing rate when attention was directed to their receptive field (Luck et al., 1997; McAdams and Maunsell, 1999; Bichot et al., 2015; Thiele et al., 2016), implying strong excitation and potentially recruiting more FS neurons. γ oscillations in the sensory cortex are often considered as a proxy for the encoding of sensory evidence during perceptual decision-making because of visual γ band activity and its dependence on stimulus strength and various stimulus features (Frien et al., 2000; Siegel and König, 2003; Kayser and König, 2004; Henrie and Shapley, 2005; Liu and Newsome, 2006; Berens et al., 2008). Based on these, we hypothesize that sensory input would evoke a stronger excitation in the sensory cortex, which could recruit more FS and increase the γ power.

Worth noting that the synapses from pyramidal cells to the nFS are dynamic synapses with short-term potentiation (STP). STP plays an essential role in the modulation of γ and β bursts. In the biological brain, there is short-term synaptic facilitation between pyramidal cells and nFS and short-term synaptic depression between pyramidal cells and FS. These short-term plasticities allow the brain to respond differently to afferents of different durations. Brief currents cause γ oscillations, while prolonged currents cause β oscillations. Based on this mechanism, these observations are reproduced in previously published work (Feng et al., 2021). In our simulation, we only briefly present the stimulus to the model, resulting in oscillations at the γ frequency, and short-term synaptic facilitation between pyramidal cells and nFS helps control excitability in our model. The synapses equipped with STP can temporarily increase inhibition by increasing the effective excitation from the pyramidal cells to the nFS and expanding the range of parameters to keep excitation and inhibition in balance. Too much excitation can lead to the merging of different information (Wei et al., 2012) or false memories (Edin et al., 2009), while too much inhibition can prevent the network from maintaining information. Suppose the synapses from pyramidal cells to the nFS are fixed without STP; the network is prone to produce spurious bursts of activity because of too much excitation or to forget the information because of too much inhibition. At present, we have not considered the effects of prolonged stimulus. It would be interesting to investigate how the brain adapts and responds to prolonged stimulus in the future.

Acknowledgments

Acknowledgements: We thank Prof. XiaoJing Wang for his helpful suggestions.

Footnotes

  • The authors declare no competing financial interests.

  • This work was supported by the National Natural Science Foundation of China Grant No. 32171094.

This is an open-access article distributed under the terms of the Creative Commons Attribution 4.0 International license, which permits unrestricted use, distribution and reproduction in any medium provided that the original work is properly attributed.

References

  1. ↵
    Baddeley A (2003) Working memory: looking back and looking forward. Nat Rev Neurosci 4:829–839.
    OpenUrlCrossRefPubMed
  2. ↵
    Bays PM, Husain M (2008) Dynamic shifts of limited working memory resources in human vision. Science 321:851–854. https://doi.org/10.1126/science.1158023 pmid:18687968
    OpenUrlAbstract/FREE Full Text
  3. ↵
    Berens P, Keliris GA, Ecker AS, Logothetis NK, Tolias AS (2008) Feature selectivity of the gamma-band of the local field potential in primate primary visual cortex. Front Neurosci 2:199–207. https://doi.org/10.3389/neuro.01.037.2008 pmid:19225593
    OpenUrlCrossRefPubMed
  4. ↵
    Bichot NP, Heard MT, DeGennaro EM, Desimone R (2015) A source for feature-based attention in the prefrontal cortex. Neuron 88:832–844. https://doi.org/10.1016/j.neuron.2015.10.001 pmid:26526392
    OpenUrlCrossRefPubMed
  5. ↵
    Braitenberg V, Schutz A (1991) Anatomy of the cortex: studies of brain function. Berlin: Springer.
  6. ↵
    Cardin JA, Carlén M, Meletis K, Knoblich U, Zhang F, Deisseroth K, Tsai LH, Moore CI (2009) Driving fast-spiking cells induces gamma rhythm and controls sensory responses. Nature 459:663–667. https://doi.org/10.1038/nature08002 pmid:19396156
    OpenUrlCrossRefPubMed
  7. ↵
    Chen G, Zhang Y, Li X, Zhao XC, Ye Q, Lin YX, Tao HW, Rasch M, Zhang XH (2017) Distinct inhibitory circuits orchestrate cortical beta and gamma band oscillations. Neuron 96:1403–1418.e6. https://doi.org/10.1016/j.neuron.2017.11.033 pmid:29268099
    OpenUrlCrossRefPubMed
  8. ↵
    Compte A, Brunel N, Goldman-Rakic PS, Wang XJ (2000) Synaptic mechanisms and network dynamics underlying spatial working memory in a cortical network model. Cereb Cortex 10:910– 923. https://doi.org/10.1093/cercor/10.9.910 pmid:10982751
    OpenUrlCrossRefPubMed
  9. ↵
    Constantinidis C, Franowicz MN, Goldman-Rakic PS (2001) Coding specificity in cortical microcircuits: a multiple-electrode analysis of primate prefrontal cortex. J Neurosci 21:3646–3655.
    OpenUrlAbstract/FREE Full Text
  10. ↵
    Constantinidis C, Funahashi S, Lee D, Murray JD, Qi XL, Wang M, Arnsten AFT (2018) Persistent spiking activity underlies working memory. J Neurosci 38:7020–7028. https://doi.org/10.1523/JNEUROSCI.2486-17.2018 pmid:30089641
    OpenUrlAbstract/FREE Full Text
  11. ↵
    Cowan N (2010) The magical mystery four: how is working memory capacity limited, and why? Curr Dir Psychol Sci 19:51–57. https://doi.org/10.1177/0963721409359277 pmid:20445769
    OpenUrlCrossRefPubMed
  12. ↵
    Domhof JW, Tiesinga PH (2021) Flexible frequency switching in adult mouse visual cortex is mediated by competition between parvalbumin and somatostatin expressing interneurons. Neural Comput 33:926–966. https://doi.org/10.1162/neco_a_01369 pmid:33513330
    OpenUrlPubMed
  13. ↵
    Edin F, Klingberg T, Johansson P, McNab F, Tegnér J, Compte A (2009) Mechanism for top-down control of working memory capacity. Proc Natl Acad Sci U S A 106:6802–6807. https://doi.org/10.1073/pnas.0901894106 pmid:19339493
    OpenUrlAbstract/FREE Full Text
  14. ↵
    Engel A, Fries P, Singer W (2001) Dynamic predictions: oscillations and synchrony in top-down processing. Nat Rev Neurosci 2:704– 716. https://doi.org/10.1038/35094565 pmid:11584308
    OpenUrlCrossRefPubMed
  15. ↵
    Feng F, Headley DB, Nair SS (2021) Model neocortical microcircuit supports beta and gamma rhythms. 2021 10th International IEEE/EMBS Conference on Neural Engineering (NER), May 4–6, pp 91–94. IEEE. https://doi.org/10.1109/ner49283.2021.9441199 pmid:35469138
  16. ↵
    Fitzpatrick D, Lund JS, Blasdel GG (1985) Intrinsic connections of macaque striate cortex: afferent and efferent connections of lamina 4C. J Neurosci 5:3329–3349.
    OpenUrlAbstract/FREE Full Text
  17. ↵
    Frien A, Eckhorn R, Bauer R, Woelbern T, Gabriel A (2000) Fast oscillations display sharper orientation tuning than slower components of the same recordings in striate cortex of the awake monkey. Eur J Neurosci 12:1453–1465. https://doi.org/10.1046/j.1460-9568.2000.00025.x pmid:10762373
    OpenUrlCrossRefPubMed
  18. ↵
    Fries P (2009) Neuronal gamma-band synchronization as a fundamental process in cortical computation. Annu Rev Neurosci 32:209–224. https://doi.org/10.1146/annurev.neuro.051508.135603 pmid:19400723
    OpenUrlCrossRefPubMed
  19. ↵
    Fries P, Reynolds J, Rorie A, Desimone R (2001) Modulation of oscillatory neuronal synchronization by selective visual attention. Science 291:1560–1563. https://doi.org/10.1126/science.1055465 pmid:11222864
    OpenUrlAbstract/FREE Full Text
  20. ↵
    Goldman-Rakic PS (1995) Cellular basis of working memory. Neuron 14:477–485.
    OpenUrlCrossRefPubMed
  21. ↵
    Gregoriou GG, Gotts SJ, Zhou H, Desimone R (2009) High-frequency, long-range coupling between prefrontal and visual cortex during attention. Science 324:1207–1210. https://doi.org/10.1126/science.1171402 pmid:19478185
    OpenUrlAbstract/FREE Full Text
  22. ↵
    Gorgoraptis N, Catalao RFG, Bays PM, Husain M (2011) Dynamic updating of working memory resources for visual objects. J Neurosci 31:8502–8511. https://doi.org/10.1523/JNEUROSCI.0208-11.2011 pmid:21653854
    OpenUrlAbstract/FREE Full Text
  23. ↵
    Hahn G, Kumar A, Schmidt H, Knösche TR, Deco G (2022) Rate and oscillatory switching dynamics of a multilayer visual microcircuit model. Elife 11:e77594. https://doi.org/10.7554/eLife.77594
    OpenUrl
  24. ↵
    Henrie JA, Shapley R (2005) LFP power spectra in V1 cortex: the graded effect of stimulus contrast. J Neurophysiol 94:479–490. https://doi.org/10.1152/jn.00919.2004 pmid:15703230
    OpenUrlCrossRefPubMed
  25. ↵
    Jahr CE, Stevens CF (1990) Voltage dependence of NMDA-activated macroscopic conductances predicted by single-channel kinetics. J Neurosci 10:3178– 3182. https://doi.org/10.1523/JNEUROSCI.10-09-03178.1990 pmid:1697902
    OpenUrlAbstract/FREE Full Text
  26. ↵
    Jensen O, Lisman JE (1998) An oscillatory short-term memory buffer model can account for data on the Sternberg task. J Neurosci 18:10688–10699. https://doi.org/10.1523/JNEUROSCI.18-24-10688.1998 pmid:9852604
    OpenUrlAbstract/FREE Full Text
  27. ↵
    Jensen O, Kaiser J, Lachaux JP (2007) Human gamma-frequency oscillations associated with attention and memory. Trends Neurosci 30:317–324. https://doi.org/10.1016/j.tins.2007.05.001 pmid:17499860
    OpenUrlCrossRefPubMed
  28. ↵
    Kaiju T, Doi K, Yokota M, Watanabe K, Inoue M, Ando H, Takahashi K, Yoshida F, Hirata M, Suzuki T (2017) High spatiotemporal resolution ECoG recording of somatosensory evoked potentials with flexible micro-electrode arrays. Front Neural Circuits 11:20.
    OpenUrl
  29. ↵
    Kayser C, König P (2004) Stimulus locking and feature selectivity prevail in complementary frequency ranges of V1 local field potentials. Eur J Neurosci 19:485–489. https://doi.org/10.1111/j.0953-816x.2003.03122.x pmid:14725644
    OpenUrlCrossRefPubMed
  30. ↵
    Kim H, Ährlund-Richter S, Wang X, Deisseroth K, Carlén M (2016) Prefrontal parvalbumin neurons in control of attention. Cell 164:208–218. https://doi.org/10.1016/j.cell.2015.11.038 pmid:26771492
    OpenUrlCrossRefPubMed
  31. ↵
    Lisman JE, Idiart MA (1995) Storage of 7+/- 2 short-term memories in oscillatory subcycles. Science 267:1512–1515. https://doi.org/10.1126/science.7878473 pmid:7878473
    OpenUrlAbstract/FREE Full Text
  32. ↵
    Liu J, Newsome WT (2006) Local field potential in cortical area MT: stimulus tuning and behavioral correlations. J Neurosci 26:7779–7790. https://doi.org/10.1523/JNEUROSCI.5052-05.2006 pmid:16870724
    OpenUrlAbstract/FREE Full Text
  33. ↵
    Luck SJ, Chelazzi L, Hillyard SA, Desimone R (1997) Neural mechanisms of spatial selective attention in areas V1, V2, and V4 of macaque visual cortex. J Neurophysiol 77:24–42. https://doi.org/10.1152/jn.1997.77.1.24 pmid:9120566
    OpenUrlCrossRefPubMed
  34. ↵
    Lundqvist M, Rose J, Herman P, Brincat SL, Buschman TJ, Miller EK (2016) Gamma and beta bursts underlie working memory. Neuron 90:152–164. https://doi.org/10.1016/j.neuron.2016.02.028 pmid:26996084
    OpenUrlCrossRefPubMed
  35. ↵
    Lundqvist M, Herman P, Warden MR, Brincat SL, Miller EK (2018) Gamma and beta bursts during working memory readout suggest roles in its volitional control. Nat Commun 9:394. https://doi.org/10.1038/s41467-017-02791-8 pmid:29374153
    OpenUrlCrossRefPubMed
  36. ↵
    Mazzoni A, Lindén H, Cuntz H, Lansner A, Panzeri S, Einevoll GT (2015) Computing the local field potential (LFP) from integrate-and-fire network models. PLoS Comput Biol 11:e1004584. https://doi.org/10.1371/journal.pcbi.1004584 pmid:26657024
    OpenUrlCrossRefPubMed
  37. ↵
    McAdams CJ, Maunsell JH (1999) Effects of attention on orientation-tuning functions of single neurons in macaque cortical area V4. J Neurosci 19:431–441. https://doi.org/10.1523/JNEUROSCI.19-01-00431.1999 pmid:9870971
    OpenUrlAbstract/FREE Full Text
  38. ↵
    Mi Y, Katkov M, Tsodyks M (2017) Synaptic correlates of working memory capacity. Neuron 93:323–330. https://doi.org/10.1016/j.neuron.2016.12.004 pmid:28041884
    OpenUrlCrossRefPubMed
  39. ↵
    Miller P, Wang XJ (2006) Power-law neuronal fluctuations in a recurrent network model of parametric working memory. J Neurophysiol 95:1099–1114. https://doi.org/10.1152/jn.00491.2005 pmid:16236788
    OpenUrlCrossRefPubMed
  40. ↵
    Mongillo G, Barak O, Tsodyks M (2008) Synaptic theory of working memory. Science 319:1543–1546. https://doi.org/10.1126/science.1150769 pmid:18339943
    OpenUrlAbstract/FREE Full Text
  41. ↵
    Pesaran B, Pezaris JS, Sahani M, Mitra PP, Andersen RA (2002) Temporal structure in neuronal activity during working memory in macaque parietal cortex. Nat Neurosci 5:805–811. https://doi.org/10.1038/nn890 pmid:12134152
    OpenUrlCrossRefPubMed
  42. ↵
    Rao SG, Williams GV, Goldman-Rakic PS (1999) Isodirectional tuning of adjacent interneurons and pyramidal cells during working memory: evidence for microcolumnar organization in PFC. J Neurophysiol 81:1903–1916.
    OpenUrlCrossRefPubMed
  43. ↵
    Siegel M, König P (2003) A functional gamma-band defined by stimulus-dependent synchronization in area 18 of awake behaving cats. J Neurosci 23:4251–4260. https://doi.org/10.1523/JNEUROSCI.23-10-04251.2003 pmid:12764113
    OpenUrlAbstract/FREE Full Text
  44. ↵
    Spiro JE, Dalva MB, Mooney R (1999) Long-range inhibition within the zebra finch song nucleus RA can coordinate the firing of multiple projection neurons. J Neurophysiol 81:3007–3020.
    OpenUrlCrossRefPubMed
  45. ↵
    Thiele A, Brandt C, Dasilva M, Gotthardt S, Chicharro D, Panzeri S, Distler C (2016) Attention induced gain stabilization in broad and narrow-spiking cells in the frontal eye-field of macaque monkeys. J Neurosci 36:7601–7612. https://doi.org/10.1523/JNEUROSCI.0872-16.2016 pmid:27445139
    OpenUrlAbstract/FREE Full Text
  46. ↵
    Troyer TW, Miller KD (1997) Physiological gain leads to high ISI variability in a simple model of a cortical regular spiking cell. Neural Comput 9:971–983. https://doi.org/10.1162/neco.1997.9.5.971 pmid:9188190
    OpenUrlCrossRefPubMed
  47. ↵
    Tuckwell HC (1988) Introduction to theoretical neurobiology. Cambridge: Cambridge University Press.
  48. ↵
    Veit J, Hakim R, Jadi MP, Sejnowski TJ, Adesnik H (2017) Cortical gamma band synchronization through somatostatin interneurons. Nat Neurosci 20:951–959. https://doi.org/10.1038/nn.4562 pmid:28481348
    OpenUrlCrossRefPubMed
  49. ↵
    Wang XJ (1998) Calcium coding and adaptive temporal computation in cortical pyramidal neurons. J Neurophysiol 79:1549–1566. https://doi.org/10.1152/jn.1998.79.3.1549 pmid:9497431
    OpenUrlCrossRefPubMed
  50. ↵
    Wang XJ (2010) Neurophysiological and computational principles of cortical rhythms in cognition. Physiol Rev 90:1195–1268. https://doi.org/10.1152/physrev.00035.2008 pmid:20664082
    OpenUrlCrossRefPubMed
  51. ↵
    Wang XJ, Yang GR (2018) A disinhibitory circuit motif and flexible information routing in the brain. Curr Opin Neurobiol 49:75–83. https://doi.org/10.1016/j.conb.2018.01.002 pmid:29414069
    OpenUrlCrossRefPubMed
  52. ↵
    Wang XJ, Tegnér J, Constantinidis C, Goldman-Rakic PS (2004) Division of labor among distinct subtypes of inhibitory neurons in a cortical microcircuit of working memory. Proc Natl Acad Sci U S A 101:1368–1373. https://doi.org/10.1073/pnas.0305337101 pmid:14742867
    OpenUrlAbstract/FREE Full Text
  53. ↵
    Wei ZQ, Wang XJ, Wang DH (2012) From distributed resources to limited slots in multiple-item working memory: a spiking network model with normalization. J Neurosci 32:11228– 11240. https://doi.org/10.1523/JNEUROSCI.0735-12.2012 pmid:22895707
    OpenUrlAbstract/FREE Full Text
  54. ↵
    Yamamoto J, Suh J, Takeuchi D, Tonegawa S (2014) Successful execution of working memory linked to synchronized high-frequency gamma oscillations. Cell 157:845–857. https://doi.org/10.1016/j.cell.2014.04.009 pmid:24768692
    OpenUrlCrossRefPubMed
  55. ↵
    Zhang W, Luck SJ (2008) Discrete fixed-resolution representations in visual working memory. Nature 453:233–235. https://doi.org/10.1038/nature06860 pmid:18385672
    OpenUrlCrossRefPubMed

Synthesis

Reviewing Editor: Niraj Desai, National Institute of Neurological Disorders and Stroke

Decisions are customarily a result of the Reviewing Editor and the peer reviewers coming together and discussing their recommendations until a consensus is reached. When revisions are invited, a fact-based synthesis statement explaining their decision and outlining what is needed to prepare a revision will be listed below. The following reviewer(s) agreed to reveal their identity: NONE.

Two reviewers and I read the revised manuscript and the rebuttal letter. While we all agree that the manuscript has considerable merit and would be valuable contribution to the literature, we also think it still needs major revision. One of the reviewers wrote a detailed list of concerns (appended below) that I would urge you to consider and address carefully should you choose to resubmit.

REVIEWER COMMENTS

This study successfully reproduces Gamma and Beta oscillations in a working memory test with concurrent cues by a spiking model with pyramidal cells and two types of interneurons. Importantly, the work contributes to the debate on whether a sustained activity or short-term plasticity is a more essential mechanism of working memory by showing that variable rhythmic activity can also be achieved through a model with bump attractors. The work is novel, and I believe it can make an impact in the field. However, this manuscript still requires detailed work to make conclusions convincing. I feel like reading a promising working progress report but not a combined polished manuscript. I have to guess what the authors mean in several places. A major revision is needed before publication.

One major concern is about the spectrum analyses. The results and the related arguments are intuitively right, but so far unconvincing. Compared to the previous version, the author claimed to reproduce an increase in the gamma range during the stimuli presentation. However, these results only appear in the rebuttal to the reviews rather than in the major figures of the main manuscript. In fact, even that figure, on the bottom left and right, the Gamma power should be comparable before the onset of the third stimulus at around time=1000ms. However, the gamma power almost doubles on the right compared to the left. Why that is the case? As Lundqvist showed in the top left and right panels, the gamma power should be consistent before S3 presents. This may be due to different random seeds. In this scenario, several runs should be combined into the analysis, but not let the readers wonder whether a cherry-picked result is shown.

Furthermore, with a closer look, the gamma power increases before the stimuli onset. This is not observed in the experiment and is not reasonable from a modeling aspect, since the statistic of the gamma rhythm should be the same before the stimuli onset. I suspect the results may come from a relatively large sliding window when generating the spectrograms and a relatively small time constant of the network dynamics. On the contrary, Figures 2, 4, 5, and 10 consistently showed that all the observable changes happen after the stimuli onset. If this is indeed an artifact, please fix it. If not, please explain the temporal differences between different analysis protocols convincingly.

Last, an accurate spectrum analysis requires subtracting 1/f^a background energy. For one example, see the discussion at Wang 2010, Physiol. Rev. Fig. 15B. I suspect the conclusion drawn from Figure 7 C to F and Figure 9 would be different if the authors did it correctly.

Another major concern is the quality of the manuscript. Little issues constantly pop up to bother me when I read the manuscript. Some examples are here: In the methods section, the parameters have three digits but without references or explanation. The authors should include the references from the Rebuttal. Even better, but not required, authors should include some recent references. Some of them argued the lateral innervation width of interneurons may be comparable with that of excitatory neurons. Further, the results in all the major figures are described with realistic units, but the parameters have no units associated when they are described in the Method section. Please include the units for all the parameters in future versions. In the figures, the direction of the y-axis labels is not consistent, and the figure quality is low. A thorough round of editing, maybe with the suggestions generated from a Large Linguistic Model, would help in the future submission.

A comparable minor thing, the authors misunderstood the meaning of a convergence test. A convergence test is generated by fitting the error line to the timestep in the log-log scale. To be specific, the author should check whether the error decay to one-fourth or one-half if decrease timestep by half. See https://math.stackexchange.com/questions/3115199/measuring-the-convergence-order-of-a-numerical-scheme-for-pde. If authors succeed in generating a 2-nd order scheme, the slope of the fitting should be 2.

Some minor things I found:

Line 65, “Fs” should be FS

Line 101 to 106, personal suggestion: “Without losing generality, a network of the size of 2^n can be utilized by the FFT”

Line 111 delete one “equation” from two “equation”. The font size is not consistent, either.

Line 134, for better visualization, the parameters should be included in a table.

Line 146, is the dendrite voltage reset to -60 as well when the neuron spikes? In realistic neurons, the reset is not that strong or even negligible. See papers from Matthew Larkum.

Line 149, is p1=p2? It should be.

Line 221, equation 20 is not mathematically correct. Rewrite 1[rand > dt ∗ rate] in a standard Poisson process with \lambda=1

Line 289, Why is it unfair? Can you modify the analysis code such that comparing “indirectly” but “fair”?

Line 341, figure 5, is not sufficiently explained. This may raise questions from readers. For example, why average gating variable NMDA_s_soma is constantly elaborated, as in A top second, but not in either preferred direction or non-preferred direction, as in B top second?

Line 423, the authors may want first summarize the model.

Line 433, I am confused about the “division of labor.” What is “labor”? I don’t see any “labor” for nFS in this model. The nFS input does not show any change before and after the stimulus onset, while the beta power and burst rate change. Actually, will model results change if replacing nFS cells with a constant inhibitory input?

Line 475 Replacing “It is not difficult to imagine that” with “Based on these, we hypothesize that”

Line 485 Should be “These observations are reproduced in Fend et al. 2021”.

Line 536 after “GABA,” maybe add “averaging.” Is that right?

Line 537, before the “FS” cell, add “from”? “NFS” should be “nFS.”

Author Response

Dear Editor,

Thanks for your comments and your patience. We also thank the reviewers for their constructive and helpful comments. Following the suggestions and comments of reviewers, we re-simulate the model using new parameters and the results are consistent with those shown in the original version of our manuscript. Attached please find the point-to-point response to reviewers’ comments.

Best

DH

Response to reviewers’ comments:

REVIEWER 1

The study succeeded in reproducing Gamma and Beta oscillations in a working memory test with concurrent cues by a spiking model with pyramidal cells and two types of interneurons. However, I am concerned about several issues not well addressed/described in the current draft.

1st, the connectivity between different pyramidal cells is way too narrow. In the original Compte et al.2000 that the author adapted, the connectivity Gaussian kernel has a width sigma_EE= 18 degrees, but in the current setup, it is just sigma_SS= 1 degree. The original parameter is chosen based on the tuning curve from experimental data but not chosen arbitrarily. The modification in the current model significantly reduces the interaction between the pyramidal cells that prefer different orientations. Essentially, the concurrent cues are stored by independent E populations in this model, which is not biologically plausible. In addition, the model losses the ability to explain other cognitive functions, such as decision-making.

Reply: Thank you for your advice. This parameter was chosen as 18 degree in original. Compte et. al.2000. In the revised version, we set parameter to 12.76 degree which is smaller than the original Compte et al.’s model, but is large enough for the mutliptle items working memory and is the same as in a multiple item working memory model by Wei et al. (2012) “From Distributed Resources to Limited Slots in Multiple-Item Working Memory: A Spiking Network Model with Normalization”. The results of our model have insignificant change and shown in the revised version of the manuscript. The new parameters are biologically more plausible than old version of our manuscript.

2nd, though the rhythmic behaviors in a WM task have not been well studied, several studies on the PV/SST and slow/fast oscillations have been established recently. However, these studies are not referred to or discussed in the manuscript. Some examples include, but are not limited to: Chen et al. 2017 (listed in the references but not referred to in the main context); Veit. J, .. Adesnik. H, 2017, Nature Neuroscience; Domhof. J, and Tiesinga. P, 2021, Neural computation; Hahn. G, .., Deco. G, 2022, eLife.

Reply: Thank you for the recommendation. We have incorporated these literatures in our revised manuscript, particularly, please refer to paragraph in INTRODUCTION section.

3rd, the authors claim that they reproduced the observed gamma rhythms during a WM task. However, in Lundqvist et al. 2016, fig. 6, the gamma power is strongest during the stimuli presentation but not during the delay period. This crucial mismatch is not explained in the paper.

Reply: thanks for the comments. We using the revised parameters and simulated the tasks. The new results (the bottom panels in the follow figure) we obtained are very similar to those reported by Lundqvist et al. 2016, with the highest gamma power during stimulus presentation. The gamma power during delay period increases with the items.

4th, the authors use the synaptic input to calculate the spectrum of LFP, but no voltage traces of pyramidal cells are shown. Surprisingly, even though only a small population of pyramidal cells (about 3 degrees) are storing the cue, their impact on the LFP is significant. I suspect this is because the inhibition current is so strong that it is no longer in a physical region. To illustrate this, plotting sample traces of voltage and gating variables is helpful.

Reply: We have shown the average membrane potential of neurons and that of example neuron and gating variables in Figures 4 and 5. Figure. 4 shows the average membrane potential for pyramidal dendrites, somas, FS cells, and NFS cells (Fig 4.A), with examples of neurons in preferred and non-preferred directions (Fig 4.B). Figure 5 shows mean gating variables over time for AMPA, NMDA, and GABA synapses (Fig 5.A), with examples of neurons in preferred and non-preferred directions (Fig 5.B).

Figure 4. Voltage traces from neurons. A) Average membrane potential for pyramidal dendrites (top left), somas (top right), FS cells (bottom left) and NFS cells (bottom right) in response to one cue. Dotted lines indicate the presence of stimulus. B) Example of single-cell membrane potential for pyramidal dendrites (top left), somas (top right), FS cells (bottom left) and NFS cells (bottom right) in the preferred and non-preferred direction in response to a cue. The voltage trace of neurons in the preferred direction is represented by the black line, and those in the non-preferred direction by the blue line.

Figure 5. Gating variables form neurons. A) AMPA, NMDA, and GABA gating variables for pyramidal dendrites, somas, FS cells and NFS cells in response to one cue. Dash lines indicate the presence of stimulation. B) Example of single cell gating variable in response to one cue. The gating variables of the individual neurons are shown at the same positions as the corresponding average gating variables in A). The gating variables of neurons in the preferred direction are represented by the black line, and those in the non-preferred direction by the blue line.

While neuronal and synaptic activity in the preferred direction is very pronounced, all neuronal activity and gating variables are within reasonable ranges.

5th, the model is not well explained in the method section. Many parameter values are missing in the whole manuscript. Further, how the noise is implemented is not described in the manuscript.

Reply:

Thanks for your advice! The Methods section has been updated, the model has been explained in detail, parameter descriptions have been added to the formula and code, and the noise implementation method has been explained.

6th, Since I have a question about how the model is simulated, I checked the code. The authors may disagree with me, but I don’t think the RK2 scheme is correctly implemented. At line 1130 of the WM_Oscil_git.m, the C_coeff and D_vec within the function RK2_4uNMDA are not evaluated at the half timestep but at the beginning of the timestep. This will introduce O(t) error to the scheme. Further, I noticed that the authors confused the value at the current step and the updated timestep. A few examples are on lines 600, 694, and 750. These mistakes can be avoided by storing these variables into a time-dependent vector or using the label _new, _old systematically in the code. Within those lines (694, 750), the scheme is forward Euler but not RK2, which makes the whole scheme 1st order but not 2nd order. I cannot debug the whole code, and I understand that some modification is reasonable in certain cases and can improve efficiency. If the authors believe what they did is correct, please do a convergent test and show the results. Or, the whole simulation should use the Euler scheme instead.

Reply:

Thank you for your advice, we have checked the code thoroughly.

1) In the function RK2_4uNMDA, C_coeff is equal to 1/syn_tau_f_NMDA, while D_vec is equal to 0. They are constant and take the same constant value at the beginning of the time step and at the half time step. The slopes (K1 and K2) at the half time step are C_coeff*(U_init-y_vec_temp) and U.*D_vec.*(1-y_vec_temp). They are updated when y_vec is updated to y_vec_temp.

2) On line 600, pyr_s_s_NMDA is incorrectly calculated with updated pyr_s_x_NMDA, we have corrected the error by saving the updated pyr_s_x_NMDA to pyr_s_x_NMDA_new.

3) On lines 694 and 750 pyr_s_Vm_new (pyr_d1_Vm_new) is correctly calculated with the current pyr_s_I_total. This is because pyr_s_I_total is calculated with pyr_s_I_leak, pyr_s_I_NMDA, and pyr_s_I_GABA, etc.,and they are all updated in the current step with pyr_s_Vm_prev from the last timestep.

We used the classical Foward-Euler LIF model to calculate the membrane potential and the second-order Range-Kutta method to integrate the gating variables. This choice is because changes in gating variables during neuronal firing are more complex and usually non-linear, whereas changes in the subthreshold membrane potential are approximately linear. We conducted a convergence test, providing 2 mV constant stimulation for 2000 ms. We found that for sufficiently small dt, the results obtained by the two methods are almost identical in 2000 ms of membrane potential computation. See the figures below.

Figure. Convergence test: The neurons were consistently stimulated with a 2mV stimulus for 2000 ms. A) The membrane potential of the neurons integrated with Euler (dt=0.02) and RK (dt=0.02) was almost identical, even in the last 40ms of the 2000ms stimulation (t=0.353, p=0.724). B) As the integration timestep is reduced, the error of the Euler method, measured as the difference in spike time, decreases. C) Similarly, as the integration timestep is reduced, the Euler method’s error, measured as the absolute difference in membrane potential, decreases.

Therefore, when the integration step equals 0.02, as in the experiment, there is almost no difference in accuracy between Euler and RK2 when calculating membrane potential, while Euler is more computationally efficient (de Alteriis, 2021). A combination of Euler and RK2 is commonly used to integrate membrane potentials. See Chou (2018).

[1] de Alteriis, G., & Oddo, C. M. (2021, May). Tradeoff between accuracy and computational cost of Euler and Runge Kutta ODE solvers for the Izhikevich spiking neuron model. In 2021 10th International IEEE/EMBS Conference on Neural Engineering (NER) (pp. 730-733). IEEE.

[2] Chou, T. S., Kashyap, H. J., Xing, J., Listopad, S., Rounds, E. L., Beyeler, M., ... & Krichmar, J. L. (2018, July). CARLsim 4: An open source library for large scale, biologically detailed spiking neural network simulation using heterogeneous clusters. In 2018 International joint conference on neural networks (IJCNN) (pp. 1-8). IEEE.

A few minor points in the Introduction section: 1st, “fast spiking interneurons” should be “fast-spiking interneurons.” “Non-fast spiking” should be “non-fast-spiking.” And use the abbreviation consistently in the latter text. 2. In the 2nd paragraph: .. updating, IS (not are) not well understood. 3rd, ..(ADP) with A time constant .. 4th, Last paragraph. Pyramidal cells EXCITE two types of interneurons, while nFS INs inhibit the (apical) dendritic compartment of pyramidal cells and FS INs inhibit the soma compartment.

Reply: Thank you! We rechecked the grammar and corrected those mistakes.

REVIEWER 2

In this study, authors constructed a spiking network model to investigate the underly circuitry mechanism of working memory task with sequentially or concurrently presented items. The study is overall interesting. However, I found the impact of the study is somewhat limited and the manuscript itself needs significant improvements, especially for figures part.

Reply: thanks for your comments.

1. In the first part of results, authors made several claims in several places, but they were not clear in the plots. For instances, “and FS almost do not discharge...”, “and form a localized oscillatory population activity (oscillatory activity bump)...” and “FS are suppressed after cue presentation...”. All of them is not clearly identifiable in plots. Theerfore, Figure 2 for this section needs revision for a better demonstration what authors wanted to claim. I suggest that authors should have some companying plots showing the average FR for each group as time evolved.

Reply: We apologies for the confusing descriptions.

As shown in Figure 2A, pyramidal cells showed relatively sparse and regular discharges in the absence of stimulation from 0-500 ms. At 500-700 ms, the stimulus elicited an intense firing in the cued direction. When the cue was withdrawn, this intense discharge in the preferred direction persisted and showed oscillatory activity throughout the 700-2000 ms delay period. We call this sustained oscillatory activity in the cued direction “oscillatory activity bump” and believe that it maintains information about the cue. Pyramidal cells in other directions were unaffected during and after the cue presentation and continued to fire sparsely and regularly. The inhibitory interneurons nFS (blue in Figure 2B) and FS (red in Figure 2B) showed opposite trends. In the pre-stimulus period (0-500 ms), the nFS discharged as sparsely and regularly as the pyramidal cell, whereas the FS barely discharged. Within 500-700 ms of stimulus presentation, the nFS stopped firing in the preferred direction, while the FS began to fire strongly and intensely in the preferred direction. This firing pattern continued throughout the 700-2000 ms delay period after stimulus withdrawal. The stimulus barely affected the non-preferred nFS and FS, and their spikes remained sparse or absent. Figure 4B provides additional information on the average membrane potential of the cells. It also gives examples of pyramidal cells, FS, and nFS in the preferred and non-preferred directions.

Hope everything is clear now.

2. For Figure 4C, I don’t see a clear drop in the total synaptic current for nFS after cue presents and delay period, as authors claimed in text.

Reply:

There was a relatively significant decrease in the total synaptic current of the NFS after the stimulus onset. Prior to the stimulus onset, the mean total synaptic current of the NFS was greater than 0, whereas after the stimulus presentation, the value was less than 0. We performed an independent samples t-test, which showed a significant result (Mpre=0.076, Maft=-0.241, t=100.96, p<0.001, Fig.6C).

3. In “WM load enhances the Gamma power during the delay period” section, when authors said “more cues elicit more distinct activity bumps with gamma oscillation”, is it because of the FR changes due to higher excitation brought by more cues?

Reply:

In our model, the increase in gamma power is likely due to more neurons involved in intensely firing oscillatory activity rather than a higher FR due to the increased excitation. After presenting one or three simultaneous cues, we compared the mean firing rates of 100 neurons in a preferred direction (from -4.39 to 4.39 degrees). We found that the firing rate in the three-cue condition was not higher (significantly lower) than in the one-cue condition (Mone=2.53, Mthree=2.17, t=21.72, p<0.001). In physiological experiments, the mechanism of gamma power increases with WM load may be more complex.

4. Regarding the last para of discussion, there was a recent modeling study by Feng et al., 2021 at 10th NER (International IEEE/EMBS Conference on Neural Engineering), which hypothesized one neocortical microcircuit can support both beta and gamma oscillations, but from different underlying mechanisms. This is very similar to what authors found here (Feng et al used two types of STP). Can authors discuss on it?

Reply:

Thanks for the suggestion. STP plays an essential role in the modulation of gamma and beta bursts. In the biological brain, there is short-term synaptic facilitation between pyramidal cells and nFS and short-term synaptic depression between pyramidal cells and FS. This short-term plasticity allows the brain to respond differently to afferents of different durations. Brief currents cause gamma oscillations, while prolonged currents cause beta oscillations. Based on this mechanism, Feng et al. developed a model in 2021 to reproduce these phenomena (Feng et al., 2021). In our research, we only briefly present the stimulus to the model, resulting in oscillations at the gamma frequency, and short-term synaptic facilitation between pyramidal cells and nFS helps control excitability in our model. The synapses equipped with STP can temporarily increase inhibition by increasing the effective excitation from the pyramidal cells to the nFS and expanding the range of parameters to keep excitation and inhibition in balance. Too much excitation can lead to the merging of different information (Wei et al,2012) or false memories (Edin et al,2008), while too much inhibition can prevent the network from maintaining information. Suppose the synapses from pyramidal cells to the nFS are fixed without STP; the network is prone to produce spurious bursts of activity due to too much excitation or to forget the information due to too much inhibition. At present, we have not considered the effects of prolonged stimulation. In future research, it would be interesting to investigate how the brain adapts and responds to prolonged stimulation.

5. Was there any particular reason for choosing number of 4098 pyramidal neurons, 512 FS and 512 nFS? And all the parameters authors chose were for some specific purposes or were they used for calibration of the network model biophysically?

Reply:

First, according to Marom and Shahaf (2002), the ratio of excitatory to inhibitory neurons in the brain is about 4:1. Considering only three cell types, the ratio of pyramidal neurons (4096) to inhibitory interneurons (FS+nFS, 1024) should approach 4:1. Second, when running the Fast Fourier Transform (FFT) algorithm, it is often desirable for the number of neurons to be a power of 2. The FFT algorithm recursively divides the input data into smaller subsets, performs Fourier transform computations on these subsets, and then combines the results in a particular way. If the number of neurons is a power of 2, it is well suited to the requirements of the FFT algorithm. The implementation of the model should be the same for other configurations of the number of neurons.

6. Authors should provide statistical tests where there were any comparisons.

Reply:

Thanks, we added statistical test results in the text. By the way, our model generated substantial data from multiple simulations with high power and narrow confidence intervals. We did not perform power and effect size tests because comparing these powers and confidence intervals directly with physiological standards would be unfair.

7. Minor: Please specify what band ranges are of interests in this study for Beta and Gamma in Abstract, Significance Statement and Introduction

Reply:

The range of interest is 10-30Hz for the beta oscillation and 30-100Hz for the gamma oscillation. This is now demonstrated in the abstract, significance statement and introduction.

8. Minor: “Beta” and “Gamma” should be lower case if they were not in the beginning of the sentences. STP and STF need be extended in the text.

Reply:

Thanks, these problems are fixed in the new manuscript.

Back to top

In this issue

eneuro: 10 (11)
eNeuro
Vol. 10, Issue 11
November 2023
Email
Print
View Full Page PDF
Citation Tools
Respond to this article
γ And β Band Oscillation in Working Memory Given Sequential or Concurrent Multiple Items: A Spiking Network Model
Shukuo Zhao, Jinpu Zhou, Yongwen Zhang, Da-Hui Wang
eNeuro 30 October 2023, 10 (11) ENEURO.0373-22.2023; DOI: 10.1523/ENEURO.0373-22.2023
Twitter logo Facebook logo Mendeley logo

Jump to section

Keywords

Responses to this article

Respond to this article

Jump to comment:

No eLetters have been published for this article.

Subjects