The effect of heterogeneity on decorrelation mechanisms in spiking neural networks: a neuromorphic-hardware study

High-level brain function such as memory, classification or reasoning can be realized by means of recurrent networks of simplified model neurons. Analog neuromorphic hardware constitutes a fast and energy efficient substrate for the implementation of such neural computing architectures in technical applications and neuroscientific research. The functional performance of neural networks is often critically dependent on the level of correlations in the neural activity. In finite networks, correlations are typically inevitable due to shared presynaptic input. Recent theoretical studies have shown that inhibitory feedback, abundant in biological neural networks, can actively suppress these shared-input correlations and thereby enable neurons to fire nearly independently. For networks of spiking neurons, the decorrelating effect of inhibitory feedback has so far been explicitly demonstrated only for homogeneous networks of neurons with linear sub-threshold dynamics. Theory, however, suggests that the effect is a general phenomenon, present in any system with sufficient inhibitory feedback, irrespective of the details of the network structure or the neuronal and synaptic properties. Here, we investigate the effect of network heterogeneity on correlations in sparse, random networks of inhibitory neurons with non-linear, conductance-based synapses. Emulations of these networks on the analog neuromorphic hardware system Spikey allow us to test the efficiency of decorrelation by inhibitory feedback in the presence of hardware-specific heterogeneities. The configurability of the hardware substrate enables us to modulate the extent of heterogeneity in a systematic manner. We selectively study the effects of shared input and recurrent connections on correlations in membrane potentials and spike trains. Our results confirm ...


