<figure><img src="image_1.png" /></figure>

HAL
open science

# Theory of the asynchronous state of structured rotator networks and its application to recurrent networks of excitatory and inhibitory units

Jonas Ranft, Benjamin Lindner

<figure><img src="image_2.png" /></figure>

## To cite this version:

Jonas Ranft, Benjamin Lindner. Theory of the asynchronous state of structured rotator networks and its application to recurrent networks of excitatory and inhibitory units. 2022. hal-04289572

HAL Id: hal-04289572

https://cnrs.hal.science/hal-04289572v1

Preprint submitted on 16 Nov 2023

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_3.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
# Theory of the asynchronous state of structured rotator networks and its application to recurrent networks of excitatory and inhibitory units

Jonas Ranft*

Institut de Biologie de l'ENS, Ecole Normale Supérieure, CNRS,
Inserm, Université PSL, 46 rue d'Ulm, 75005 Paris, France

Benjamin Lindner

Bernstein Center for Computational Neuroscience, Berlin, Philippstraße 13, Haus 2, 10115 Berlin,
Germany and Physics Department of Humboldt University Berlin, Newtonstraße 15, 12489 Berlin, Germany

(Dated: November 22, 2022)

Recurrently coupled oscillators that are sufficiently heterogeneous and/or randomly coupled can show an asynchronous activity in which there are no significant correlations among the units of the network. The asynchronous state can nevertheless exhibit a rich temporal correlation statistics, that is generally difficult to capture theoretically. For randomly coupled rotator networks, it is possible to derive differential equations that determine the autocorrelation functions of the network noise and of the single elements in the network. So far, the theory has been restricted to statistically homogeneous networks, making it difficult to apply this framework to real-world networks, which are structured with respect to the properties of the single units and their connectivity. A particularly striking case are neural networks for which one has to distinguish between excitatory and inhibitory neurons, which drive their target neurons towards or away from firing threshold. To take into account network structures like that, here we extend the theory for rotator networks to the case of multiple populations. Specifically, we derive a system of differential equations that govern the self-consistent autocorrelation functions of the network fluctuations in the respective populations. We then apply this general theory to the special but important case of recurrent networks of excitatory and inhibitory units in the balanced case and compare our theory to numerical simulations. We inspect the effect of the network structure on the noise statistics by comparing our results to the case of an equivalent homogeneous network devoid of internal structure. Our results show that structured connectivity and heterogeneity of the oscillator type can both enhance or reduce the overall strength of the generated network noise and shape its temporal correlations.

## I. INTRODUCTION

Groups of interacting oscillators are capable of the prominent phenomenon of *synchronization* – famous text books have been devoted to the topic [1, 2]. A broad literature identified conditions under which oscillators rotate to some degree in synchrony and explored how this is affected by the distribution of eigenfrequencies, by disorder in the coupling coefficients, or by dynamical noise (see e.g. [3–7]). However, interacting oscillators can also be in an *asynchronous state* which is far from being a trivial counterpart of the synchronized case – it is both important in many applications as well as mathematically challenging to describe, and it is what we focus on in this article.

In neuroscience, networks of interest consist of excitatory and inhibitory neurons that are sparsely coupled to each other via chemical synapses [8, 9]. In theoretical studies, networks with sparse random recurrent connectivity display distinct states of activity, among them prominently the asynchronous state, in which the firing of one neuron is not or only very weakly correlated to the firing of other neurons [10–13]. Interestingly, this type of state is often observed in cortical networks of animals in the attentive state [14, 15].

In the asynchronous state of neural networks, the main statistics of interest is the mean activity and the temporal correlation of spikes within the spike train of a single cell –

as different neurons are not correlated by definition of this state, their cross-correlations can be neglected. We note that early theories of recurrent networks often focussed predominantly on the mean activity (the firing rate), employing a coarse Poisson (white-noise) approximation for the network-generated fluctuations seen by the single cell [16–18]. However, it has been shown that the generated network noise is white only in the limit case of very weak activity [19, 20] and displays in general pronounced correlations – in particular if realistic synaptic coupling strengths are taken into account [13, 21, 22].

What makes the second problem of the output’s autocorrelations mathematically tricky, is their self-consistence [19, 20]: The correlations in the output of one cell (output correlation) constitute the input-noise correlation for another cell that receives input from this very neuron. One may naively hope that the input correlations will wash out into a temporally uncorrelated (Poissonian) total input when we sum over the many statistically independent input spike trains arriving at the synapses of the considered neuron, but it is relatively easy to show that this is not the case [23]. Hence, when considering spike generation for a neuron embedded in the recurrent network, we do not know the statistics of the network noise that this neuron is subject to; we only know that this input noise should have the same correlation function, or equivalently the same power spectrum as the spike train (the output) generated by the very same neuron. This problem of determining self-consistently the correlation statistics in recurrent networks of neural oscillators has been pursued in three different ways in

<footnote>* jonas.ranft@ens.psl.eu</footnote>
the past.

First, we can turn the condition of self consistence into an iterative single-neuron simulation scheme [19, 20, 24, 25], in which we replace the recurrent input by a surrogate Gaussian noise. In contrast to the diffusion approximation, which employs white Gaussian noise [16–18], here temporal correlations of the input noise are taken into account by changing iteratively the power spectrum of the input noise according to the output spectrum of the previous iteration. This can be repeated until input and output spectra match to some degree of accuracy. The simplest version of this procedure works only for a homogeneous network [19, 20] but it can be generalized to heterogeneous setups as well [22]. This approach provides efficient numerical schemes to generate spike trains with power spectra that are in many cases indistinguishable from those measured in large recurrent networks. However, it does not give us a full theory of the self-consistent autocorrelations as the scheme ultimately requires stochastic simulations.

Secondly, it is possible to extend the classical diffusion approximation via the Fokker-Planck approach for recurrent networks of integrate-and-fire neurons. In a more accurate description of the network noise, additional variables via a Markovian embedding account for the self-consistent correlations of the input noise [26]. The white noise is here replaced by a multidimensional Ornstein-Uhlenbeck process the coefficients of which have to be found from the self-consistency condition for the autocorrelations of the spike train. The resulting high-dimensional Fokker-Planck equation is difficult to treat both analytically and numerically.

A third approach to the problem of the self-consistent correlations of neural oscillators is a radical simplification of the neuron model. Instead of an integrate-and-fire model as used in most of the aforementioned papers, one can e.g. use rate units as done in the pioneering study by Sompolinsky, Crisanti and Sommer [27] and derive a single differential equation giving access to the self-consistent autocorrelation function of the units (as comprehensively explained in [28]). This ansatz has more recently been generalized to take into account heterogeneity in the network [29–31].

