Nonlinear response characteristics of neural networks and single neurons undergoing optogenetic excitation

Over the last several years, the field of optogenetics has led to the development of extremely useful experimental tools that can be employed to stimulate single cells or entire neuronal populations. Essentially, optogenetics enables precisely targeted stimulation of specific cell typesby genetic modification and exposure to light. Optogenetic techniques have already facilitated many explorations, and likely many more will follow (Deisseroth, 2015).

An experimentally well-studied class of optogenetic tools are channelrhodopsins (Berndt et al., 2011; Lin, 2011; Schneider, Grimm, & Hegemann, 2015). Being a light-gated ion channel that leads to excitation of its host cell, wild-type channelrhodopsin-2 (ChR2) was the first microbial rhodopsin successfully employed as an optogenetic tool (Boyden, Zhang, Bamberg, Nagel, & Deisseroth, 2005; Nagel et al., 2003). Nowadays, a multitude of engineered variants of ChR2 exist. The work we present here is based on a model of the ChR2/H134R variant, which is popular for its enhanced photocurrents as well as for having a good peak/steady-state ratio and increased light sensitivity (Lin, 2011; Tchumatchenko, Newman, Fong, & Potter, 2013; Williams et al., 2013; Yawo, Asano, Sakai, & Ishizuka, 2013). The ChR2/H134R variant, however, is only an example of an optogenetic channel that can be used with our model. Parameter modifications and customization of the simulational tools which we provide will facilitate analogous network studies for other optogenetic variants.

Making quantitative activity predictions for neuronal networks in vivo is complicated by the interaction of the different timescales of rhodopsins, neurons, and synapses, as well by the voltage dependence of conductances. Determining the stimulus-response relationship in a network by simultaneously measuring the activities of thousands of neurons is experimentally challenging, and only recently pioneering experiments (Chettih & Harvey, 2019; Daie, Svoboda, & Druckmann, 2019; Russell et al., 2019) started tackling this question. Here, we present theoretical results that aid in the prediction of experimental outcomes as experiments push in the direction of understanding how populations of neurons respond to specific stimuli (Humphries, 2017). Hence, a reliable prediction of the network effects evoked by optogenetic excitation of a defined cell population would be an extremely helpful tool. Theoretical studies on the network response to pulsed current stimulation are not sufficient to describe optogenetic stimulation because the resonance properties of the rhodopsins have a complex effect on the properties of single neurons and whole networks. By including a detailed channel model in a network setting, however, optogenetically evoked network responses can be predicted more accurately. Proposing such a model, we aim to bridge the gap between channel dynamics and spiking, in addition to providing predictions for network and single-neuron dynamics following optogenetic stimulation.

We constructed our computational model by using a network of Leaky Integrate-and-Fire neurons, which are biologically realistic and have been used in a number of studies addressing the response of cortical populations (Fourcaud-Trocmé, Hansel, van Vreeswijk, & Brunel, 2003; Gerstner & Naud, 2009; Tchumatchenko, Malyshev, Wolf, & Volgushev, 2011). We augmented the Leaky Integrate-and-Fire neuron model by including in its input a population of ChR2 channels, simulated with a three-state Markov model previously investigated by others (Nikolic et al., 2009; Tchumatchenko et al., 2013; Williams et al., 2013). Given its few parameters and short runtime, the three-state model is a prime candidate for large network simulations. In addition, the model’s good response match to continuously varying stimuli (Tchumatchenko et al., 2013) made us choose the three-state model over a competing four-state model (Nikolic et al., 2009; Williams et al., 2013).