I. INTRODUCTION
Dynamical systems in nature often exhibit a remarkable degree of diversity, specialization or anticorrelation across their components, despite equalizing factors such as common input or homogeneity in component and interaction parameters. In many cases, these observations can be explained by the effect of negative feedback. Cell differentiation caused by lateral inhibition [1], formation of new species driven by competition [2] or antiferromagnetism [3] constitute just a few examples. In recurrent neuronal networks, inhibitory feedback constitutes a powerful decorrelation mechanism which allows different neurons to respond nearly independently, even if they are driven by largely overlapping local or external inputs [4][5][6]. Decorrelation by negative feedback hence implements an efficient form of redundancy reduction. In biological systems, it may serve similar purposes as decorrelation in technical applications, where it is used in data compression (e.g., principal-component analysis [7]), crosstalk reduction (e.g., in digital signal processing [8]), echo suppression (e.g., in acoustics [9]) or random-number generation in hardware [10]. Moreover, inhibitory feedback suppresses "quantization noise" at low frequencies and can thereby increase the dynamical range and signal-to-noise ratio for the encoding of analog signals in the spiking activity of recurrent neural networks [11]. It is tempting to exploit these mechanisms in synthetic, neurally inspired architectures such as analog neuromorphic hardware.
Analog neuromorphic hardware mimics properties of biological neural systems using physical models of neurons and synapses (capacitors, for example, emulate insulating cell membranes) [12,13]. The temporal evolution of the analog circuits represents a solution to the corresponding model equations. In consequence, neuralnetwork emulations on analog neuromorphic hardware are massively parallel, extremely fast and energy efficient. Analog neuromorphic devices are therefore highly attractive as tools for neuroscientific research, e.g., for the investigation of learning on long time scales, and technical applications [14][15][16][17]. A biologically inspired neural network (olfactory system of insects) performing rapid online data (odor) classification, for example, has recently been successfully implemented on the analog neuromorphic hardware system Spikey [18,19]. In this application, decorrelation by inhibition is an essential ingredient to guarantee high classification performance. The suppression of quantization noise by inhibitory feedback [11] has been used as a means of noise shaping in several neuromorphic-hardware applications, aiming at the construction of biologically inspired ultra-low power analogto-digital converters [? ? ].
For the functional performance of neuronal architectures, the level of correlations between the activities of individual neurons is often pivotal. Whether such correlations are beneficial or not is context dependent. A number of previous studies emphasize a functional benefit of certain types of correlation for encoding/decoding of information in/from populations of neurons [20][21][22], information transmission [23][24][25], robustness against noise [26], or gain control of postsynaptic neurons [27]. Other studies argue that positive cross-correlations are detrimental as they decrease the precision or sparseness of population codes [19,[28][29][30][31]. Cohen & Maunsell [32], for example, have shown that decreased spike-train correlations in macaque visual area V4 are accompanied by increased behavioral performance in an orientation changedetection task. Depending on the similarity between the trial-averaged responses of different neurons to external stimuli (signal correlation), noise correlations (correlations not explained by signal correlations) can either increase or decrease the amount of information that can be encoded in or decoded from a population of neurons. In populations of neurons with high signal correlation, vanishing or even negative noise correlations are desirable to improve the population code [21].
In finite neural networks, an inevitable source of correlated neural activity is common presynaptic input, shared by multiple postsynaptic neurons. In network models and in-vivo recordings, however, pairwise correlations in the activity of neighboring neurons have been found to be substantially smaller than expected given the amount of shared input [4, 5, [33][34][35][36]. In several studies, this observation has been explained by inhibitory coupling. While Ly et al. [37] and Middleton et al. [38] primarily focused on the effect of feedforward inhibition, Renart et al.
[4], Wiechert et al. [39], and Tetzlaff et al. [5] attributed the smallness of correlations to an active decorrelation of neural activity by inhibitory feedback. The mechanism underlying this active decorrelation has already been described by Mar et al. [11]. In this study, the authors focused on the suppression of low-frequency fluctuations of the population firing rate by recurrent dynamics. As the amplitude of population-rate fluctuations is directly linked to pair-wise correlations (see, e.g., [40]), the effect described in [11] corresponds to a suppression of pairwise correlations in the spiking activity. The theory underlying decorrelation by inhibitory feedback suggests the ef-fect to be general: Decorrelation should be observable in any system with sufficiently strong inhibitory feedback, irrespective of the details of the network structure and the cell and synapse properties. For networks of spiking neurons, however, the effect has so far been explicitly demonstrated only for the homogeneous case, where all neurons have identical properties, receive (approximately) the same number of inputs, and, hence, fire at about the same rate [4,5]. Moreover, the sub-threshold dynamics of individual neurons was assumed to be linear.
Biological neuronal networks typically exhibit broad, heavy-tailed firing-rate distributions [41][42][43][44][45][46][47], indicating a high degree of heterogeneity, e.g., in synaptic weights [48][49][50][51][52], in-degrees [53] or time constants [46,54]. The same holds for neural networks implemented on analog neuromorphic hardware. All analog circuits suffer from device variations caused by unavoidable variability in the manufacturing process. Neurons and synapses implemented in analog neuromorphic hardware therefore exhibit heterogeneous response properties, similar to their biological counterparts [55,56]. To understand the dynamics and function of recurrent neural networks in both biological and synthetic substrates, it is therefore essential to account for such heterogeneities.
Previous work on recurrent neural networks has shown that heterogeneity in single-neuron properties or connectivity broadens the distribution of firing rates [46,57] and affects the stability of asynchronous or oscillatory states [53,[58][59][60][61][62]. A number of studies pointed at a potential benefit of heterogeneity for the informationprocessing capabilities of neural networks [62][63][64][65][66][67][68][69][70][71][72][73]. The effect of heterogeneity on correlations in the activity of recurrent networks of spiking neurons, however, remains unclear. Padmanabhan & Urban [67] have shown that the responses of a population of unconnected neurons are decorrelated by heterogeneity in the neuronal response properties. These results are supported by the subsequent theoretical analysis in [70]. In the following, we refer to this type of decorrelation by heterogeneity as feedforward decorrelation. It does not account for the effect of the recurrent network dynamics. Active decorrelation due to inhibitory feedback [see above; 4, 5], in contrast, constitutes a very different mechanism. The effect of heterogeneity on this feedback decorrelation has lately been studied by Bernacchia & Wang [74] in the framework of a recurrent network of linear firing-rate neurons. In this setup, correlations are suppressed by heterogeneity in the network connectivity (distributions of coupling strengths or random dilution of connectivity). It remains unclear, however, whether this holds true for networks of (nonlinear) spiking neurons.
In this study, we investigate the impact of heterogeneity on input and output correlations in the asynchronous regime of sparse networks of leaky integrateand-fire (LIF) neurons with conductance-based synapses. Emulation of the networks on the analog neuromorphic hardware system Spikey ( Figure 1) [18,75] enable us to investigate the impact of substrate specific properties on a b FIG. 1. The neuromorphic hardware system Spikey. (a) Photograph of the Spikey chip (size 5 × 5 mm 2 ). It comprises analog circuits of 384 neurons and 98304 synapses, is highly configurable and emulates neural-network dynamics with a speed-up of 10 4 with respect to biological real-time. (b) Photograph of the partly cased Spikey system, carrying the Spikey chip (covered by a black round seal) and conventional memory. The system is connected to the host computer via USB 2.0, consumes 6 W of power in total and less than 1 nJ per synaptic transmission (see Supplements 1). the network dynamics. Insights about the interplay between features of the computing substrate and network dynamics are a necessary prerequisite for the development of algorithms that exploit the benefits of analog neuromorphic systems at best. The configurability of this system [18] enables us to systematically vary the level of heterogeneity, and to disentangle the effects of heterogeneity on feedforward and feedback decorrelation (see above). For simplicity, we focus on purely inhibitory networks, thereby emphasizing that active decorrelation by inhibitory feedback does not rely on a dynamical balance between excitation and inhibition [5,6]. We show that decorrelation by inhibitory feedback is effective even in highly heterogeneous networks with broad distributions of firing rates (Section III A). Increasing the level of heterogeneity has two effects: Feedforward decorrelation is enhanced, feedback decorrelation is impaired. Due to the latter, the overall input and output correlations do not necessarily become smaller with increasing heterogeneity. They can even increase (Section III B).
Note that results from specific network emulations on hardware do not directly translate to those obtained by simulations on conventional computers, because the dynamics, parametrization and interplay of analog circuits is very complex and difficult to reproduce with classical simulations. If simplified models for spatial and temporal variability are considered in software simulations, however, emulation results can be reproduced qualitatively, thereby verifying the design of the hardware system. While our hardware system is designed to physically implement biologically realistic neural algorithms in a fast and energy-efficient way, software simulations are used as a complementary tool to isolate, verify and investigate different hardware features, such as spatial and temporal parameter variations. Due to the limited access and configurability of network parameters, this would be difficult to achieve with hardware studies alone. In analogy to the necessity of performing experiments on biological neural systems to verify assumptions made in Computational Neuroscience, actual emulations on neuromorphic hardware are essential to understand its properties and develop efficient neural algorithms for these devices. The fact that our main findings hold true for both emulations on hardware and simulations with software, and that they can be distilled to simple linear models, support their broad relevance and robustness.

A. Network model
Details on the network, neuron and synapse model are provided in Table I. Parameter values are given in Table II. Briefly: We consider a purely inhibitory, sparse network of N (N = 192, unless stated otherwise) LIF neurons with conductance-based synapses. Each neuron receives input from a fixed number K = 15 of randomly chosen presynaptic sources, independently of the network size N . Self-connections and multiple connections between neurons are excluded. Resting potentials E l are set above the firing thresholds Θ (equivalent to applying a constant supra-threshold input current). We thereby ensure autonomous firing in the absence of any further external input. Due to temporal noise, the initial conditions are essentially random.
B. Network emulations on the neuromorphic-hardware system Spikey The Spikey chip (Figure 1) consists of physical models of LIF neurons and conductance-based synapses with exponentially decaying dynamics (for details, see Table I). The emergent dynamics of these physical models represents a solution for the model equations of neurons and synapses in continuous time, in parallel for all units. In contrast, in classical simulations on von-Neumann architectures, model equations are solved by step-wise numerical integration, where parallelization is limited by the available number of processor cores. To emphasize the difference between simulations using software and simulations using physical models, the term emulation is used for the latter [18].
The response properties of physical neurons and synapses vary across the chip due to unavoidable variations in the production process that manifest in a spatially disordered pattern (fixed-pattern noise). In contrast to the approximately static fixed-pattern noise, temporal noise, including electronic noise and transient experiment conditions (e.g., chip temperature), impairs the reproducibility of emulations. In general, two network emulations with identical configuration and stimulation do not result in identical network activity. Both fixed-pattern and temporal noise need to be taken into account when developing models for analog neuromorphic hardware.
The key features of the Spikey chip are the high acceleration and configurability of the analog network implementation. Some network parameters, e.g., synaptic weights and leak conductances, are configurable for each unit, while other parameters are shared for several units (for details see [18]). The hardware system is optimized for spike in-and output and allows to record the membrane potential of one (arbitrarily chosen) neuron with a sampling frequency of 96 MHz in hardware time. On the Spikey chip, capacitances are smaller and conductances are much higher than in biological nervous systems. In consequence, networks on the Spikey chip are emulated with a speed-up of approximately 10 4 with respect to biological real-time. Due to this high acceleration of the neuromorphic chip, the data bandwidth of the connection between the neuromorphic system and the host computer is not sufficient to communicate with the chip in real time. Consequently, input and output spikes (for stimulation and from recordings, respectively) are buffered in a local memory next to the chip. The high acceleration of the Spikey chip allows most of the transistors to operate outside of weak inversion, thereby reducing the effect of transistor variations and minimizing fixed-pattern noise.
Access to the Spikey system is encapsulated by the simulator-independent language PyNN [85,86], providing a stable and user-friendly interface. PyNN integrates the hardware into the computational neuroscience tool chain and has facilitated the implementation of several network models on the Spikey chip [18,19,[87][88][89]. On the Spikey system, a spiking neural network is emulated as follows ( Figure 2a): First, the network described in PyNN is mapped to the Spikey chip, i.e., neurons and synapses are allocated and parametrized. Second, input spikes, if available, are prepared on the host computer and transferred to the local memory on the hardware system. Third, the emulation is triggered and available input spikes are generated. Output spikes and membrane data are recorded to local memory. Last, spike and membrane data are transferred to the host computer and scaled back into the biological domain of the PyNN model description.
For consistency with the model description and simplified comparison to the existing literature, all hardware times and all hardware voltages are expressed in terms of the quantities they represent in the neurobiological model, throughout this study.

C. Experimental setup
To differentiate and compare the effects of shared inputs and feedback connections on correlations, we investigate two different emulation scenarios: First, we emulate networks with intact feedback (FB, Figure 2b), and second, the contribution of shared input is isolated by randomizing the temporal order of this feedback (RAND, Figure 2d).
In the RAND scenario, the inputs of neurons are decoupled from their outputs. Spatio-temporal correlations in presynaptic spike trains are removed by randomizing the presynaptic spike times.
Input correlations between neurons are measured via their free membrane potential, i.e., the membrane potential with disabled spiking mechanism (technically, the threshold is set very high). Because membrane potential traces can be recorded in the hardware only one at a time, traces are obtained consecutively, while repeatedly replaying the previously recorded activity of the FB network to a population of unconnected neurons of equal size. We keep the connectivity the same, and hence each neuron receives the same number of spikes as in the recurrent network during the whole emulation, either without (FB replay , Figure 2c) or with randomization of presynaptic spike times (RAND, Figure 2d), respectively. To preserve the fixed pattern of variability of synaptic weights in hardware, the same hardware synapses are used for each connection in both scenarios. If network dynamics were reproduced perfectly, membrane potential traces and spike times would be identical in the FB and FB replay cases (see also Section II D).
Drawing two different network realizations (i.e., the connectivity matrix) results in the allocation of different hardware synapses, and, due to fixed-pattern noise, in different values of synaptic weights. To average over this variability, throughout this study, emulation results are averaged over M = 100 network realizations, if not stated otherwise.

D. Reproducibility of hardware emulations
Since the initial conditions of the recurrent network on hardware are undefined, consecutive emulations of the FB network result in different network activities. In the RAND and FB replay case, however, the input of neurons is decoupled from their output. Although unavoidable temporal noise is present, the system's state-space trajectory returns to the trajectory of the previously recorded FB case. A certain degree of reproducibility is required for two reasons: First, the investigated effect of decorrelation by inhibitory feedback requires a precise relation between spike input and output. Thus our method of replacing the feedback loop by replay is only valid if temporal noise does not substantially corrupt this relationship. Second, to record the membrane potentials of all neurons, as if recorded at once, neuron dynamics have to be reasonably similar in consecutive emulations.
We measure the reproducibility of neuron dynamics by comparing consecutive emulations with identical configuration, i.e., connectivity and stimulation. For this purpose the spiking activity of a FB network is first recorded ( Figure 2b) and then repeatedly replayed ( Figure 2c). Reproducibility is quantified by the correlations (κ X in Table III) of free membrane potential traces and output spike trains obtained for individual neurons in L = 25 different trials.
Free membrane potentials are reproduced quite well, while spike trains show larger deviations across trials ( Figure 3). Small deviations in the membrane potential ( Figure 3b) are amplified by the thresholding procedure [39,90,91] and can lead to large differences between spike trains (Figure 3c). Consequently, measures based on data of several consecutive replays are more precise for membrane potentials than for spike trains. Nevertheless, results have to be interpreted with care in both cases.

E. Calibration
The heterogeneity of the Spikey hardware is adjusted by calibrating the leak conductance 1 for each individual neuron, compensating for fixed-pattern noise of neuron parameters. To this end, a population of unconnected neurons is driven by a constant supra-threshold current and the time-averaged population activityr is measured. Then, we applied the bisection method [92] to adjust the leak conductance g l of each neuron, such that the neuron's firing rate matches the target rater. This results in calibration values b for the leak conductance g l = g l,0 (1 + b), where g l,0 is the leak conductance before calibration. Because emulations on hardware are not perfectly reproducible, more precise calibration was achieved by evaluating the median over 25 identically configured trials instead of single trials. Furthermore, the bisection method was modified for noisy systems (for details, see Supplements 2).
Intermediate calibration states are obtained by linearly scaling the full calibration: (1) The heterogeneity a is chosen in [0, 1] for calibrations between the uncalibrated (a = 1) and calibrated state (a = 0). In the following, the fully calibrated chip (a = 0) is used, if not stated otherwise. This calibration substantially narrows the distribution of firing rates compared to the uncalibrated state (Figure 4). With respect to the stationary firing rate, variability on the neuron level is reduced from 35.1 to 0.9 s −1 .
Even in the fully calibrated state, leak conductances can still be widely distributed. Due to the chosen calibration procedure, they are likely to be correlated to other parameters that influence the neurons' response to a constant supra-threshold current after calibration. This mutual compensation can lead to similar phenomenology (here: firing rates) despite disparate parameter values, similar to what is observed in biology [93]. In addition to remaining variations in neuron parameters, synaptic parameters are significantly distributed [19,94].

F. Correlation measures
In the following, we introduce definitions used to analyze the recorded data. For clarity, all relevant equations and their parametrization are listed in Table III and IV, respectively.
We quantify correlations of membrane potentials v i (t) and spike trains s i (t) by the population-averaged low- frequency coherence κ V and κ S , respectively. At frequency zero, the coherence corresponds to the normalized integral of the cross-covariance function, i.e., it measures correlations on all time scales. We define the lowfrequency coherence κ X , with X ∈ {S, V }, to be the average coherence over a frequency interval from 0.1 to 20 Hz. In this interval, the suppression of population-rate fluctuations in recurrent networks due to inhibitory feedback is most pronounced, and the coherence is approximately constant. Before calculating the coherence, we convolve the power-and cross-spectra with a rectangular window to average out random fluctuations. This measure, or a variant of it, is commonly used in the neuroscientific literature [4,5,70,91,[95][96][97]. We use the terms lowfrequency coherence and correlation interchangeably.
Throughout this study, the term input correlations is used for correlations between free membrane poten- tials, and output correlations for correlations between spike trains. Shared-input correlations are membranepotential correlations that are exclusively caused by overlapping presynaptic sources, ignoring possible correlations in the presynaptic activity. In homogeneous networks, the average pairwise shared-input correlation is given by the connectivity K/N [5]. In heterogeneous networks, shared-input correlations can be reduced. In the presence of heterogeneous synaptic weights, for example, the shared-input correlation is decreased by a factor 1/(1 + CV 2 J ), where CV J denotes the coefficient of variation of the (non-zero) synaptic weights. Note, however, that heterogeneities which affect only the spike generation but not the integration of synaptic inputs, e.g. distributions of firing thresholds, have no effect on the shared-input correlation.
We assess the significance of correlations by comparing the results from emulations to correlations in surrogate data, in which we removed spatial correlations. For every neuron, we randomly shuffled bins of the membrane potential trace, and assigned a new timestamp uniformly drawn from the emulation interval to every spike, respectively. We thereby remove all spatio-temporal correlations between neurons recorded in parallel. By this procedure we create 100 surrogate trials, across which we calculate the average correlations and the standard error.
To quantify fluctuations in the population activitys ( ure 5e), which we scale with the duration T of the emulation. Consequently, the population power spectrum A(f ), scaled by the population size, coincides with the time-averaged population activityr for high frequencies: [35]. As a measure of pairwise correlations in the time domain ( Figure 5d), we compute the population-averaged cross-correlation function c(τ ) by Fourier transforming the population-averaged cross-spectrum C(f ) to time do-main.

III. RESULTS
In this study, we investigate the roles of shared input, feedback and heterogeneity on input and output correlations in random, sparse networks of inhibitory LIF neurons with conductance-based synapses (Table I), im-plemented on the analog neuromorphic hardware chip Spikey (Figure 1). Similarly to [5], we separate the contributions of shared input and feedback by studying different network scenarios ( Figure 2): In the FB case, we emulate the recurrent network with intact feedback loop ( Figure 2b) and record its spiking activity (Figure 5a). In the FB replay case (Figure 2c), the feedback loop is cut and replaced by the activity recorded in the FB network. Ideally, the input to each neuron in the FB replay case should be identical to the input of the corresponding neuron in the FB network. As the replay of spikes and the resulting postsynaptic currents and membrane potentials are not perfectly reproducible on the Spikey chip, the neural responses in the FB and in the FB replay scenario are slightly different (compare Figures 5a and  5b). In the RAND case (Figures 2d and 5c), we use the same setup as in the FB replay case. However, the spike times in each presynaptic spike train are randomized. While the average presynaptic firing rates and the shared-input structure are exactly preserved in this scenario, the spatio-temporal correlations in the presynaptic spiking activity are destroyed.
Using this setup, we first demonstrate in Section III A that active decorrelation by inhibitory feedback [4,5] is effective in heterogeneous networks with conductancebase synapses over a range of different network sizes. In Section III B, we show that decreasing the level of heterogeneity by calibration of hardware neurons leads to an enhancement of this active decorrelation.

A. Decorrelation by inhibitory feedback
The time-averaged population activities in the FB, FB replay and RAND scenarios are roughly identical (vertical histograms in Figures 5a-c; see also high-frequency power in Figure 5e). In the FB and FB replay scenario, fluctuations in the population-averaged activity are small (horizontal histograms in Figure 5a and b). The removal of spatial and temporal correlations in the presynaptic spike trains in the RAND case leads to a significant increase in the fluctuations of the population-averaged response activity (horizontal histogram in Figure 5c). At low frequencies (≤ 20 Hz), the population-rate power in the FB and in the RAND case differs by about two orders of magnitudes (dark and light gray curves in Figure 5e). This increase in low-frequency fluctuations in the RAND case is mainly caused by an increase in pairwise correlations in the spiking activity ( Figure 5d; the power spectra of individual spike trains [inset in Figure 5e] are only marginally affected by a randomization of presynaptic spike times) [5]. In other words, shared-input correlations, i.e., those leading to large spike-train correlations in the RAND scenario, are efficiently suppressed by the feedback loop in the FB case.
On the neuromorphic hardware, the replay of network activity is not perfectly reproducible (Section II D). While the across-trial variability in membrane potentials is small, postsynaptic spikes are dithered by few milliseconds ( Figure 3). In the FB replay case, the suppression of shared-input correlations by correlations in presynaptic spike trains is slightly less efficient as compared to the intact network (FB). The differences in the population-rate power spectra and in the spike-train correlations between the FB replay and RAND case, respectively, are nevertheless substantial (solid black and light gray curves in Figure 5d and e; note the logarithmic scale; for a detailed investigation of spike dither see Supplements 4 and Supplements Figure 7).
Note that the suppression of correlations and, hence, population-rate fluctuations by inhibitory feedback is restricted to low frequencies (here, to frequencies < 50 Hz; see Figure 5e). In the remainder of this study, we will quantify pairwise correlations by the low-frequency coherence in the range 0.1-20 Hz (see Section II F). At higher frequencies, the population-rate power spectra in the FB, FB replay and RAND case are similar. In Figure 5e, the peaks at ∼ 50 Hz and higher harmonics result from the single-cell spike-train statistics (they are also visible in the single-cell spectra; see inset): A large fraction of cells, in particular those firing at higher rates, generate regular spike trains with low inter-spike-interval (ISI) variability (cf. Figure 6d-f). These (fast spiking) cells contribute maxima to the spike-train spectra at frequencies close to their firing rates (and higher harmonics). The structure of the population-rate spectra at higher frequencies (≥ 50 Hz) is reproduced using surrogate data where the ISI distributions of the individual neurons (and, hence, their firing rates and ISI variability) are preserved, but serial ISI correlations and crosscorrelations between spike trains are destroyed (data not shown).
In the RAND case, presynaptic spike-train correlations were removed, and hence input (i.e., free-membranepotential) correlations are exclusively determined by the number of shared presynaptic sources (Equation 2). If the in-degree K is fixed, input correlations will decrease with network size N (Equation 2, light gray curve and symbols in Figure 7a). In purely inhibitory networks with intact feedback loop (FB scenario), correlations in presynaptic spike trains are on average significantly smaller than zero (dark gray diamonds in Figure 7b, [5]), and largely cancel the positive contribution from shared-input correlations. Average input correlations are therefore significantly reduced (black symbols in Figure 7a). As both shared-input and spike-train correlations scale with the inverse of the network size (N −1 ; light gray curve in Figure 7a and inset in Figure 7b, respectively) [91], this suppression of correlations in the FB (and FB replay ) case is observed for all investigated network sizes N . Note that output correlations are negative even though input correlations are positive. This effect is predicted by theory and also observed in linear network models as well as LIF-network simulations on conventional computers (see Figure 9, Supplements 4, Supplements 5 and Section IV). ). The inset in (b) shows a magnified view of the spike-train correlations in the FB case (dark gray diamonds) with a power-law fit ∼ N −1 (dark gray curve). The light gray horizontal band represents mean ± three standard deviations of (spurious) correlations in surrogate data where correlations were removed. Note that free membrane potentials cannot be recorded in the FB case (see Section II). Hence, there are no gray diamonds in (a).

B. Effect of heterogeneity on decorrelation
In neural networks implemented in analog neuromorphic hardware, neuron and synapse parameters vary significantly across the population of cells (fixed-pattern noise; see Section II B). For a population of mutually unconnected neurons with distributed parameters, injection of a constant (supra-threshold) input current leads to a distribution of response firing rates (Figure 4). In this study, we consider the width of this firing-rate distribution as a representation of neuron heterogeneity. It is systematically varied by calibration of leak conductances. The extent of heterogeneity is quantified by the calibration parameter a (a = 1 and a = 0 correspond to the uncalibrated and the fully calibrated system, respectively; for details, see Section II E). For an unconnected population of neurons subject to constant input, the width of the firing-rate distribution increases monotonically with a.
As shown in Figure 6, the level of heterogeneity (i.e., the calibration state a) is clearly reflected in the activity of the recurrent network (FB case). Both the width of the distribution of mean free membrane potentials (Figure 6a-c) as well as the width of the firing-rate distribution increases with a (Figure 6d-f; horizontal histograms). In the uncalibrated system (a = 1), a substantial fraction of neurons is predominantly driven by constant supra-threshold input currents and therefore generates highly regular spike trains (CV ISI ≈ 0) with high firing rates (r > 60 s −1 ). Simultaneously, about 37% of the neurons are silent (r = 0 s −1 ). Neurons with intermediate firing rates (0 s −1 < r < 15 s −1 ), however, show quite irregular activity (CV ISI > 0.5). After calibration, the firing-rate distribution is narrowed. For a = 0, the fraction of silent neurons is reduced to about 8%. Maximum rates are limited to < 66 s −1 . Note that our calibration routine compensates only for the distribution of neuron parameters, but not for the heterogeneity in synapse properties (synaptic weights, synaptic time constants; see Section IV). Even for the fully calibrated network (a = 0), the firing-rate distribution is therefore still broad. In the RAND case, we obtain similar firing-rate and ISI statistics as in the FB case (Supplements 6).
For all levels of heterogeneity attainable by our calibration procedure (a ∈ [0, 1]), input and output correlations are significantly suppressed by the recurrent-network dynamics (cf. black and dark gray vs. light gray symbols in Figure 8). In a homogeneous, random (Erdős-Rényi) network with fixed in-degree K and linear sub-threshold dynamics, the contribution of shared input to the input (free-membrane-potential) correlation is given by the network connectivity K/N [Equation 2; reference 5] (thin light gray curves in Figure 7a and Figure 8a). Nonlinearities in synaptic and/or spike-generation dynamics [91] as well as heterogeneity in neuron (and synapse) parameters lead to a suppression of this contribution (Equation 3, [67]). Here, we refer to this type of decorrelation as feedforward decorrelation. In fact, in our setup the input and spike-train correlations in the RAND case decrease with increasing heterogeneity (light gray symbols in Figure 8). Even in the fully calibrated case input correlations are slightly smaller than K/N (gray symbols vs. thin gray curve in Figure 8a) due to remaining heterogeneities. Spike-train correlations decrease slightly faster with increasing heterogeneity than input correlations. These observations indicate that both the synaptic integration and spike generation are affected by heterogeneities on hardware. To illustrate the different effects of heterogeneity in synaptic integration and spike generation, we performed network simulations on conventional computers where we distributed either firing thresholds (Figure 9) or synaptic weights (Supplements Figure 4). In the RAND case, the decrease of input correlations on the hardware can be attributed to an increase in heterogeneity in parameters affecting synaptic integration (compare light gray symbols in Figure 8a and Supplements Figure 4a). Heterogeneity in spike thresholds, in contrast, does not affect input correlations in the RAND scenario, but strongly reduces spike-train correlations (light gray symbols in Figure 9). Overall, feedforward decorrelation, i.e. the suppression of correlations in the RAND case, becomes more effective in networks with heterogeneous cell parameters.
Despite this enhancement of feedforward decorrelation, input and output correlations increase with the level of heterogeneity in the presence of an intact feedback loop (FB and FB replay scenarios; black and dark gray symbols in Figure 8). We attribute this effect to a weakening of the effective feedback loop in the recurrent circuit: In heterogeneous networks with broad firing-rate distributions, neurons firing with low or high rates, corresponding to mean inputs far below or far above firing threshold (see Figure 6a-c), are less sensitive to input fluctuations than moderately active neurons (see Supplements Figure 2). Hence, they contribute less to the overall feedback. In consequence, feedback decorrelation is impaired by heterogeneity (see also Section IV).

IV. DISCUSSION
We have shown that inhibitory feedback effectively suppresses correlations in heterogeneous recurrent neural networks of leaky integrate-and-fire (LIF) neurons with nonlinear subthreshold dynamics, emulated on analog neuromorphic hardware (Spikey; [18,75]). Both input and output correlations are substantially smaller in networks with intact feedback loop (FB), as compared to the case where the feedback is replaced by randomized input while preserving the connectivity structure and presynaptic firing rates (RAND). Our results hence show that active decorrelation of network activity by inhibitory feedback [4, 5] is a general phenomenon which can be observed in realistic, highly heterogeneous networks with nonlinear interaction and sufficiently strong negative feedback. Moreover, the study serves as a proofof-principle that network activity can be efficiently decor- related even on heterogeneous hardware, which can be exploited in functional applications, e.g., in the neuromorphic algorithms developed by Pfeil et al. [18] and Schmuker et al. [19].
Functional neural architectures often rely on stochastic dynamics of its constituents or on some form of background noise [see, e.g., 18,19,98]. Deterministic recurrent neural networks with inhibitory feedback could provide decorrelated noise to such functional networks, both in artificial as well as in biological substrates. In neuromorphic hardware applications, these "noise networks" could thereby replace conventional random-number generators and avoid a costly transmission of background noise from a host computer to the hardware substrate (which may be particularly relevant for mobile applications with low power consumption; see Supplements 1). It needs to be investigated, however, how well functional stochastic circuits perform in the presence of such network-generated noise.
Partial calibration of hardware neurons allowed us to modulate the level of network heterogeneity and, therefore, to systematically study its effect on correlations in the network activity. The analysis revealed two counteracting contributions: As shown in previous studies [e.g., 67], neuron heterogeneity decorrelates (shared) feedforward input (feedforward decorrelation). On the other hand, however, heterogeneity impairs feedback decorrelation (see next paragraph). In our network model, this weakening of feedback decorrelation is the dominating factor. Overall, we observed a slight increase in correlations with increasing level of heterogeneity. We cannot exclude that feedforward decorrelation may play a more significant role for different network configurations (e.g., different connection strengths or network topologies, different structure of external inputs, different types of heterogeneity). Our study demonstrates, however, that heterogeneity is not necessarily suppressing correlations in recurrent systems. In this context, it would be interesting to investigate the interplay of signal and noise correlations in the presence of network heterogeneities in recurrent systems. We leave this intriguing topic to future studies.
As shown in [5], feedback decorrelation in recurrent networks becomes more (less) efficient with increasing (decreasing) strength of the effective negative feedback. For networks of spiking neurons, the effective connection strength w ij between two neurons j and i corresponds to the total number of extra spikes emitted by neuron i in response to an additional input spike generated by neuron j (see, e.g., [99]). Assuming that the effect of a single additional input spike is small, the effective connectivity can be obtained by linear-response theory [91]. Note that the effective weights w ij depend on the working point, i.e., the average firing rates of all pre-and postsynaptic neurons (mathematically, w ij is given by the derivative of the stationary response firing rate r i = φ i (r 1 , . . . , r j , . . . , r N ) of neuron i with respect to the input firing rate r j , evaluated at the working point; for details, see [5]). Neurons firing at very low or very high rates are typically less sensitive to input fluctuations than neurons firing at intermediate rates (due to the shape of the response function φ i (r 1 , . . . , r N )). Their dynamical range is reduced. In consequence, they hardly mediate feedback in a recurrent network. In heterogeneous networks with broad distributions of firing rates, the number of these insensitive neurons is increased. Hence, the effective feedback is weak-ened (see Supplements 3). We can qualitatively reproduce this effect of heterogeneity on correlations in recurrent networks (FB case) by means of a simplified linear rate model where increasing heterogeneity is described as a decrease in the effective-weight amplitudes (see Supplements 5). A more quantitative analysis requires an explicit mapping of the synaptic weights in the LIF-neuron network to the effective weights of the linear model (as in, e.g., [100]) in the presence of distributed firing rates. We commit this task to future studies. Note that the rate dependence of the effective weights and the resulting effects on correlations are consistent with our observation that neuron pairs with very low firing rates exhibit spiketrain correlations close to zero, whereas pairs with high firing rates are positively correlated (see Supplements 7). Pairs with one neuron firing at an intermediate rate often exhibit negative spike-train correlations. As shown in [4, 5], these negative spike-train correlations are essential for compensating the positive contribution of shared inputs to the total input correlation (at least in purely inhibitory networks). Narrowing the firing rate distribution (e.g., by calibration of hardware neurons) increases the number of neurons contributing to the negative feedback, which, in turn, leads to more neuron pairs with negative spike-train correlations and, therefore, to smaller overall correlations.
Seemingly contrary to our findings, Bernacchia & Wang [74] report a decrease in correlations with increasing level of heterogeneity. The results of their study are obtained for a linear network model, which can be considered the outcome of the linearization procedure described above. Hence, the connectivity of their model corresponds to an effective connectivity (see above). Their study neglects the rate (working-point) dependence of the effective weights and can therefore not account for the effect of firing-rate heterogeneity. In [74], heterogeneity is quantified by the variance of the (effective) weight matrix (Equations 2.2 and 2.4 in [74]). For sparse connectivity matrices (with a large number of zero elements), the variance of the weight matrix reflects not only the width of the non-zero-weight distribution, but also its mean (Equation 2.4 in [74]). For networks of nonlinear spiking neurons, heterogeneities in neuron and/or synapse parameters broadens the distribution of non-zero effective weights, but may simultaneously reduce its mean (see above, Supplements 3, and [46,100]). Hence, the variance of the full weight matrix may decrease (for illustration, see Supplements Figure 9). In other words, increasing heterogeneity in the nonlinear system may correspond to decreasing heterogeneity in the linearized system. A direct test of this hypothesis requires an explicit linearization of the nonlinear heterogeneous system.
The results of this study were obtained by network emulations on analog neuromorphic hardware. We reproduced the main findings by means of conventional computer simulations of LIF-neuron networks with distributed firing thresholds (see Figure 9). The focus on threshold heterogeneity allows us to isolate the effect of firing-rate distributions on correlations. It does not affect shared-input correlations (see Section II F). Although networks simulated on conventional computers and those emulated on the neuromorphic hardware differ in several respects (e.g., in the exact implementation of heterogeneity or the synapse model; compare Supplements Table I and II to Table I and II, respectively), the qualitative results are very similar: In networks with intact feedback loop, input and output correlations are substantially reduced (as compared to the case where the feedback is replaced by randomized input), but increase with the extent of heterogeneity. As predicted by the theory for homogeneous inhibitory networks, we observe positive input correlations and negative output correlations (see Equation 21 in [5] and in the paragraph which follows; see also [103] and Supplements 5). Further, note that heterogeneity in neuron parameters does not "average out" in larger networks. Upscaling the network size by a factor of 25 (N = 4800, in-degree K = 375) yields smaller spike-train correlations, but the qualitative results are similar to those obtained for the smaller network (N = 192, K = 15) emulated on the Spikey chip (compare Figure 8 to Supplements Figure 3).
In networks with intact feedback loop (FB and FB replay scenarios), the precise spatio-temporal structure of spike trains arranges such that the self-consistent input and output correlations are suppressed. Perturbations of this structure in the local input typically lead to an increase in correlations [5]. In this study, we demonstrate this by replaying spiking activity after randomization of spike times, i.e., by replacing the time of each input spike by a random number uniformly drawn from the full emulation time interval [0, T ) (RAND case). However, even subtle modifications of input spike trains, such as random dither of spike times by few milliseconds, lead to an increase of correlations. On the neuromorphic hardware, replay of spike trains is not entirely reproducible (see Section II D and Supplements 9). Hence, spike-train correlations measured in the FB replay mode are slightly larger than in the FB case. We would expect the same effect on the input side (free membrane potentials). Due to hardware limitations, however, we can measure input correlations only in replay mode (FB replay or RAND), but not in the fully connected network (FB). Therefore, all reported input correlations are likely to be slightly overestimated. In conventional network simulations, we mimicked the effect of unreliable replay by input-spike dithering and, indeed, find a gradual increase in input and output correlations (see Supplements Figure 6). These results seem to be contrary to the study by Rosenbaum & Josic [104], in which synaptic noise leads to a decrease of output correlations in a feedforward scenario. In our case, spiketrain correlations, which suppress shared-input correlations, are removed by dithering spikes, thereby increasing correlations on the output side. In [104], in contrast, spike-train correlations are always zero, and shared-input correlations are decreased by synaptic failure, explaining the decreased output correlations. We attribute this con- tradiction to the missing feedback loop in their system, and expect correlations to increase in recurrent networks subject to similar perturbations. Despite the imperfect replay of input spikes, the decorrelation effect is clearly visible in hardware emulations, both on the input and on the output side. The reproducibility of emulations on neuromorphic hardware could be improved by stabilizing the environment of the system, e.g., the chip temperature or the support electronics (under development). Analog hardware, however, will never reach the level of reproducibility of digital computers. But note that, similar to analog hardware, biological neurons exhibit a considerable amount of trialto-trial variability, even under controlled in-vitro conditions [90]. So far, the details of how neuronal noise, for example, stochastic synapses (spontaneous postsynaptic events, stochastic spike transmission, synaptic failure [105]), affects correlations in recurrent neural circuits remain unclear.
Although different Spikey chips exhibit different realizations of fixed-pattern noise, they show a comparable extent of heterogeneity and yield results which are qualitatively similar to those presented in this article (Supplements 8).
We have shown that negative feedback in recurrent circuits can efficiently suppress correlations, even in highly heterogeneous systems such as the analog neuromorphic architecture Spikey. Correlations can be further reduced by minimizing the level of network heterogeneity. In this study, we reduced the level of heterogeneity through calibration of neuron parameters in the unconnected case (see Section II E). The calibration could, in principle, be improved by calibrating neuron (and possibly synapse) parameters in the full recurrent network. Such calibration procedures are however time consuming and cumber-some. In biological substrates, homeostasis mechanisms [56,106] keep neurons in a responsive regime and reduce the level of firing-rate heterogeneity in a self-regulating manner. Future neuromorphic devices could mimic this behavior, thereby reducing the necessity of time consuming calibration procedures. Alternatively, the analog circuits could be optimized to reduce fixed-pattern noise. This would likely require the allocation of more chip resources, hence reducing the network size per chip area.
For simplicity, this work focuses on purely inhibitory networks (as in [11]). This demonstrates that decorrelation by inhibitory feedback does not rely on a dynamical balance between excitation and inhibition (note that the external "excitatory" drive is constant in our model) [5,6]. Previous studies have shown that, for the homogeneous case, decorrelation by inhibitory feedback is a general phenomenon, which also occurs in excitatoryinhibitory networks, provided the overall inhibition is sufficiently strong (which is typically the case to ensure stability) [4-6, 74]. For the heterogeneous case, computer simulations of excitatory-inhibitory networks show qualitatively the same results as purely inhibitory networks (compare Figure 8 to Supplements Figure 5), confirming that our results generalize to the case of mixed excitatory-inhibitory coupling.
Similar to our study, Giulioni et al. [107] use a theoryguided approach to implement, verify and investigate network dynamics on analog neuromorphic hardware. In their study, an attractor network is implemented that is inspired by a mean-field model. Due to heterogeneities in synaptic efficacies on the hardware, stability analysis of attractor states requires the authors to measure effective response functions of populations of hardware neurons. To this end, they replace recurrent connections in one population of neurons by external input. This al- The acceleration factor is defined as the ratio between the emulated network time T (in biological time) and the execution time (wall clock time). In the record case, a network realization is generated on the host computer and uploaded to the chip. During the subsequent emulation, spike trains are recorded. In the replay case, spikes are replayed and the membrane potential of one neuron is recorded with full sampling frequency (9.6 kHz). The execution time covers the full data flow from a network description in PyNN to the emulation on the Spikey system and back to the network representation in PyNN. The time-averaged population firing rate isr = (23.7 ± 0.2) s −1 . The vertical dashed line depicts the runtime used in this study. The hardware system has to be initialized once before usage (< 1 s), which is not considered here.
lows them to measure the firing rate of the population as a function of the external input, while the activity of the population is in equilibrium with that of other recurrently connected populations. This study represents another example illustrating that investigations of actual hardware emulations are a prerequisite for successful application of analog neuromorphic hardware. This study demonstrates that the Spikey system has matured to a level that permits its use as a tool for neuroscientific research. For the results presented in this study, we recorded in total 10 11 membrane-potential and spike-train samples, representing more than 100 days of biological time. Due to the 10 4 -fold acceleration of the Spikey chip, this corresponds to less than 15 minutes in the hardware-time domain. Interfacing the hardware system, however, reduces the acceleration to an approximately 50-fold speed-up ( Figure 10). The translation between the network description and its hardware representation claims the majority of execution time, more than the network emulation and the transfer of data to and from the hardware system together. Encoding and decoding spike times on the host computer is particularly expensive. Obviously, the system could be opti-mized by processing the data directly on the hardware, or by choosing a data representation which is closer to the format used on the Spikey chip, but this would impair user-friendliness, and hence, the effectiveness of prototyping. While the Spikey system permits the monitoring of the spiking activity of all neurons simultaneously, access to the membrane potentials is limited to a single (albeit arbitrary) neuron in each emulation run. Monitoring of membrane potentials of a population of n neurons therefore requires n repetitions of the same emulation. Extending the hardware system to enable access to the membrane potentials of at least two neurons simultaneously would allow for a direct observation of input correlations in the intact network (and thereby avoid problems with replay reproducibility; see above) and reduce execution time (the Spikey chip itself permits recording of up to eight neurons in parallel, the support electronics, however, does not). While the Spikey system does not significantly outperform conventional computers in terms of computational power, emulations on this system are much more energy efficient (Supplements 1). A substantial increase of computational power is expected for large systems exploiting the scalability of this technology without slow-down [108].

ACKNOWLEDGMENTS
We would like to thank Matthias Hock, Andreas Hartel, Eric Müller and Christoph Koke for their support in developing and building the Spikey system, as well as Mihai Petrovici and Sebastian Schmitt for fruitful discussions. This research is partially supported by the Helmholtz Association portfolio theme SMHB, the Jülich Aachen Research Alliance (JARA), EU Grant 269921 (BrainScaleS), and EU Grant 604102 (Human Brain Project, HBP). We acknowledge financial support by Deutsche Forschungsgemeinschaft and Ruprecht-Karls-Universität Heidelberg within the funding programme Open Access Publishing. All network simulations carried out with NEST (http://www.nest-simulator.org).

Appendix A: Network description
Details on the network model as well as parameter values are provided in Table I and II, respectively.

Appendix B: Description of data analysis
A summary of the quantities and the data analysis used in this study is provided in Table III. Parameter values for the data analysis are given in Table IV v(t) = v reset This model is emulated by analog circuitry on the Spikey chip [13].

Conductance dynamics
For each presynaptic spike at time t * (t > t * + d): g syn (t) ≈ J exp(− t−t * −d τsyn )Θ(t), with J = w hw g max and Heaviside function Θ(t). This model is emulated by analog circuitry on the Spikey chip [18].

Spiking
If v(t * −) < Θ ∧ v(t * +) ≥ Θ: emit spike with time stamp t * TABLE I. Description of the network model (according to [109]).     The Spikey system consumes approximately 6 W of power, and the chip itself less than 0.6 W. On the chip most power is consumed by digital communication infrastructure, which is not part of the neuromorphic network. In the following, we estimate the power consumption for a single synaptic event using the data set partly shown in Figure 5a. This emulation lasts T = 10 s in biological time and generates approximately 45 · 10 3 spikes. Considering the acceleration of the hardware network (10 4 ) and the synapse count per neuron (K = 15), the system generates 7 · 10 8 synaptic events per second in hardware time. If we consider the total power consumption of the Spikey chip, the upper bound of energy consumed by each synaptic transmission will be approximately 1 nJ. Because these measurements include the communication infrastructure and other support electronics to observe spike times and membrane traces, the real energy consumption for synaptic transmissions is estimated to be approximately ten times smaller. Network simulations on conventional supercomputers a far less energy efficient and consume tens of µJ for each synaptic transmission [110].

Supplements 2: Modification of the bisection method
In each iteration of the bisection method that is used to calibrate the leak conductances of hardware neurons (Section II E), we evaluated the firing rate for each neuron by the median over L = 25 identical trials. However, if this measure is compared between consecutive identical iterations, temporal noise on time scales longer than the duration of one iteration may still lead to variability. In the original bisection method, the interval of possible solutions is halved after each iteration step [92]. To improve the convergence of this method in the context of our calibration we expanded the halved interval by 20% at both ends after each iteration. This prevents the algorithm to get stuck in an interval wrongly chosen by random fluctuations of the firing rates.

Supplements 3: Effective weights
We quantify the effect of a single spike of neuron j on the firing rate of a postsynaptic neuron i by the effective weight w ij of the connection i ← j. Assuming that the activity of neuron i does not affect the activity of neuron j (i.e., the RAND case), we define w ij as the cross-correlation between the spike trains s j (t) and s i (t + τ ), where, due to causality, τ is positive. Then, we average the effective weight over the emulated network time T , and subtract the baseline determined by the average correlation for negative τ : (S1) Here, we chose τ max = 50 ms and τ min = 50 ms. τ max was determined by measuring the average duration in which a spike from neuron j has an influence on neuron i (data not shown). τ min was then chosen symmetrically. The mass of the effective weights density shifts towards less negative effective weights for increasing heterogeneity (Supplements Figure 1a). We obtainw by averaging over all possible connections: For increasing heterogeneity the absolute value of the average effective weightw decreases (Supplements Figure 1b). This can be explained by the dependence of the effective weight on the firing rate of the postsynaptic neuron (Supplements Figure 2). In the regime of small rates (< 10 s −1 ), incoming spikes hardly affect the neuron's firing and the effective weight is small. Similarly, in the case for large rates (> 40 s −1 ). Neurons with intermediate firing rates are sensitive to input and hence have a more negative effective weight. Heterogeneity increases the number of neurons with small and large firing rates, and hence the absolute value of the average effective weight decreases, which in turn weakens the effective negative feedback of the network.

Supplements 4: Simulations with software
We validate our results by comparing them to simulations with software (NEST [101], PyNN [102]). In these simulations we modulated the degree of heterogeneity by distributing the firing thresholds of all neurons according to a normal distribution with mean Θ and standard deviation σ Θ . Details about the network, neuron and synapse models and their parameters can be found in Supplements Table I and II, respectively. The results are qualitatively similar to network emulations on the Spikey chip (compare Figure 8 to 9) and also hold for larger network sizes (Supplements Figure 3). In the FB case, input correlations and output correlations increase with network heterogeneity. In contrast to hardware emulations, input correlations stay approximately constant in the RAND case. Output correlations strongly decrease with the standard deviation σ Θ . Here, heterogeneity only affects the output spike times, and not the integrative properties of the neurons (see also Section III B and [70]).
In addition, we compare the effect of heterogeneity in firing thresholds on correlations to that of heterogeneity in synaptic weights. Details about the network, neuron and synapse models and their parameters can be found in Supplements Table I and III, respectively. If we distribute synaptic weights, shared-input correlations decrease with the width of the weight distribution (Supplements Figure 4). Output correlations are also reduced, however, only proportionally to the reduction of input correlations (compare insets in Supplements Figure 4). While output correlations are overall smaller than input correlations due to the non-linearity in spike generation, we do not observe a boost of this decrease for large heterogeneities. Overall, the dynamics of the recurrent system is more sensitive to heterogeneities in firing thresholds than in synaptic weights (compare scale of abscissas of Figure 9 to Supplements Figure 4).
To investigate the generalization of our results to mixed excitatory-inhibitory networks, correlations were measured ). The insets show correlations normalized to unity at σJ = 0, with the same abscissa as in the main plot. Note that in simulations the FB replay is identical to the FB case, and is hence not shown. Networks simulated with NEST [101] and PyNN [102].
in a network of N = 192 neurons consisting of half excitatory and half inhibitory neurons. Details about the network, neuron and synapse models and their parameters can be found in Supplements Table IV and V, respectively. As described above, we distribute the firing thresholds of all neurons according to a normal distribution. The results are consistent with those obtained from purely inhibitory networks on hardware demonstrating the generality of our findings (compare Figure 8 to Supplements Figure 5). By modulating the level of heterogeneity separately for the excitatory or the inhibitory population, we observe that network dynamics are more sensitive to heterogeneities in the inhibitory than in the excitatory population (data not shown).
We investigate the effect of temporal noise on correlations by dithering spikes in the FB replay case before replaying them to the network (for network, neuron and synapse models see Supplements Table I and II). Each spike time t is replaced by a spike time t randomly drawn from a normal distribution with width ϑ: t ∼ N (t, ϑ). Even for small ϑ, correlations increase on the input as well as on the output side, demonstrating the sensitivity of correlations to perturbations in the feedback loop (see Supplements Figure 6 and [5]). This effect is also reflected in an increase of the power of the population activity at small frequencies (see Supplements Figure 7). These results suggest that on hardware temporal noise is responsible for the increase of correlations and population power in the FB replay case (compare FB replay to FB in Figure 5, 7 and 8). ). Correlations in the FB case correspond to correlations in the FB replay case for ϑ = 0 ms. Networks simulated with NEST [101] and PyNN [102].

Supplements 5: Linear model
We investigate the consistency of our results with a linear rate model that allows us to numerically calculate the average correlations from a given connectivity matrix W. The model is defined as (according to, e.g., [100]) Here, r(t) denotes the rate of the individual neurons and x(t) a Gaussian white noise input that is independent for each neuron. The linear filter kernel h(t) depends on the details of the model, is not relevant in our calculation, and hence is not further specified here. Equation S3 can be transformed to Fourier domain, where the input and output spectral matrices can be expressed by with T (ω) = (1 − H(ω)W) −1 [100]. In the RAND case, the linear equation for the rate of the (unconnected) neurons reads wherer(t) has the same auto-correlations as r(t) but zero cross correlations, i.e., CRR = diag(C RR ), since the randomization of spike times removes all spatio-temporal correlations. According to Tetzlaff et al. [5] spectral matrices in the RAND case are given by We calculate the population-averaged power-and cross-spectra from the full matrices: (S10) Here, X ∈ {R, Q} denotes the FB and RAND case, respectively. The low frequency coherence is the cross-spectra normalized by the power spectra: Note that in the linear model we are actually taking the zero frequency coherence.
As in the spiking model, we consider a sparse network, i.e., we randomly choose for each neuron i ∈ [1, N ] an identical number of presynaptic partners (K = 15). In the linear model we do not consider a distribution of non-zero effective weights. Instead, each realized connection is assigned the same weight value −w. To mimic the effect of calibration we vary the absolute value of the effective weight by scaling the weights of the non-zero connections with (S12) This procedure changes the variance of the weight matrix [74] and henceã is denoted the heterogeneity of the network. More homogeneous (heterogeneous) networks have larger (smaller) effective weights and hence stronger (weaker) feedback. We obtain qualitatively similar results as we observe on the Spikey chip (compare Supplements Figure 8 to Figure 8). Output correlations in the RAND case decrease, while both input and output correlations in the FB case increase with network heterogeneity, i.e., with the variance of the effective weight matrix. Note that, in this simplified model, input correlations in the RAND case are not affected by the level of "heterogeneity" because the non-zero weights are homogeneous, i.e. have zero variance, irrespectively ofã (cf. Equation 3). In Supplements Figure 9 we illustrate, how in a sparse network the variance of the weight matrix can increase, although the distribution of non-zero weights becomes narrower. The standard deviation σ w of a distribution of non-zero weights with mean µ w is (in this example) smaller than the standard deviation σ W of the full effective weight matrix, due to the sparseness of the matrix (Supplements Figure 9a; here, we chose = 0.8). If we, at the same time, increase the mean µ w and decrease the standard deviation σ w of non-zero weights, the standard deviation of the weight matrix σ W can increase significantly (Supplements Figure 9). While the distribution of effective weights is broadened, the mean is decreased, which has a greater impact on the size of correlations. This observation could explain the decrease of correlations with increased calibration of the neuromorphic chip.

Supplements 6: Firing statistics in the RAND case
The firing rate distributions in the RAND case are similar to those in the FB scenario (compare Supplements Figure 10 to Figure 6d-f). As in the FB case they become broader for increasing heterogeneity. Due to a larger variance of the membrane potentials, inactive neurons in the FB case have a higher probability to fire in the RAND case (compare gray bars in horizontal histograms of Supplements Figure 10a-c and Figure 6d-f). Neurons with firing rates larger than the average firing rate are strongly driven by constant current influx, and hence show similar firing rates than in the FB scenario. The irregularity of firing increases for the RAND compared to the FB case. This can also be traced back to an increased variance of the membrane potential for the RAND case (data not shown). neurons (Supplements Figure 11). This reveals a dependence of the pairwise correlation on the firing rate of the respective neurons. If both neurons fire at low rate (here < 5 s −1 ), correlations will be close to zero similar to the results in [91]. For high rates (here > 25 s −1 ) we find mostly positive correlations. However, for neurons firing at intermediate rates, the activity of the pair is often anti-correlated, and can hence effectively suppress (positive) shared-input correlations. After calibration, the amount of neurons firing at intermediate rates and showing negative correlations increases (Supplements Figure 11b, Figure 6). Shared-input correlations are hence suppressed by more neurons leading to smaller input correlations than in the uncalibrated case ( Figure 8).

Supplements 8: Results for different Spikey chips
We have performed the same experiments as described in the main text using two additional Spikey chips of the current version (5). Different chips show different realizations of fixed-pattern noise, and hence calibration was repeated for each chip separately (Figure 4, Supplements Figure 12A). For all chips we find qualitatively the same results as described above: in the FB case, free membrane potentials (and spike trains) are decorrelated by inhibitory feedback and correlations increase with an increasing level of heterogeneity ( Figure 8, Supplements Figure 12B).
This also applies to the results of identical experiments on three previous-generation chips (version 4). These old chips show more pronounced heterogeneity in their parameters. In addition, we observe larger quantitative differences between old chips for the uncalibrated case, which is likely to be caused by different extents of intrinsic fixed-pattern noise. Unfortunately, the old systems show an artifact affecting the population power at around 100 Hz (Supplements Figure 13), consistently across all chips of identical revision. We can exclude that this effect can be traced back to network effects (it occurs across the FB, FB replay and RAND cases) or to single-cell power spectra (data not shown). We hence have to consider it an artificial spurious synchronization. As the artifact lies outside the range of frequencies that are relevant for our measures of correlation it does not affect the main results and experiments on the old chips qualitatively confirm our results in the main manuscript (compare Supplements Figure 14 to Figure 4, 8 and Supplements Figure 12).

Supplements 9: Reproducibility of networks with intact feedback
We measured the standard deviation of free membrane potential and spike train correlations over several trials for a single network realization (Supplements Figure 15). This standard deviation is smaller than the standard deviation we observe over different network realizations (compare to Figure 8 and Supplements Figure 12), which indicates that the latter is mostly caused by the different realizations of the connectivity, not by trial-to-trial variability. Note that the variability between trials of networks with intact feedback is likely to be larger than between replays of network activity as shown in Figure 3, because network dynamics may be chaotic. The data shown in Supplements Figure 15 has to be interpreted with care, because the reproducibility of only a single network realization is considered. For different realizations, the mean and standard deviation will in general be different.

Supplements 10: Extraction of effective PSP amplitudes, time constants and delays
In order to provide an estimation for the order of magnitude of various neuron parameters, we fitted the postsynaptic potentials (PSPs) of current-based LIF neurons with exponential postsynaptic currents (PSCs) to the spike-triggered average (STA) of free membrane potentials obtained from hardware emulations. The differential equations describing the dynamics of membrane potentials V j (t) and synaptic currents I j (t) for this simplified model are given by τ eff mVj (t) = − V j (t) + RI j (t) (S13) τ eff synİj (t) = − I j (t) + τ eff syn k J jk s k (t − d eff ) . (S14) where R is the resistance of the membrane, τ eff m the membrane time constant, τ eff syn the synaptic time constant and J jk the synaptic weight of the connection from neuron k to j. Here the superscript eff indicates that these are not the parameters realized on hardware, but of an abstract model that we fit to hardware measurements. Data were obtained from 100 network realizations of the RAND case (same dataset as for Figure 8). For each neuron we computed the STA of its free membrane potential individually for each spike train of its presynaptic partners. The PSP for the current based model described above is given by where we introduced V eff max = RJ to simplify the fitting procedure. For each pair of neurons i, j we obtain a set of parameters: V eff max , τ eff syn , τ eff m and d eff . The fitted curves (using the function scipy.optimize.curve_fit with default parameters from Python's scipy module v0.17.0) match fairly well to the empirical data to allow an estimation of the effective values of the neuron parameters governing the shape of the postsynaptic potentials (Supplements Figure 16a). However, note that the resulting parameter distributions do not precisely reflect the actual parameter distributions on the chip.
The distributions of neuron and synapse parameters change slightly with calibration. For the calibrated network the absolute values of V eff max , τ eff syn , τ eff m and d eff are larger than for the uncalibrated network (see Table II and Supplements Figure 16b).  Table II.    If v(t * −) < Θ ∧ v(t * +) ≥ Θ: emit spike with time stamp t *

E Measurements
Spike trains recorded from P E s excitatory and P I s inhibitory neurons Membrane potentials recorded from P E v excitatory and P I v inhibitory neurons SUP.