Another approach employing a simplified model has been suggested by van Meegen and Lindner [32] who studied a randomly connected network of phase rotators and derived equations for the autocorrelation function of the network noise and of the single units. We have, quite recently, extended their theory to the case of additional (common or individual) dynamic noise [33]. In both studies, the model was restricted to have a homogeneous random connectivity, a setting in which it is not possible to strictly define excitatory and inhibitory units. Here we present another extension of the theory, namely, for networks comprising populations of phase rotators with distinct properties. We test the theory for the important special case of recurrent networks of excitatory and inhibitory units, a network of two populations that differ in their statistics of eigenfrequencies and connectivity. We discuss how the heterogeneity of the network affects the temporal correlations of the network noise and the correlations of the single units and compare specifically to the case of a homogeneous network.

Our paper is organized as follows. In the next section, we

introduce the model for a network of multiple populations of rotators and develop the general theory for the network noise correlation functions. In sec. III we apply this general theory to the case of EI networks in the balanced state and highlight the differences for the correlation statistics (correlation functions and power spectra) to the statistics that would be observed in a homogeneous network with statistically identical but randomized (unstructured) single-unit properties and connectivities. Finally, in sec. IV we summarize our results and discuss open problems for future research.

## II. GENERAL THEORY FOR MULTIPLE POPULATIONS

We consider the following dynamical equation for phase variables $\theta_m^\alpha$ of P distinct populations $\alpha \in \mathcal{P} = \{A, B, \dots\}$, $m \in \mathcal{N}_\alpha = \{1, \dots, N_\alpha\}$:

$$ \dot{\theta}_m^\alpha = \Omega_m^\alpha + \sum_{\beta \in \mathcal{P}} \sum_{n \in \mathcal{N}_\beta} K_{mn}^{\alpha\beta} F_{\alpha\beta}(\theta_n^\beta(t)) \quad (1) $$

where in all generality also the coupling functions $f_{\alpha\beta}$ may depend both on pre- and postsynaptic populations. This amounts to a total number of $N = \sum_\alpha N_\alpha$ equations for the $N$ neurons in the heterogeneous network.

We assume that all values of intrinsic frequencies $\Omega_m^\alpha$ and the coupling coefficients $K_{mn}^{\alpha\beta}$ are independently drawn for all $N$ units and $N(N-1)$ connections between distinct neurons, respectively, excluding self-connections. The heterogeneity of the network is captured by potentially distinct statistics of the underlying distributions, for each population $\alpha$ in the case of the intrinsic frequencies and for pairs of populations $\alpha\beta$ in the case of the coupling coefficients. Note that the block matrices $K_{mn}^{\alpha\beta}$ together contain the $N \times N$ recurrent connections among all neurons (self-connections being set to zero). We furthermore take into account potential heterogeneity in the coupling functions $F_{\alpha\beta}$ between (not necessarily distinct) populations $\beta$ and $\alpha$.

Specifically, we consider the intrinsic frequencies to be Gaussian distributed with respective mean values $\langle \Omega_m^\alpha \rangle = \Omega_0^\alpha$ and variances $\langle (\Omega_m^\alpha - \Omega_0^\alpha)^2 \rangle = \tilde{\sigma}_\alpha^2$. For the theory, we do not have to assume a particular distribution for the coupling coefficients $K_{mn}^{\alpha\beta}$ but characterize them exclusively by their first- and second-order statistics:

$$ \langle K_{mn}^{\alpha\beta} \rangle = \kappa_1^{\alpha\beta} \quad (2a) $$

$$ \langle K_{mn}^{\alpha\beta} K_{op}^{\gamma\delta} \rangle = \kappa_1^{\alpha\beta} \kappa_1^{\gamma\delta} + \kappa_2^{\alpha\beta} \delta_{\alpha\gamma} \delta_{\beta\delta} \delta_{mo} \delta_{np} \quad (2b) $$

$$ \langle (K_{mn}^{\alpha\beta})^2 \rangle = (\kappa_1^{\alpha\beta})^2 + (\kappa_2^{\alpha\beta})^2 \quad (2c) $$

Here and in the following, the angular brackets represent averages over different realizations of the network, the intrinsic frequencies, and initial conditions. Specific choices for the coupling matrices are introduced and discussed further below.

In the following, we will assume that the presynaptic phases are uniformly distributed, as can be expected for an asynchronous state. We can then introduce shifted intrinsic frequencies $\omega_m^\alpha$ and coupling functions $f_{\alpha\beta}$ to absorb any finite
mean input. With

$$F^{\alpha\beta}(\theta) = \sum_{l} A_{l}^{\alpha\beta} e^{il\theta} = A_{0}^{\alpha\beta} + f_{\alpha\beta}(\theta) \quad (3a)$$

$$f_{\alpha\beta}(\theta) = \sum_{l \neq 0} A_{l}^{\alpha\beta} e^{il\theta} \quad (3b)$$

we obtain

$$\dot{\theta}_{m}^{\alpha} = \omega_{m}^{\alpha} + \sum_{\beta \in \mathcal{P}} \sum_{n \in \mathcal{N}_{\beta}} K_{mn}^{\alpha\beta} f_{\alpha\beta}(\theta_{n}^{\beta}(t)), \quad (4)$$

where

$$\langle \omega_m^\alpha \rangle = \Omega_0^\alpha + \sum_{\beta \in \mathcal{P}} \sum_{n \in \mathcal{N}_\beta} \kappa_1^{\alpha\beta} A_0^{\alpha\beta} \equiv \omega_0^\alpha \quad (5a)$$

$$\langle (\omega_m^\alpha - \omega_0^\alpha)^2 \rangle = \tilde{\sigma}_\alpha^2 + \sum_{\beta \in \mathcal{P}} \sum_{n \in \mathcal{N}_\beta} \kappa_2^{\alpha\beta} (A_0^{\alpha\beta})^2 \equiv \sigma_\alpha^2 \quad (5b)$$

Of particular importance is the so-called balanced state in which the recurrent input to the single unit does not diverge in the limit of an infinite network []. A simple choice for such a scenario is a vanishing mean value. From the above it can be seen that the recurrent input vanishes when

$$\sum_{\beta \in \mathcal{P}} \sum_{n \in \mathcal{N}_\beta} \kappa_1^{\alpha\beta} A_0^{\alpha\beta} = \sum_{\beta \in \mathcal{P}} N_\beta \kappa_1^{\alpha\beta} A_0^{\alpha\beta} \equiv 0. \quad (6)$$

These P equations impose constraints on the parameters $N_\beta$, $\kappa_1^{\alpha\beta}$, and $A_0^{\alpha\beta}$. A typical choice is to scale the $\kappa_1^{\alpha\beta}$ with inverse powers of $N_\beta$; the above constraints allow to have a finite (in fact vanishing) mean value for infinite populations even when the individual terms $N_\beta \kappa_1^{\alpha\beta} A_0^{\alpha\beta}$ would diverge [].

In order to develop a self-consistent theory for coupled populations of neurons, we rewrite the recurrent input and the governing equations for the phases as follows:

$$\dot{\theta}_{m}^{\alpha} = \omega_{m}^{\alpha} + \xi_{m}^{\alpha} \quad (7a)$$