In the channel model we consider, each individual ChR2 channel can only be in one of the three states at a time. For the thousands of channels that are typically found in the membrane of ChR2-transfected neurons, a stochastic description via three continuous probability variables is sufficient. The expression level of ChR2 (i.e., the number of channels in a neuron) is difficult to determine in experiments, so that only rough estimations of its magnitude can be made. Here, we study two expression levels: a low expression level (60,000 channels per neuron) and a high expression level (300,000 channels per neuron). We chose these values because we found through simulations that 60,000 channels per neuron by themselves, without external or recurrent input, are not sufficient to trigger spikes, while the stimulation of a neuron with 300,000 channels reliably correlates with spiking. We estimated these numbers from conductance values in experiments, using them as representative values for relatively low and high expression of ChR2 (Nagel et al., 2003; Zhang, Wang, Boyden, & Deisseroth, 2006). The lower expression level yielded a small, but measurable, correlation between light pulses and spikes, and can therefore be considered a subthreshold modulation in in vivo networks. The higher level evoked spikes that were nearly phase-locked to the light pulses.

To evoke network responses, we chose a periodic stimulation with a frequency of 50 Hz, which is commonly used in experiments (Lin, Lin, Steinbach, & Tsien, 2009; Zhang et al., 2006), and we investigated the response across a broad range of light intensities. To monitor the impact of the connection probability on the network response, we also varied the connection probability across a broad range of values. Furthermore, we aimed to obtain a detailed characterization of the single neuron response by considering a large number of frequencies and light intensities.

If the network response were linear, we could expect a constant spatial profile at all stimulation amplitudes. We found, however, that the spatial distribution of activities in realistic, nonlinear neural networks is much broader than the spatial light stimulation, and that the width of this activity distribution depends on the light intensity, as well as on the number of channels expressed in neurons. Thus, our main finding is that a narrow stimulation profile evokes a broad response profile. We further explored the nonlinear dependence of the peak and the baseline activity as well as the mean activity in the network as a function of the light intensity.

Interestingly, for single neurons under pulsed light stimulation we discovered that the magnitude of the firing rate response minimum increases with higher stimulation frequency, while the absolute height of the evoked response pulses decreases. Moreover, we found the duration of the firing rate response pulses to exhibit a nonmonotonic relation to the stimulation frequency. With respect to the light intensity, we found a monotonic saturation of the firing rate response.

Our results show that the stimulus-response relationship of networks and single neurons exhibiting ChR2 is highly nonlinear but can, nevertheless, be described in a tractable manner. We provide a characterization of this nonlinear behavior, along with our customizable simulation tools, which can be used to predict the ChR2-mediated response of networks and of single neurons to pulsatile stimulation protocols. Our model and our results contribute a quantitative perspective and offer tractable in silico predictions that can be used to design experiments.

For fitting, we used gnuplot (version 5.0.3) and GNU Scientific Library (version 2.4). We also used gnuplot to create plots. We implemented all of our simulations in C++ (ISO 2011 standard) and compiled them with g++ in version 7.5.0. The source code and ready-to-run binary versions of our simulation tools can be retrieved from https://github.com/jlubo/nn-lightchannels-sim or from https://jlubo.net/nn-lightchannels-sim.

Our model assumes a channel to be in one of the following three states:

• 1. 

  closed (with probability C): no current is flowing, activation by effective photon flux ϵϕ(t) is possible;

• 2. 

  open (with probability O): current is flowing, desensitization will occur stochastically at rate Γ_d;

• 3. 

  desensitized (with probability D): no current is flowing, no activation is possible, recovery will occur stochastically at rat Γ_r.

A list of the most important parameters used in our simulations is given in Table 1.

The external input causes spiking dynamics to occur even in the absence of light stimulation. We adjusted the mean current I₀ for all simulations, such that unstimulated neurons were firing at a mean rate of 5 Hz. For single neurons that are not simulated in a network, the synaptic input equals the external input: I_syn(t) = I_ext(t).

In order to characterize the long-term activity dynamics of a single neuron, we examined three measures of the neuronal activity in the steady state. These are the response minimum, which is the lowest value the steady-state firing rate takes under persistent stimulation, the response maximum, which is the total height of the firing rate response pulses evoked by stimulus pulses, and the duration of those pulses, measured by the full width at half maximum. We averaged the pulses across the whole steady state and over 900 trials. To draw a comparison between the impact of the neuron dynamics and the channel dynamics, we determined the same three measures from the pulses in the steady-state open-state probability of a channel. Examples of the pulses in the steady-state open-state probability and the steady-state activity are shown in Figure 2C, D.

The onset phase of the light stimulation causes a strong but very brief spiking response lasting approximately 100 ms (see Figure 2A and B). This onset response vanishes quickly and its features are different from the steady-state response that follows. Therefore, to measure the experimentally more representative long-term spiking response, we considered the spikes which followed the onset phase and imposed a wait time in our simulations of approximately 100 ms (see Figure 2).

We used stimulation protocols consisting of temporally rectangular light pulses, delivered with different frequencies and light intensities. The duration of the pulses was kept constant at 4 ms, which is a value that has been used in experiments when employing moderately short light pulses (Boyden et al., 2005; Gunaydin et al., 2010). We used a sliding window/boxcar kernel approach to compute the neuronal activity from spike trains (Dayan & Abbott, 2001). This method is particularly useful because it prevents the resulting activity from being dependent on the placement of the time window, as would be the case with peristimulus time histogram spike densities.

Our network was represented by a square grid of excitatory neurons, and a population of inhibitory neurons. All neurons were coupled via random recurrent connections (Figure 3A). The probability of a connection p_c (i.e., the connectivity) was the same for all types of connections (E → E, E → I, I → E, I → I) and was varied across simulations.

As mentioned previously, the neurons receive input accounting for projections from outside the network, modeled by an Ornstein–Uhlenbeck process (Equation 8). In our network simulations, the synaptic current I_syn(t) further contains contributions from synapses within the network: I_syn(t) = I_ext(t) + ∑_iI_int,i(t).

First, we investigated how a single neuron containing a given number of ChR2 channels responds to the application of pulsed light stimuli. We used the Leaky integrate-and-fire model to describe the neuron. We selected this neuron model because it reproduces various properties of neocortical pyramidal cells (Gerstner & Naud, 2009; Tchumatchenko et al., 2011) and because of its outstanding numerical speed in network simulations. For the ChR2 channels, we used a probabilistic three-state Markov model, based on previous studies (Nikolic et al., 2009; Tchumatchenko et al., 2013).

As a first test for our model, we evaluated the impact of the expression level on the membrane potential dynamics. The results showed that at a low expression level of 60,000 channels, the number of channels was not sufficient to depolarize the neuron strongly enough to make it fire in synchrony with the frequency of light stimulation. At a higher expression level of 300,000 channels, synchronized firing occurred. Membrane potential traces for these two cases are shown in Figure 4.

Under sustained stimulation with a certain frequency and intensity, the dynamics of the open-state probability as well as the firing rate dynamics reached a steady state. In this steady state, we characterized the activity response by measuring the minimum firing rate, the pulse duration, and the maximum of the firing rate pulses in response to light pulses (cf. Figure 1B, Methods). As we varied the frequency and intensity of the light pulses, we recorded the response characteristics for a broad range of stimulus protocols.

We found that the activity response of a single neuron under optogenetic excitation exhibited characteristics that might be considered surprising or counterintuitive. In a pulsed stimulation protocol, the minimum of the evoked spiking activity ν_s,min increased with higher stimulation frequency (Figure 5A), while the maximum of the evoked activity pulses ν_s,max decreased (Figure 5B). Furthermore, it is remarkable that the duration of the firing rate pulses $tFWHMν$ exhibited a nonmonotonic relation in response to varying stimulation frequency (Figure 5C). In Figure 5 we show results for neurons with 300,000 channels while we found qualitatively similar behavior also for 60,000 channels per neuron.

One could expect that the activity response increases its minimum with increasing stimulation frequency, since a higher frequency should lead to more depolarization because it raises the level of continuous activity (Figure 5A). This is also what we found. Considering the dependence on the light intensity, there is a slight decrease in the activity response minimum for lower light intensity. The response minimum of the open-state probability showed similar behavior (Figure 6A).