$$\xi_m^\alpha = \sum_{\beta \in \mathcal{P}} \xi_m^{\alpha\beta}, \quad \xi_m^{\alpha\beta} = \sum_{n \in \mathcal{N}_\beta} K_{mn}^{\alpha\beta} f_{\alpha\beta}(\theta_n^\beta(t)) \quad (7b)$$

In line with previous theories [27, 32] (but see [33]), we will now assume that in the asynchronous state the $\xi_m^{\alpha\beta}$ can be approximated as stochastic noise processes with Gaussian statistics and given temporal correlations. Exactly these statistics, captured by the correlation functions of the $\xi_m^{\alpha\beta}$, we want to determine self-consistently, which eventually will allow us to describe the statistics of the network without solving the full (deterministic) N equations of the network. The following derivations, shown for completeness and to make our presentation self-contained, are similar in spirit to what was done in refs. [32, 33] and we refer the interested reader to these works for more detailed discussions of the steps involved.

We first note that the auto-correlation function $C_\xi^\alpha(\tau) = \langle \xi_m^\alpha(t+\tau) \xi_m^\alpha(t) \rangle$ of the combined input $\xi_m^\alpha$ to population $\alpha$ is given by a simple sum over the auto-correlation functions

$C_\xi^{\alpha\beta}(\tau) = \langle \xi_m^{\alpha\beta}(t+\tau) \xi_m^{\alpha\beta}(t) \rangle$ of the individual inputs $\xi_m^{\alpha\beta}$. According to Eq. (7b),

$$C_{\xi}^{\alpha}(\tau) = \sum_{\beta \in \mathcal{P}} \sum_{\gamma \in \mathcal{P}} \langle \xi_{m}^{\alpha\beta}(t+\tau) \xi_{m}^{\alpha\gamma}(t) \rangle \\ = \sum_{\beta \in \mathcal{P}} \langle \xi_{m}^{\alpha\beta}(t+\tau) \xi_{m}^{\alpha\beta}(t) \rangle = \sum_{\beta \in \mathcal{P}} C_{\xi}^{\alpha\beta}(\tau); \quad (8)$$

in the second line, we used additionally that the sum over $\gamma$ reduces to one term because for $\beta \neq \gamma$

$$\begin{aligned} & \langle \xi_m^{\alpha\beta}(t+\tau) \xi_m^{\alpha\gamma}(t) \rangle \\ &= \sum_{n \in \mathcal{N}_\beta} \sum_{n' \in \mathcal{N}_\gamma} \langle K_{mn}^{\alpha\beta} K_{mn'}^{\alpha\gamma} f_{\alpha\beta}(\theta_n^\beta(t+\tau)) f_{\alpha\gamma}(\theta_{n'}^\gamma(t)) \rangle \\ &= \sum_{n \in \mathcal{N}_\beta} \sum_{n' \in \mathcal{N}_\gamma} \langle K_{mn}^{\alpha\beta} K_{mn'}^{\alpha\gamma} \rangle \langle f_{\alpha\beta}(\theta_n^\beta(t+\tau)) \rangle \langle f_{\alpha\gamma}(\theta_{n'}^\gamma(t)) \rangle \\ &= \sum_{n \in \mathcal{N}_\beta} \sum_{n' \in \mathcal{N}_\gamma} \kappa_1^{\alpha\beta} \kappa_1^{\gamma\delta} \cdot 0 \cdot 0 = 0, \end{aligned} \quad (9)$$

where $\langle f_{\alpha\beta}(\theta_n^\beta) \rangle = 0$ follows from averaging over random, uniformly distributed $\theta_n^\beta$ (note that $f_{\alpha\beta}(\theta) = \sum_{l \neq 0} A_l^{\alpha\beta} e^{il\theta}$).

Turning to the correlation functions $C_\xi^{\alpha\beta}(\tau)$, we find

$$\begin{aligned} & \langle \xi_m^{\alpha\beta}(t+\tau) \xi_m^{\alpha\beta}(t) \rangle \\ &= \sum_{n \in \mathcal{N}_\beta} \sum_{n' \in \mathcal{N}_\beta} \langle K_{mn}^{\alpha\beta} K_{mn'}^{\alpha\beta} f_{\alpha\beta}(\theta_n^\beta(t+\tau)) f_{\alpha\beta}(\theta_{n'}^\beta(t)) \rangle \\ &= \sum_{n \in \mathcal{N}_\beta} \sum_{n' \in \mathcal{N}_\beta} \langle K_{mn}^{\alpha\beta} K_{mn'}^{\alpha\beta} \rangle \langle f_{\alpha\beta}(\theta_n^\beta(t+\tau)) f_{\alpha\beta}(\theta_{n'}^\beta(t)) \rangle \\ &= \sum_{n \in \mathcal{N}_\beta} \langle K_{mn}^{\alpha\beta} \rangle^2 \langle f_{\alpha\beta}(\theta_n^\beta(t+\tau)) f_{\alpha\beta}(\theta_n^\beta(t)) \rangle \\ &\stackrel{(2c)}{=} N_\beta [(\kappa_1^{\alpha\beta})^2 + (\kappa_2^{\alpha\beta})] \langle f_{\alpha\beta}(\theta_n^\beta(t+\tau)) f_{\alpha\beta}(\theta_n^\beta(t)) \rangle. \end{aligned} \quad (10)$$

As above, terms with $n' \neq n$ drop out, assuming that $\theta_n^\beta(t_1)$, $\theta_{n'}^\beta(t_2)$ are uncorrelated since we consider an asynchronous state.

Using the Fourier representation of $f_{\alpha\beta}$ and expressing the phases as integrals of the respective inputs,

$$\begin{aligned} \theta_n^\beta(t) &= \theta_n^\beta(t_0) + \int_{t_0}^t dt' \dot{\theta}_n^\beta(t') \\ &= \theta_n^\beta(t_0) + \omega_n^\beta(t-t_0) + \int_{t_0}^t dt' \xi_n^\beta(t'), \end{aligned} \quad (11)$$
one eventually obtains for the remaining average in Eq. (10)