The detailed behavior of the response maxima, however, is more intriguing. Our simulations show that as the stimulation frequency increases, the absolute height of the response pulses (i.e., the response maximum) decreases as it approaches the response minimum (Figure 5B). For very high stimulation frequencies, the response maximum and minimum have to be equal because above a certain frequency, which depends on the light pulse duration (250 Hz for 4 ms), pulsed stimulation would become in fact constant and the response pulses would vanish as well. Studying the dependence on the light intensity, we found a monotonic increase of the response maximum, leading into saturation. The response maximum of the open-state probability showed similar behavior (Figure 6B).

Studying the activity pulse duration (Figure 5C), we first found that it was longer than the duration of the stimulus pulses, which had a constant duration of 4 ms. In addition, we made a surprising discovery. The relationship of the activity pulse duration to the frequency was not monotonic, unlike the monotonic increase of the response minimum and decrease of the response maximum with increasing frequency. Instead, the pulse duration exhibited a maximum whose width and height depend on the stimulation frequency and light intensity. The comparison with the pulse duration of the open-state probability (Figure 6C), which does not exhibit such a maximum, suggests that the occurrence of the maximum is due to the restricted firing capability of neurons.

The comparison of the steady-state dynamics of the neuronal firing rate (Figure 5) and the open-state probability of a ChR2 channel (Figure 6) shows that the behavior of the response minimum and maximum of the activity seems to be determined predominantly by the channel dynamics. Since the pulses in the firing rate as well as in the open-state probability are longer than the stimulus pulses, some filtering has to occur already at the channel level. Nevertheless, the pulse duration of the firing rate shows nonmonotonic behavior in relation to the stimulus frequency, while the pulse duration of the open-state probability has a monotonic relationship to the stimulus frequency. Hence, the pulses of the firing rate reflect nonlinear neuronal processing. In fact, the pulse duration of the firing rate is smoothed out by the processes that occur in Leaky integrate-and-fire neurons (cf. Equations 6 and 7).

At the network level, we investigated the spatial distribution of activity following excitation through a spatially Gaussian-distributed light stimulus. As for the investigation of a single neuron, we used pulsed light stimuli at different light intensities, but now maintained their frequency at 50 Hz. Additionally, we spatially modulated the intensity with a Gaussian distribution. The center of the Gaussian spatial distribution conformed with the center of the square grid of excitatory neurons. We observed the resulting activity levels of the excitatory neurons following this stimulation paradigm. By averaging over the time-averaged activities of all individual neurons, we obtained the population activity 〈ν〉. To further analyze the network response, we sought a measure for the spatial distribution of evoked activity. To this end, we fitted different functions to the evoked activity and compared the R² values of these functional fits. We found that the Gaussian had the highest R² values (sample values for p_c = 1%, Ê = 5 mW/mm², and 60,000 channels: Gaussian 0.373, Lorentzian 0.363, Logistic distribution 0.368). Thus, we used a Gaussian fit function to describe the spatial distribution of the evoked activity (see Equation 11). We used the standard deviation σ_FR of this fit function to measure the width, the amplitude ν_max to measure the maximum, and the vertical shift ν_base to measure the baseline of the activity distribution. Light distributions, resulting activity distributions and Gaussian fit functions to the activity distributions are shown in Figures 7D and Figure 8D. These figures provide a glimpse on our main finding that the spatial width of the evoked activity distribution is considerably larger than the width of the light stimulus.

We found that the population activity ν〉 exhibited the same monotonically increasing behavior as a function of the maximum light intensity for both expression levels (Figures 7B, C and 8B, C). Increasing the connection probability p_c affected the magnitude but not the shape of the activity profile. We obtained a similar response for the height of the spatial activity distribution (ν_max). The baseline in the activity distributions (ν_base) declined as the peak light intensity increased, which could be explained by increased inhibition in the network evoked by the increased excitatory population activity. The width of the spatial activity distribution (σ_FR) increased monotonically with growing light intensity (Figures 7A and 8A). Varying the connection probability had a weak effect on the peak amplitude of the response, but not the shape of the response function. In the considered parameter space, the response width was between 1.5 and 2.5 times the width of the light stimulus (σ_light = 8 grid units). To summarize, our results revealed that for a large set of parameters, light stimulation evokes substantially broader firing activity responses.

Optogenetic manipulation of neural network activity has become a widely used method to modulate neuronal activity in vivo (Masseck, 2018). One of the most commonly employed optogenetic tools is still channelrhodopsin-2 (ChR2), a blue-light-gated cation-selective ion channel from a species of green algae (Chlamydomonas reinhardtii) which conducts H⁺, Na⁺, K⁺, and Ca²⁺ (Nagel et al., 2003). Historically, ChR2 was the first optogenetic channel used in neurons, and it has become the basis for other ChR variants (Masseck, 2018). Anticipating the effectiveness of an optogenetic stimulation and interpreting its outcome in vivo is often difficult because of the influences of channel activation, the single-neuron excitability, and the recurrent network dynamics. These effects interact with one another and can lead to a complex set of outcomes. Analyzing the stimulus-response relationship of neurons is of outstanding importance to theoreticians who aim to model and understand network dynamics following optogenetic stimulation, as well as for experimentalists who wish to alter neuronal activity in brain tissue in a desired way. Here, we present a spiking model and its analysis at the network level in order to help predict the outcome of excitatory optogenetic manipulation. Our computational model considers recurrent networks of excitatory and inhibitory neurons, in which excitatory neurons express light-sensitive ChR2/H134R channels.

To accurately represent the physiological properties of light-stimulated pyramidal neurons, we chose the Leaky integrate-and-fire (LIF) model. This choice is due to the observation that LIF neurons can capture the broad range of frequencies that are reliably encoded by cortical pyramidal neurons (0–200 Hz frequency range; Tchumatchenko et al., 2011), while other more detailed, conductance-based models, including the Hodgkin–Huxley models, fail to represent the dynamical response of cortical neurons (Fourcaud-Trocmé et al., 2003). Furthermore, the integrate-and-fire type models have been shown to accurately reproduce the experimentally recorded spike pattern (Gerstner & Naud, 2009) and thereby offer various computational advantages compared to more detailed conductance-based models. To model synapses, we used a monoexponential model, introducing postsynaptic currents with finite decay time. Because this model facilitates discretization in numerical computations, it is much faster than slightly more detailed alpha-function-based models. The LIF model parameters which we used correspond to the irregular asynchronous regime of cortical neurons: synaptic strength and time constant are chosen so that they account for the amplitude and width of postsynaptic potentials of AMPA, NMDA, and GABA synaptic currents (Gerstner & Kistler, 2002; London, Roth, Beeren, Häusser, & Latham, 2010; Roth & van Rossum, 2009). Moreover, the physiological parameters such as membrane time constant and targeted coefficient of variation range (0.5–1) were chosen to match experimental reports for the cortex (Stevens & Zador, 1998). The LIF model has, as well as in vivo neurons, two important parameter regimes that shape the spiking response of neurons: the noise-driven and the mean-driven regime (Petersen & Berg, 2016). In our study, neurons operate in the noise-driven regime that is typical for cortical neurons in vivo. Therefore, as long as the parameters stay within this biologically plausible regime, our prior work (Herfurth & Tchumatchenko, 2019; Tchumatchenko & Wolf, 2011) and the work of others (Brunel, Chance, Fourcaud, & Abbott, 2001; Gerstner & Kistler, 2002) suggests that one can expect similar results because the shape of the frequency response function is largely preserved across a broad range of membrane time constants, firing rates, and noise levels. To facilitate the exploration of other spiking models or parameter regimes in combination with optogenetic stimulation, we release our program code, so that our readers can consider the effects of any other custom neuron model on the recurrent network activity and, if necessary, adapt the parameters to any particular value of interest (on instructions see the file CUSTOM_MODELS.txt from the code linked in the Methods section). To model optogenetic light stimulation, we used the three-state model of ChR2 because it is computationally efficient, it describes the monoexponential photocurrent decay of ChR2/H134R (Williams et al., 2013), and it has been shown to reliably reproduce experimentally measured responses to continuously varying stimuli (Tchumatchenko et al., 2013). However, our program code can also be used to implement alternative ChR2 dynamics. Extensions of our model could include a fourth state that can account for biexponential decay of the photocurrent, which is exhibited by some ChR2 variants (Nikolic et al., 2009; Williams et al., 2013), light-assisted recovery, which has been observed for some ChR2 variants (Nagel et al., 2003; Nikolic et al., 2009), or a nonlinear voltage dependence factor to the photocurrent, which could account for specific inward rectification effects (Gradmann, Berndt, Schneider, & Hegemann, 2011; Grossman, Nikolic, Toumazou, & Degenaar, 2011; Lin et al., 2009).