$$
\begin{align}
& \langle f_{\alpha\beta}(\theta_n^\beta(t+\tau)) f_{\alpha\beta}(\theta_n^\beta(t)) \rangle \nonumber \\
&= \sum_{l \neq 0} \sum_{l' \neq 0} A_l^{\alpha\beta} A_{l'}^{\alpha\beta} \langle e^{il\theta_n^\beta(t+\tau)} e^{il'\theta_n^\beta(t)} \rangle \nonumber \\
&= \sum_{l \neq 0} \sum_{l' \neq 0} A_l^{\alpha\beta} A_{l'}^{\alpha\beta} \underbrace{\langle e^{il(l+l')\theta_n^\beta(t_0)} \rangle}_{\delta_{l,-l'}} \langle e^{il(l_0^{l+\tau}t)\theta_n^\beta(t')} e^{il(l_0^{l+\tau}t')\theta_n^\beta(t')} \rangle \nonumber \\
&= \sum_{l \neq 0} |A_l^{\alpha\beta}|^2 \langle e^{il\omega_n^\beta l\tau} \rangle \langle e^{il(l_0^{l+\tau}t)\tau'\xi_n^\beta(t')} \rangle \nonumber \\
&= \sum_{l \neq 0} |A_l^{\alpha\beta}|^2 \Phi_\beta(l\tau) e^{-\frac{l^2}{2} \langle |\xi_n^\beta(\tau)|^2 \rangle} \tag{12}
\end{align}
$$

In the last step, we introduced the characteristic function of
the distribution over the (effective) intrinsic frequencies $\omega_n^\beta$
of population $\beta$,

$$
\Phi_{\beta}(l\tau) = \langle e^{i\omega_{n}^{\beta}l\tau} \rangle, \quad (13)
$$

and assumed that the integrated recurrent input (or integrated
network noise)

$$
y_n^\beta(\tau; t) = \int_t^{l+\tau} dt' \xi_n^\beta(t') \quad (14)
$$

is Gaussian-distributed with zero mean, which implies
$\langle e^{ily_n^\beta(\tau)} \rangle = e^{-\frac{l^2}{2} \langle |y_n^\beta(\tau)|^2 \rangle}.$

Using furthermore

$$
\begin{align}
& \int_t^{l+\tau} dt' \int_t^{l+\tau} dt'' \langle \xi_n^\beta(t') \xi_n^\beta(t'') \rangle \nonumber \\
&= 2 \int_0^\tau dt (\tau - t) \langle \xi_n^\beta(t+\tau) \xi_n^\beta(t) \rangle \nonumber \\
&= 2 \int_0^\tau dt (\tau - t) C_\xi^\beta(t) \tag{15}
\end{align}
$$

and the definition

$$
\Lambda_{\alpha}(\tau) = \int_{0}^{\tau} dt (\tau - t) C_{\xi}^{\alpha}(t) \quad (16)
$$

we eventually get the following system of ordinary differential
equations for the $\Lambda_\alpha(\tau)$:

$$
\ddot{\Lambda}_{\alpha}(\tau) = \sum_{\beta \in \mathcal{P}} N_{\beta} [(\kappa_{1}^{\alpha\beta})^2 + \kappa_{2}^{\alpha\beta}] \sum_{l \neq 0} |A_{l}^{\alpha\beta}|^2 \Phi_{\beta}(l\tau) e^{-l^2 \Lambda_{\beta}(\tau)}, \quad (17)
$$

the solutions of which give the network-noise autocorrelation
according to

$$
C_{\xi}^{\alpha}(\tau) = \ddot{\Lambda}_{\alpha}(\tau). \qquad (18)
$$

Eq. (17) is the key result of our paper. We have reduced
the *N* equations for all the neurons in the nework to *P* cou-
pled second-order ordinary differential equations for the func-
tions Λ<sub>α</sub>(τ) that provide the self-consistent network statistics.
In the previously considered case of a homogeneous network

with $P = 1$, the set of equations reduces to a single one that
has been derived by van Meegen and Lindner [32]. The same
reduction applies when all statistics for the different subpop-
ulations are identical, i.e., all parameters become independent
of $\alpha$ and $\beta$ ($\kappa_1^{\alpha\beta} = \kappa_1, ...$).

Besides the correlation function of the network noise, we
can also determine the correlation statistics of the single rota-
tors. Using the phase evolution given by Eq. (7a), it is rel-
atively straightforward to calculate the correlation function
for the pointer $x_m^\alpha(t) = e^{i\theta_m^\alpha(t)}$ for rotator *m* of population *α*
with effective intrinsic frequency $\omega_m^\alpha$. (Note that the dynamics
$\dot{\theta} = \omega + \xi$, where $\xi$ is a noise, corresponds to the Kubo oscil-
lator [34], with the added complication that the noise is deter-
mined self-consistently here.) Along the lines of the deriva-
tion of $C_\xi^\alpha$, one eventually obtains

$$
C_{x_m^\alpha}(\tau) = e^{i\omega_m^\alpha\tau - \Lambda_\alpha(\tau)}. \quad (19)
$$

Often, it makes sense to pool the data of the rotators of one
type. This corresponds to an average over the autocorrelation
functions of the pointers $x_m^\alpha(t)$ of the specific subpopulation.
In the limit of large $N_\alpha$, this approaches an average over the
frequencies $\omega_m^\alpha$,

$$
C_x^\alpha(\tau) = \langle e^{i\omega_m^\alpha \tau - \Lambda_\alpha(\tau)} \rangle_{\omega^\alpha} = \Phi_\alpha(\tau) e^{-\Lambda_\alpha(\tau)}, \quad (20)
$$

where $\Phi_\alpha(\tau)$ is the characteristic function of the effective in-
trinsic frequency distribution of population $\alpha$. We note that
this averaged correlation function is also proportional to the
autocorrelation function of a subpopulation activity

$$
X^{\alpha}(t) = \sum_{n \in \mathcal{N}_{\alpha}} x_{n}^{\alpha}(t), \qquad (21)
$$

because for uncorrelated rotators, one finds readily

$$
C_X^\alpha(\tau) = \sum_{n \in \mathcal{N}_\alpha} C_{x_n}^\alpha(\tau) \approx N_\alpha C_X^\alpha(\tau). \quad (22)
$$

Here, we used again that for large $N_\alpha$ the sampling over distinct rotators corresponds to an ensemble average over the effective frequencies.

Instead of correlation functions with pronounced oscilla-
tory components, it can be instructive to consider power spec-
tra of the observables. These are here obtained by numer-
ical Fourier transformation of the autocorrelation functions,
according to the Wiener-Khinchin theorem [34], $S(\omega) =$
$2\Re \int_0^\infty d\tau C(\tau) e^{-i\omega\tau}$.

### III. APPLICATION TO AN EXCITATORY-INHIBITORY NETWORK IN THE BALANCED STATE

As an example application, we consider an excitatory-
inhibitory (E-I) network comprising two populations $\mathcal{P} =$
$\{E, I\}$. Units are connected with probability $p$ and—if
connected—possess a fixed synaptic connection strength $j_{\alpha\beta}$
that depends on pre- and postsynaptic neuron identity, where
$j_{\alpha E} > 0$ and $j_{\alpha I} < 0$ to reflect the excitatory and inhibitory
character of populations E and I, respectively. For unconnected neurons, the coupling coefficients are set to zero. We can then give the network statistics introduced above in terms of the $j_{\alpha\beta}$ and $p$:

$$ \kappa_1^{\alpha\beta} = p j_{\alpha\beta} \qquad (23a) $$

$$ \kappa_2^{\alpha\beta} = p(1-p)j_{\alpha\beta}^2 \qquad (23b) $$

$$ \langle (K_{mn}^{\alpha\beta})^2 \rangle = p j_{\alpha\beta}^2 \qquad (23c) $$

The self-consistent equations (17) for the network-noise statistics thus become

$$ \ddot{\Lambda}_{\alpha}(\tau) = \sum_{\beta \in \mathcal{P}} p N_{\beta} j_{\alpha\beta}^2 \sum_{l \neq 0} |A_l^{\alpha\beta}|^2 \Phi_{\beta}(l\tau) e^{-l^2 \Lambda_{\beta}(\tau)}. \quad (24) $$

We remind the reader that the characteristic functions $\Phi_{\beta}$ reflect the effective distributions of the intrinsic frequencies taking into account the mean recurrent input, see the discussion surrounding Eqs. (5). In particular, we see that the sign of the synaptic coupling does not influence the dynamics directly as only $j_{\alpha\beta}^2$ enters Eq. (24). To distinguish excitation and inhibition in the framework of our model, we thus need coupling functions $f_{\alpha\beta}$ with non-vanishing Fourier components $A_0^{\alpha\beta}$, giving rise to non-vanishing mean synaptic inputs. For the effective intrinsic frequencies, we consider in the following without loss of generality that $A_0^{\alpha\beta} = 1$.

### A. Simplification of the two-population theory resulting from the balance condition

As can be seen from the products $pN_{\beta}j_{\alpha\beta}^2$ appearing in each term of Eq. (24), the network noise can be expected to be finite and non-zero in the limit of large N if the synaptic strengths $j_{\alpha\beta}$ scale with the average number of connections according to

$$ j_{\alpha\beta} = \frac{J_{\alpha\beta}}{\sqrt{pN_{\beta}}}. \qquad (25) $$

We then have

$$ \omega_0^{\alpha} = \Omega_0^{\alpha} + \sqrt{pN_E} J_{\alpha E} + \sqrt{pN_I} J_{\alpha I} \qquad (26) $$

for the mean effective intrinsic frequency of population $\alpha$. In order to ensure that the intrinsic frequencies remain finite in the limit of large N and constant p, we choose the case of balanced input, i.e.,

$$ J_{\alpha I} = -J_{\alpha E} \sqrt{N_E/N_I}. \qquad (27) $$

The spread of the intrinsic frequencies due to the random connectivity remains finite and is in this case given by

$$ \sigma_{\alpha}^2 = \tilde{\sigma}_{\alpha}^2 + (1-p)J_{\alpha E}^2(1+N_E/N_I). \qquad (28) $$

Note that the effective spread of the frequencies depends on our choice of the connectivity, here of the Erdos-Rényi type.

If for instance we chose a balanced networks with a fixed in-degree, the intrinsic frequencies would not be affected by the (otherwise random) recurrent connections, in contrast to Eq. (28).

Let us furthermore assume that the coupling functions only depend on the presynaptic identity, i.e.,

$$ f_{\alpha\beta}(\theta) \equiv f_{\beta}(\theta) \text{ or } A_l^{\alpha E} \equiv A_l^E, \quad A_l^{\alpha I} \equiv A_l^I \qquad (29) $$

For the balanced network, we can then further simplify the equations for $\Lambda_{\alpha}$ with the substitutions

$$ \Lambda_{\alpha}(\tau) = J_{\alpha E}^2 \lambda(\tau), \qquad (30) $$

giving rise to a single ordinary differential equation for $\lambda(\tau)$:

$$ \begin{aligned} \ddot{\lambda}(\tau) &= \frac{\ddot{\Lambda}_{\alpha}(\tau)}{J_{\alpha E}^2} = \sum_{\beta \in \mathcal{P}} \frac{J_{\alpha\beta}^2}{J_{\alpha E}^2} \sum_{l \neq 0} |A_l^{\alpha\beta}|^2 \Phi_{\beta}(l\tau) e^{-l^2 \Lambda_{\beta}(\tau)} \\ &= \sum_{l \neq 0} |A_l^E|^2 \Phi_E(l\tau) e^{-l^2 J_{EE}^2 \lambda(\tau)} + \frac{N_E}{N_I} \sum_{l \neq 0} |A_l^I|^2 \Phi_I(l\tau) e^{-l^2 J_{II}^2 \lambda(\tau)}. \end{aligned} \qquad (31) $$

Remarkably, this means that for a balanced network, the theory for a two-population system is almost as simple as the theory for a single homogeneous population. Notably, the correlation functions $C_{\xi}^{E}(\tau)$ and $C_{\xi}^{I}(\tau)$ of the network input to populations E and I, respectively, are simply scaled versions of each other.

### B. Comparison to an unstructured network and the corresponding single-population theory

The theory for the balanced E-I network derived above takes into account the particular structure of the network, which is comprised of two distinct populations with specific input noise statistics. We can now ask how well such a network would be described by the original theory of van Meeugen and Lindner developed for networks devoid of internal structure, i.e. for a homogeneous network. This begs the question of how such a comparison can be meaningfully conceived. We consider here for comparison a randomized version of the original, structured network where all connection strengths $K_{mn}$, with $m, n = 1, ..., N$, are drawn from the combined and properly weighted distributions of the $K_{mn}^{\alpha\beta}$, with $m = 1, ..., N_{\alpha}$, $n = 1, ..., N_{\beta}$ and $\alpha, \beta \in \mathcal{P}$. Since in the theory for homogeneous networks the coupling functions between units were assumed to be the same for all connections, we will restrict the comparison further to the case where $f_E(\theta) = f_I(\theta) \equiv f(\theta)$ for all $\theta$.

To apply the single-population theory and assess the effect of network structure by comparison, we thus consider the effective “single-population” network

$$ \dot{\theta}_m = \omega_m + \xi_m \qquad (32a) $$

$$ \xi_m = \sum_n K_{mn} f(\theta_m) \qquad (32b) $$
where *m*, *n* are now indices that run over all $N = N_E + N_I$ units. We consider the effective distribution of the intrinsic frequencies of the homogeneous network to be given by $p_{\text{eff}}(\omega) = \frac{N_E}{N} \mathcal{N}(\omega_0^E, \sigma_E) + \frac{N_I}{N} \mathcal{N}(\omega_0^I, \sigma_I)$, that is, by the weighted combination of the frequency distributions of the original populations *E* and *I*, where we take into account the additional broadening due to the recurrent input expected for each population, see Eq. (28). Note that one could also consider instead the broadening that would be expected from the unstructured recurrent input $\sum_n K_{mn}F(\theta_n)$ on rotators with intrinsic frequencies $\Omega_m$; however, we choose here to restrict our comparison to the effects related to the self-generated network noise.

For an unstructured network, the network-noise correlation function $C_{\xi}(\tau) = \langle \xi_m(t + \tau) \xi_m(t) \rangle$ is given by

$$C_{\xi}(\tau) = \ddot{\Lambda}(\tau) \qquad (33a)$$

$$\ddot{\Lambda}(\tau) = K^2 \sum_{l \neq 0} |A_l|^2 \Phi(l\tau) e^{-l^2 \Lambda(\tau)} \qquad (33b)$$

according to van Meegen and Lindner [32]. Here, $K^2$ is related to the variance of the connectivity matrix $K_{mn}$ that contains the recurrent synaptic weights of all *N* units, and $\Phi(t)$ is the characteristic function of the distribution of effective intrinsic frequencies for both populations combined. With the synaptic weights between populations given above, one finds for the variance of the connectivity matrix

$$\langle K_{mn}^2 \rangle = \frac{1}{N} \left( J_{EE}^2 \frac{N_E}{N_I} + J_{IE}^2 \right) . \qquad (34)$$

Our approximation of a two-population network would thus become

$$\ddot{\Lambda}(\tau) = K^2 \sum_{l \neq 0} |A_l|^2 \left( \frac{N_E}{N} \Phi_E(l\tau) + \frac{N_I}{N} \Phi_I(l\tau) \right) e^{-l^2 \Lambda(\tau)}, \qquad (35)$$

with $K^2 = J_{EE}^2 N_E / N_I + J_{IE}^2$.

From the solution of Eq. (35), one obtains the correlation function of the pointer $x_m(t) = e^{i\theta_m(t)}$ of rotator *m* with effective intrinsic frequency $\omega_m$ according to

$$C_{x_m}(\tau) = e^{i\omega_m \tau - \Lambda(\tau)}. \qquad (36)$$

(Note that the network-noise autocorrelation function $C_{\xi}(t)$ is directly given by Eq. (33a).) In a strict sense, the theory for an unstructured network does not permit to consider correlation functions at the level of distinct populations. To compare the population-averaged rotator autocorrelation functions predicted by the two-population theory to the 'unstructured' network, we can, however, average over the effective intrinsic frequency distributions of populations *E* and *I*, respectively, such that

$$C_{x}^{E,I}(\tau) = \Phi_{E,I}(\tau) e^{-\Lambda(\tau)}. \qquad (37)$$

Thus comparing the two-population theory to the equivalent 'unstructured' network allows one to assess the importance of having a correct description of the self-generated network noise at the level of the individual populations, whereas the

first-order effect due to the spread of intrinsic frequencies is taken into account also in the one-population theory.

Let us finally consider the specific case where the connection strengths $J_{\alpha\beta}$ do not depend on the postsynaptic identity, i.e., $J_{EE} = J_{IE} \equiv J$ (and consequently $J_{EI} = J_{II} = -J\sqrt{N_E/N_I}$). For the structured network, with $\Lambda(\tau) = J^2\lambda(\tau)$, one has

$$C_{\xi}^{E}(\tau) = C_{\xi}^{I}(\tau) = \ddot{\Lambda}(\tau) \qquad (38)$$

and with $K^2 = J^2 N / N_I$ the two-population theory reduces to

$$\ddot{\Lambda}(\tau) = K^2 \sum_{l \neq 0} |A_l|^2 \left( \frac{N_I}{N} \Phi_E(l\tau) + \frac{N_E}{N} \Phi_I(l\tau) \right) e^{-l^2 \Lambda(\tau)} \qquad (39)$$

which indeed is very similar to the single-population theory Eq. (35) but interestingly does not exactly correspond to that result. This difference can be understood as follows. Whereas in the unstructured network (Eq. (35)), the characteristic functions of each population are weighted with the population size, they contribute in the actual two-population network (Eq. (39)) with a weight that scales with the synaptic strength. For a typical value of $g \equiv N_E / N_I = 4$, the balance condition implies that albeit the inhibitory population is smaller in size, its characteristic function influences the dynamics of the self-generated noise more strongly.

### C. Comparison to network simulations

Even though we consider the balanced case of an E-I network, we still have some freedom in choosing the synaptic coupling strength for the distinct connections. In the following, we separately study three such choices that differ in the overall scaling of the synaptic inputs to the inhibitory population compared to the inputs to the excitatory population.

#### 1. Strong input to the inhibitory population

We start with the case illustrated in Fig. 1 on the left: both excitatory and inhibitory inputs to the inhibitory population are chosen to be very strong compared to the inputs to the excitatory population, i.e. $J_{IE} > J_{EE}$ and $|J_{II}| > |J_{EI}|$, where the latter follows directly from the first with the balance condition (27). First of all, we remark that all measured correlation functions and power spectra (blue and orange for E and I populations, respectively) are well described by our two-population theory (dashed and dotted black lines), see Fig. 1a-f.

As expected for our choice of input amplitudes, the network noise for neurons in the I population (orange) is stronger than for neurons in the E population (blue), indicated by a larger network-noise correlation function (Fig. 1a) and power spectrum (Fig. 1b); in fact, according to our theory and as confirmed by our simulations, the correlation function of the network noise for the I units is just an upscaled version of that for the E units. In particular, this implies that $C_{\xi}^{E}(\tau)$ and $C_{\xi}^{I}(\tau)$
<figure><img src="image_4.png" /><figcaption>Figure 1. Network statistics for a balanced E-I network with strong input to the inhibitory population. Left: sketch of the network with excitatory and inhibitory populations and corresponding connection strengths. Right: measured correlation functions (a,c,e) and power spectra (b,d,f) for network noise (a,b), single rotators (c,d) and population-averaged rotators (e,f). Network parameters: $N_E = 800$, $N_I = 200$, $p = 0.2$, $J_{EE} = 0.5$, $J_{IE} = 2$, $\Omega_0^E = 1$, $\Omega_0^I = 3$, $\bar{\sigma}_E = \bar{\sigma}_I = 0$; simulation parameters: time step $dt = 0.01$, simulation length $T = n_0 T_0$ with time window $T_0 = 1000$ and $n_0 = 10$, repeated for $n_r = 12$ network realizations with newly drawn connectivity matrices; numerical integration parameters: time step $dt = 0.01$, $T = 1000$.</figcaption></figure>

share the same time dependence. We can also compare the network noise correlation functions to the prediction of the single-population theory (grey line), which does *not* share this time course. We note that in the considered case the single-population theory provides correlation functions and power spectra that are in between the corresponding statistics for the E and I populations. This seems to be plausible but has not necessarily to be the case (see below).

Turning to single-rotator statistics, we show the correlation functions (Fig. 1c) and power spectra (Fig. 1d) for specific rotators, one from each population (blue and orange lines for E and I unit, respectively), and use their eigenfrequencies $\omega_i$ to calculate the corresponding theoretical predictions (dashed and dotted black lines, respectively). While these statistics clearly indicate the presence of intrinsic oscillations, a pronounced damping of the correlation function (Fig. 1c) and, equivalently, a considerable width of the spectral peak (Fig. 1d) is apparent. The damping of the correlation function is related to the phase diffusion of the stochastic oscillator, which in turn is entirely driven by the network noise (see e.g. Eq. (7a)). Consistent with our choice of strong synaptic couplings onto the I population and, in consequence, a stronger network noise for the units in that population, this damping is considerably stronger for those units. Of note, this difference in damping cannot be captured by the equivalent one-population theory (dashed and dotted grey lines for E and I unit, respectively), as only a single, “aggregate” network noise is described by this theory, see also Fig. 1a,b.

We finally show the population-averaged correlation func-

tions (Fig. 1e) and power spectra (Fig. 1f) of the rotators. The intrinsic eigenfrequencies of the rotators differ significantly even within populations—according to Eq. (28), for excitatory units the standard deviation $\sigma_E = 1$ is as large as mean frequency $\Omega_0^E = 1$; for inhibitory units, the spread is even more larger ($\sigma_I = 4$, $\Omega_0^I = 3$). For this reason, the intrinsic oscillations are much less apparent in the correlation functions, see Fig. 1e. The power spectra are peaked at the respective mean frequencies of both populations, see Fig. 1f, where the spectral width is now a combined effect of the spread of the intrinsic frequencies and the network noise. For the population-averaged rotator statistics, the single-population theory for the equivalent homogeneous network (dashed and dotted grey lines) is closer to the network simulation and two-population theory than for the single rotators (cf. Fig. 1c,d). Notably the spectral widths of the peaks in the power spectra given by the one-population theory are not identical for both populations (Fig. 1f). We remind the reader that we incorporate an explicit average over the correct intrinsic frequency distribution of each population, see Eq. (37). Thus, the improved correspondence of single-population theory for this specific statistics indicates that these spectral widths are to a significant degree determined by the spread of the effective eigenfrequencies.
<figure><img src="image_5.png" /><figcaption>Figure 2. Network statistics for a balanced E-I network with weak input to the inhibitory population. Left: sketch of the network with excitatory and inhibitory populations and corresponding connection strengths. Right: measured correlation functions (a,c,e) and power spectra (b,d,f) for network noise (a,b), single rotators (c,d) and population-averaged rotators (e,f). Parameters as in Fig. 1 except $J_{IE} = 0.2$.</figcaption></figure>

## 2. Weak input to the inhibitory population

We now discuss the opposite case, where the inputs to the I population are chosen to be smaller than the inputs to the E population, see sketch of Fig. 2, left. Again, the two-population theory (dashed and dotted black lines) nicely captures all measured statistics of the network simulations (blue and orange lines). Since we present the same statistics (correlation functions and power spectra of network simulations, the two-population theory, and the single-population theory) as in the previous case, we will focus in the following on those observations that highlight the difference of the present case from the former.

First of all, the network noise is now clearly stronger for the excitatory population than for the inhibitory population, see Fig. 2a,b, consistent with a reversal of ‘strong’ and ‘weak’ synaptic input compared to the situation in Fig. 1. More interestingly, the network noise predicted by the one-population theory exhibits a strong low-frequency component and its statistics does not appear to be a simple interpolation between the statistics of the network noise of excitatory and inhibitory populations. A discrepancy between the one-population theory and the network simulations is again found in the statistics for single rotators (Fig. 2c,d), where the distinct spectral widths of both populations cannot be captured by the approximation of the true network by an unstructured network. We note that the two-population theory describes very well the observed spectra including the weak side peaks, see inset of Fig. 2d. In line with the weak network noise for the I population, the intrinsic oscillation of the I unit is very weakly damped as indicated by the small spectral width (Fig. 2d) and

the persistent oscillations in the correlation function (Fig. 2c). For the population averages of the single rotator statistics, the one-population theory again fares better than for the statistics of an individual rotator. However, in the present case of weak input to the inhibitory population, the spread of the intrinsic frequencies of the I population $\sigma_I = 0.4$ is considerably smaller than the mean frequency $\Omega_0^I = 3$, and the effect of the network noise is more important than in the previous case. This explains the stronger deviation of the population-averaged correlation function and power spectrum predicted by the one-population theory for the I population.

## 3. Inputs of similar strength for E and I units

As a last example, we consider the case where the synaptic input strengths are the same for all units, i.e. $J_{EE} = J_{IE}$ and $J_{EI} = J_{II}$, see the sketch of Fig. 3, left. As predicted by our two-population theory and confirmed by simulations, the network noise statistics are then identical for E and I populations (Fig. 3a,b).

In light of the identical network noise, one might be tempted to think that for the considered situation the theory for the equivalent ‘unstructured’ network should give the same results as the full theory. It is thus instructive to compare the remaining differences between the full heterogeneous network and the one-population theory. As can be seen from Fig. 3 (grey lines), the one-population theory fails to describe the observed network noise and rotator statistics correctly: (i) In the network noise spectrum (Fig. 3b), we observe increased power at low frequencies and the absence of a spectral peak at a non-vanishing frequency, in marked contrast to
<figure><img src="image_6.png" /><figcaption>Figure 3. Network statistics for a balanced E-I network with equal input to the excitatory and the inhibitory populations. Left: sketch of the network with excitatory and inhibitory populations and corresponding connection strengths. Right: measured correlation functions (a,c,e) and power spectra (b,d,f) for network noise (a,b), single rotators (c,d) and population-averaged rotators (e,f). Parameters as in Fig. 1 except $J_{IE} = 0.5$.</figcaption></figure>

the two-population theory and the network simulations; (ii) in the single-rotator correlation functions and power spectra (Fig. 3c,d), the one-population theory predicts a significantly stronger damping of the correlations, an effect which is also reflected in the population-averaged correlation functions and spectra (Fig. 3e,f).

As already briefly mentioned (cf. discussion following Eq. (39)), the remaining discrepancy between the ‘structured’ and the ‘unstructured’ network is due to the difference in connection weights $J_{\alpha E}$ and $J_{\alpha I}$ of excitatory and inhibitory units, respectively, combined with the difference in intrinsic frequencies. The remaining structure of the considered network is that the faster oscillators (I units) still have much stronger synaptic weight than the slower oscillators (E units)—an effect that cannot be taken into account in the one-population theory.

We further remark that in the case of identical couplings for inhibitory and excitatory units, the contribution of recurrent inputs to the spread of the intrinsic frequencies is the same for both populations (Eq. (28)). For our parameters, the total frequency spread is thus identical ($\sigma_E = \sigma_I = 1$) as can also be observed in the single-rotator statistics (Fig. 3c-f). We finally note that while in the present case the one-population theory overestimate the phase diffusion of rotators, this is not a general feature of our setup when coupling strengths $J_{EE} = J_{IE}$ are identical but other parameters can differ. If we, for instance, switch the roles of E and I units with respect to their intrinsic frequencies, see Fig. 4, the predictions for correlation functions and power spectra of the two-population theory become exactly the ones of the one-population theory of Fig. 3 and vice versa. In particular, the one-population theory

now *underestimates* the broadening of the single-rotator power spectra (Fig. 4d,f).

## IV. SUMMARY & OUTLOOK

In this paper we have developed a general theory for the correlation functions of the self-generated network noise and single units in a heterogeneous network of rotators. The heterogeneity was introduced by means of populations with distinct properties of the single elements and distinct connectivities. Our theory, which builds on the work by van Meegen and Lindner, culminated in the system of equations (17), meaning that in the end we have reduced the problem of calculating the correlations from simulating $N = \sum_{\alpha=1}^{P} N_{\alpha}$ elements to solving $P$ ordinary differential equations.

Our results demonstrate, first of all, that the theory works. Despite being strictly valid only in the limit of $N_{\alpha} \to \infty$ for all $\alpha$ (due to the assumption of Gaussianity for network fluctuations), it works already for comparatively small numbers of units; in our test we used only $N_I = 200$ inhibitory units. This means that it can be applied for realistic numbers of oscillators and is useful to capture finite-size effects of the network noise even in cases in which a balance condition is not obeyed.

The theory introduced here could be further explored and generalized in various directions. One could try to solve the resulting differential equations analytically (for examples in a homogeneous situation, see [32]), using simplifying assumptions for the coupling statistics and frequency distributions; explicit formulas for the correlation functions and spectra may provide further insights into the effects of network hetero-
<figure><img src="image_7.png" /><figcaption>Figure 4. Network statistics for a balanced E-I network with equal input to the excitatory and the inhibitory populations, with inverted intrinsic frequencies of excitatory and inhibitory units. Left: sketch of the network with excitatory and inhibitory populations and corresponding connection strengths. Right: measured correlation functions (a,c,e) and power spectra (b,d,f) for network noise (a,b), single rotators (c,d) and population-averaged rotators (e,f). Parameters as in Fig. 1 except $J_{IE} = 0.5$, $\Omega_0^E = 3$, and $\Omega_0^I = 1$.</figcaption></figure>

geneity on the generated network noise.

Beyond the case of two distinct populations, one should explore networks with three or more populations. In particular the diversity of inhibitory cells in the cortex (see e.g. [35]) suggests to study networks with one population of excitatory (pyramidal) cells but two or more types of interneurons (see e.g. [36] for a recent computational study of such a network).

Another extension is to make the rotator coupling more similar to the pulse coupling in recurrent networks of integrate-and-fire models. For homogeneous networks such an approach has been suggested already in [32] and it was shown that the self-consistent power spectra of single spike trains could be well approximated if the single integrate-and-

fire neurons of the recurrent network operated in a strongly mean-driven regime. We cannot think of a reason why the same approach should not work for a heterogeneous network.

Furthermore, the inclusion of intrinsic noise sources as in [33] and of external signals is certainly worth the exploration but may require also complicated corrections to the Gaussian theory put forward here. Inclusion of common stimuli would be interesting also from an information-theoretic perspective. Besides the autocorrelation functions of the autonomous network activity, then the correlation between information-carrying stimuli and the network activity may be calculated permitting to quantify the information transmission by the network. All of these are exciting topics for future studies.

[1] Y. Kuramoto, *Chemical Oscillations Waves and Turbulence* (Dover, 2003).

[2] A. Pikovsky, M. Rosenblum, and J. Kurths, *Synchronization: A universal concept in nonlinear sciences* (Cambridge Univ. Press, U. K., 2001).

[3] V. S. Anishchenko, V. V. Astakhov, A. B. Neiman, T. E. Vadivasova, and L. Schimansky-Geier, *Nonlinear Dynamics of Chaotic and Stochastic Systems* (Springer, Berlin, 2002).

[4] J. A. Freund, L. Schimansky-Geier, and P. Hänggi, *Chaos* **13**, 225 (2003).

[5] J. A. Acebrón, C. J. P. V. L. L. Bonilla, F. Ritort, and R. Spigler, *Rev. Mod. Phys.* **77**, 137 (2005).

[6] N. Brunel and D. Hansel, *Neural Comp.* **18**, 1066 (2006).

[7] B. Doiron, J. Rinzel, and A. Reyes, *Phys. Rev. E.* **74**, 030903 (2006).

[8] V. Braitenberg and A. Schüz, *Cortex: Statistics and Geometry of Neuronal Connectivity* (Springer, Berlin London New York, 1998).

[9] O. Sporns, *Networks of the Brain* (MIT Press, Cambridge MA, 2016).

[10] C. van Vreeswijk and H. Sompolinsky, *Science* **274**, 1724 (1996).

[11] A. Renart, J. D. L. Rocha, P. Bartho, L. Hollender, N. Parga, A. Reyes, and K. D. Harris, *Science* **327**, 587 (2010).

[12] M. Helias, T. Tetzlaff, and M. Diesmann, *PLoS Comp. Biol.* **10**, e1003428 (2014).

[13] S. Ostojic, *Nat. Neurosci.* **17**, 594 (2014).

[14] J. F. A. Poulet and C. C. H. Petersen, *Nature* **454**, 881 (2008).

[15] K. D. Harris and A. Thiele, *Nat. Rev. Neurosci.* **12**, 509 (2011).

[16] D. J. Amit and N. Brunel, *Cereb. Cortex* **7**, 237 (1997).
[17] N. Brunel and V. Hakim, Neural Comput. **11**, 1621 (1999).

[18] N. Brunel, J. Comput. Neurosci. **8**, 183 (2000).

[19] A. Lerchner, G. Sterner, J. Hertz, and M. Ahmadi, Network : Computation in Neural Systems **17**, 131 (2006).

[20] B. Dummer, S. Wieland, and B. Lindner, Front. Comp. Neurosci. **8**, 104 (2014).

[21] S. Wieland, D. Bernardi, T. Schwalger, and B. Lindner, Phys. Rev. E **92**, 040901(R) (2015).

[22] R. Pena, S. Vellmer, D. Bernardi, A. Roque, and B. Lindner, Front. Comp. Neurosci. **12** (2018).

[23] B. Lindner, Phys. Rev. E. **73**, 022901 (2006).

[24] C. F. Mari, Phys. Rev. Lett. **85** (2000).

[25] A. Lerchner, C. Ursta, J. Hertz, M. Ahmadi, P. Ruffiot, and S. Enemark, Neural Comput. **18**, 634 (2006).

[26] S. Vellmer and B. Lindner, Phys. Rev. Res. **1**, 023024 (2019).

[27] H. Sompolinsky, A. Crisanti, and H. J. Sommers, Phys. Rev. Lett. **61**, 259 (1988).

[28] J. Schücker, S. Goedeke, D. Dahmen, and M. Helias, arXiv preprint arXiv:1605.06758 (2016).

[29] J. Aljadeff, M. Stern, and T. Sharpee, Phys. Rev. Lett. **114**, 088101 (2015).

[30] J. Kadmon and H. Sompolinsky, Phys. Rev. X **5**, 041030 (2015).

[31] F. Mastrogiuseppe and S. Ostojic, PLoS Comput. Biol. **13**, e1005498 (2017).

[32] A. Van Meegen and B. Lindner, Phys. Rev. Lett. **121**, 258302 (2018).

[33] J. Ranft and B. Lindner, Chaos **32**, 063131 (2022).

[34] C. W. Gardiner, *Handbook of Stochastic Methods* (Springer-Verlag, Berlin, 1985).

[35] X. Jiang, S. Shen, C. R. Cadwell, P. Berens, F. Sinz, A. S. Ecker, S. Patel, and A. S. Tolias, Science **350**, aac9462 (2015).

[36] D. Bernardi, G. Doron, M. Brecht, and B. Lindner, PLoS Comput. Biol. **17**, e1007831 (2021).