We considered the network activity in response to an optogenetic excitation by spatially Gaussian-distributed light stimulation and studied its spatial profile. We found that the spatial extent, that is, the width, of the network activation can be 1.5 to 2.5 times as large as that of the light source. This indicates that in experiments targeting a particular spatial activation profile the light beam width should be chosen smaller than the intended area of activation, and that unintended co-activation of neighboring regions should be monitored. How small the light beam width will need to be depends on the light intensity and channelrhodopsin expression levels. Interestingly, we found that the response profile and its width depended on the expression levels, while the synaptic connectivity in the stimulated region modulated only the peak height but not the width of the response.

Here, we presented results for a number of basic measures of the activity distribution such as width, height, and baseline, while readers interested in more complex measures such as peak-to-width ratios (ν_max/σ_FR) can use use the program code we provide to explore additional quantities of interest.

To clarify the temporal dynamics of the network response we examined the firing rate dynamics at the single-neuron level. We found that under pulsatile light stimulation, which is commonly used in experiments, the neuronal activity response minimum increases with the stimulation frequency, while the activity response maximum decreases. Our results revealed the existence of local maxima in the stimulus-response relationship, which means that specific stimulus parameters can evoke particularly strong and particularly narrow responses. In particular for the pulse duration of the firing rate response we found that at low stimulation frequencies, the depolarization of neurons is not strong enough to produce sustained activation. This is due to the low-pass-filtering effect of neurons, which are not capable of immediate response to fast, single pulses. Furthermore, we found that the response dynamics saturate with higher light intensities, which is caused by the slow recovery of the channels from desensitization.

Although the number of theoretical studies on optogenetic effects is increasing, no study has yet presented a systematic investigation of the nonlinear effects that are subject to our work. Recently, interesting studies were published which use field-programmable gate array (FPGA) processors to simulate networks of small size with a very detailed neuron model (Luo et al., 2017) or provide a framework which enables simulating optogenetic impact on neurons and networks (Evans, Jarvis, Schultz, & Nikolic, 2016). So far, however, no explicit and quantitative predictions were derived for the optogenetically triggered spatial response of neuronal activity in large networks. Similarly, the effects occurring in the firing rate response of single neurons to light pulses of different frequencies have not been quantified. Providing a model and analysis tools, as well as reporting quantitative predictions that can help in the design and interpretation of optogenetic experiments in recurrent networks and single neurons, was the goal of our study. Since our software package is written in standard C++, it enables very fast computation and will in most cases outperform other tools that are based on higher-level languages such as Python, while offering maximal customizability.

What experimentally relevant effects did we find? To achieve controlled optogenetic stimulation in experiments in vivo, it is helpful to be able to reverse engineer the optogenetic excitation profile from the desired neural activation profile. In our study, we identified three important aspects that are relevant for this procedure. First, the size of the area where the stimulus evokes direct or indirect excitation can be substantially larger than the area covered by light (up to 2.5 times; cf. Figures 7 and 8). For experiments, it is therefore beneficial to start with a small area of excitation and broaden it only when there is evidence that not all intended neurons are reached. Starting with a large light stimulation area bears the risk of unintended co-activation of neighboring areas, since the effective footprint of light on the firing rates of neurons will be much broader than the light profile itself. In addition to the broadening of the spiking activity profile relative to light that we have characterized, further broadening of the area of response is to be expected by light scattering in brain tissue. While this investigation was beyond the scope of this study, it should also be kept in mind when setting up light stimulation protocols. Second, our results indicate that the duration of the firing rate pulses in response to light pulses can be substantially longer than the duration of each light pulse that evokes it. Our analysis has shown that evoking very brief firing rate responses is not possible with the ChR2 variant that we considered in our study (ChR2/H134R). Experiments which rely on millisecond or submillisecond precision of few carefully introduced extra spikes need to either drastically lower the light level such that the integrated firing rate response to light stimulation carries only a few spikes, or choose another optogenetic variant that may allow for a higher precision (Lin, 2011; Masseck, 2018). Third, our results show that the amplitude of the steady-state firing rate response is not only proportional to the light intensity, it is also highly sensitive to the frequency of the light stimulation. Such behavior has also been observed in experiments (Ni et al., 2017). If the experimental aim is to elicit strong excitation pulses in neurons, then it is best to use low stimulation frequencies. Otherwise, if the goal is to elicit maximal numbers of spikes that do not need to occur in pulses, then it is best to use intermediate stimulation frequencies of around 30 Hz, which provide an optimal balance between high response minimum and high response maximum. It is important to know that the response pulses of the spiking response in this stimulation range are long, an observation which needs to be considered when measuring the resonance properties of specific neural populations or investigating the synchronization properties in particular frequency ranges (e.g., Cardin et al., 2009). In summary, our study enables a quantitative prediction of neural network activity that can help guide and interpret the outcome of optogenetic experiments.

Jannik Luboeinski: Conceptualization; Data curation; Formal analysis; Investigation; Methodology; Software; Validation; Visualization; Writing - Original Draft; Writing - Review & Editing. Tatjana Tchumatchenko: Conceptualization; Formal analysis; Funding acquisition; Methodology; Project administration; Resources; Supervision; Validation; Writing - Original Draft; Writing - Review & Editing.

Tatjana Tchumatchenko, Max Planck Society. Tatjana Tchumatchenko, German Research foundation via SPP 2041.

We would like to thank Sara Konrad for advice on network dynamics, Rogier Poorthuis and Hiroshi Ito for scientific discussions, and Leonor Fajardo Rebelo for proofreading help.

Gaussian light beam:

A light beam whose transversal intensity profile follows a Gaussian function.

Optogenetics:

Umbrella term for techniques that make use of genetic modifications to introduce light sensitivity into neurons.

Ion channel dynamics:

Ion channels undergo transitions between different states, depending on factors like time, membrane voltage or light stimulation.

Channelrhodopsin-2:

A molecule complex that naturally occurs in unicellular algae. It can be genetically inserted into neurons and depolarize host neurons upon stimulation with blue light.

Markov model:

A Markov model assumes that the transition between two states depend only on the very last state and not on the whole history of states.

Nonlinear function:

A mathematical function that cannot be expressed by only using the operations of scalar multiplication and addition.

Light stimulation:

Optogenetically modified cells can be excited or inhibited by illuminating them with light of a certain wavelength.

Ornstein–Uhlenbeck process:

A stochastic process that has an exponentially decaying autocorrelation function and can therefore be used to model cortical noise.

Population activity:

The population activity or population rate is the average firing rate of the neurons in a population.

Neuronal firing rate:

The neuronal firing rate is the rate of spikes per time fired by a neuron and can be computed from spike trains, for example, by using a sliding time window.

Monotonically increasing function:

A function of a parameter that increases or stays constant, but never decreases, as the parameter increases.