Global organization of neuronal activity only requires unstructured local connectivity

Abstract

Modern electrophysiological recordings simultaneously capture single-unit spiking activities of hundreds of neurons spread across large cortical distances. Yet, this parallel activity is often confined to relatively low-dimensional manifolds. This implies strong coordination also among neurons that are most likely not even connected. Here, we combine in vivo recordings with network models and theory to characterize the nature of mesoscopic coordination patterns in macaque motor cortex and to expose their origin: We find that heterogeneity in local connectivity supports network states with complex long-range cooperation between neurons that arises from multi-synaptic, short-range connections. Our theory explains the experimentally observed spatial organization of covariances in resting state recordings as well as the behaviorally related modulation of covariance patterns during a reach-to-grasp task. The ubiquity of heterogeneity in local cortical circuits suggests that the brain uses the described mechanism to flexibly adapt neuronal coordination to momentary demands.

Keywords: balanced state; beyond-mean-field theory; correlated activity; long-range coordination; motor cortex; neuroscience; reach-to-grasp task; rhesus macaque.

Figure 1.. Salt-and-pepper structure of covariances in motor cortex.

(A) Sketch of 10 × 10 Utah electrode array recording in motor cortex of macaque monkey during rest. (B) Spikes are sorted into putative excitatory (blue triangles) and inhibitory (red circles) single units according to widths of spike waveforms (see Appendix 1 Section 2). Resulting spike trains are binned in 1 s bins to obtain spike counts. (C) Population-resolved distribution of pairwise spike-count Pearson correlation coefficients in session E2 (E-E: excitatory-excitatory, E-I: excitatory-inhibitory, I-I: inhibitory-inhibitory). (D) Pairwise spike-count correlation coefficients with respect to the neuron marked by black triangle in one recording (session E2, see Materials and methods). Grid indicates electrodes of a Utah array, triangles and circles correspond to putative excitatory and inhibitory neurons, respectively. Size as well as color of markers represent correlation. Neurons within the same square were recorded on the same electrode.

Figure 2.. Correlations from direct and indirect connections.

(A) Positive correlation (green neuron i) follows from direct excitatory connection (top) or shared input (middle). (B) Negative correlation (magenta) between two excitatory neurons cannot be explained by direct connections: Neuronal interactions are not only mediated via direct connections ($n=1$; sign uniquely determined by presynaptic neuron type) but also via indirect paths of different length $n\S gt;1$. The latter may have any sign (positive: green; negative: purple) due to intermediate neurons of arbitrary type (triangle: excitatory, circle: inhibitory).

Figure 3.. Spatially organized E-I network model.

(A) Network model: space is divided into cells with four excitatory (triangles) and one inhibitory (circle) neuron each. Distance-dependent connection probabilities (shaded areas) are defined with respect to cell locations. (B) Eigenvalues λ of effective connectivity matrix for network in dynamically balanced critical state. Each dot shows the real part $\mathrm{Re}(\lambda )$ and imaginary part $\mathrm{Im}(\lambda )$ of one complex eigenvalue. The spectral bound (dashed vertical line) denotes the right-most edge of the eigenvalue spectrum. (C) Simulation: covariances of excitatory neurons over distance $x$ between cells (blue dots: individual pairs; cyan: mean; orange: standard deviation; sample of 150 covariances at each of 200 chosen distances). (D) Theory: variance of covariance distribution as a function of distance $x$ for different spectral bounds of the effective connectivity matrix. Inset: effective decay constant of variances diverges as the spectral bound approaches one. (E) For large spectral bounds, the variances of EE, EI, and II covariances decay on a similar length scale. Spectral bound $R=0.95$. Other parameters see Appendix 1—table 3.

Figure 4.. Long-range covariances in macaque motor cortex.

Variance of covariances as a function of distance. (A) Population-specific exponential fits (lines) to variances of covariances (dots) in session E2, with fitted decay constants indicated in the legend (I-I: putative inhibitory neuron pairs, I-E: inhibitory-excitatory, E-E: excitatory pairs). Dots show the empirical estimate of the variance of the covariance distribution for each distance. Size of the dots represents relative count of pairs per distance and was used as weighting factor for the fits to compensate for uncertainty at large distances, where variance estimates are based on fewer samples. Mean squared error 2.918. (B) Population-specific exponential fits (lines) analogous to (A), with slopes constrained to be identical. This procedure yields a single fitted decay constant of 1.029 mm. Mean squared error 2.934. (C) Table listing decay constants fitted as in (B) for all recording sessions and the ratios between mean squared errors of the fits obtained in procedures B and A.

Figure 5.. Changes in effective connectivity modify coordination patterns.

(A) Visualization of effective connectivity: A sparse random network with given structural connectivity (left network sketch) is simulated with two different input levels for each neuron (depicted by insets), resulting in different firing rates (grayscale in right network sketches) and therefore different effective connectivities (thickness of connections) in the two simulations. Parameters can be found in Appendix 1—table 4. (B1) Histogram of input currents across neurons for the two simulations (N1 and N2). (B2) Scatter plot of inputs to subset of 1500 corresponding neurons in the first and the second simulation (Pearson correlation coefficient $\rho =0.90$). (C) Correlation coefficients of rates and of covariances between the two simulations (b, black) and within two epochs of the same simulation (w, gray). (D) Correlation coefficient of rates (gray) and covariances (black) between the two simulations as a function of the spectral bound $R$. (E1) Distribution of rates in the two simulations (excluding silent neurons with $\left|\mathrm{rate}\right|\S lt;{10}^{-3}$). (E2) Scatter plot of rates in the first compared to the second simulation (Pearson correlation coefficient $\rho =0.81$). (F1) Distribution of covariances in the two simulation (excluding silent neurons). (F2) Scatter plot of sample of 5000 covariances in first compared to the second simulation (Pearson correlation coefficient $\rho =0.40$). Here silent neurons are included (accumulation of markers on the axes). Other parameters: number of neurons $N=2000$, connection probability $p=0.1$, spectral bound for panels B, C, E, F is $R=0.8$.

Figure 6.. Behavioral condition reshapes mesoscopic neuronal coordination.

(A) Trial structure of the reach-to-grasp experiment (Brochier et al., 2018). Blue segments above the time axis indicate data pieces at trial start (dark blue: S (S1+ S2)) and during the preparatory period (light blue: P (P1+ P2)). (B) Salt-and-pepper structure of covariance during two different epochs (S1 and P1) of one recording session of monkey N (151 trials, 106 single units, Figure 1 for recording setup). For some neurons, the covariance completely reverses, while in the others it does not change. Inhibitory reference neuron indicated by black circle. (C1) Distributions of firing rates during S1 and P1. (C2) Scatter plot comparing firing rates in S1 and P1 (Pearson correlation coefficient $\rho =0.69$). (D1/D2) Same as panels C1/C2, but for covariances (Pearson correlation coefficient $\rho =0.40$). (E) Correlation coefficient of firing rates across neurons in different epochs of a trial for eight recorded sessions. Correlations between sub-periods of the same epoch (S1-S2, P1-P2; within-epoch, gray) and between sub-periods of different epochs (Sx-Py; between-epochs, black). Data shown in panels B-D is from session 8. Box plots to the right of the black dashed line show distributions obtained after pooling across all analyzed recording sessions per monkey. The line in the center of each box represents the median, box’s area represents the interquartile range, and the whiskers indicate minimum and maximum of the distribution (outliers excluded). Those distributions differ significantly (Student t-test, two-sided, $p\ll 0.001$). (F) Correlation coefficient of covariances, analogous to panel E. The distributions of values pooled across sessions also differ significantly (Student t-test, two-sided, $p\ll 0.001$). For details of the statistical tests, see Materials and methods. Details on number of trials and units in each recording session are provided in Appendix 1—table 1.

Appendix 1—figure 1.. Correlations and covariances.

The shown data is taken from session E2. (E-E: excitatory-excitatory, E-I: excitatory-inhibitory, I-I: inhibitory-inhibitory). (A) Population-resolved distribution of pairwise spike-count Pearson correlation coefficients. Same data as in Figure 1C. (B) Population-resolved distribution of pairwise spike-count covariances. (C) Population-resolved distribution of variances. (D) Pairwise spike-count correlation coefficients with respect to the neuron marked by black triangle. Grid indicates electrodes of a Utah array, triangles and circles correspond to putative excitatory and inhibitory neurons, respectively. Size as well as color of markers represent correlation. Neurons within the same square were recorded on the same electrode. Same data as in Figure 1D. (E) Pairwise spike-count covariances with respect to the neuron marked by black triangle.

Appendix 1—figure 2.. Distance-resolved variance of covariance: robustness of decay constant estimation.

Exponential fits (lines) to variances of covariances (dots) analogous to Figure 4A and B in the main text (columns 1&3 and 2&4, respectively) for all analyzed resting state sessions. The two sets of plots differ in E/I separation consistency values chosen during data preprocessing. Panels A-H: default (lowest) required consistency (~0.6), used throughout the main analysis; panels I-P: $c=0.75$ The values of the obtained decay constants are listed in Appendix 1—table 1.

Appendix 1—figure 3.. Rate and covariance stationarity during a reach-to-grasp trial: monkey E.

Black line indicates population mean and gray area +/- 1 population standard deviation of single unit firing rate (left column) and pairwise zero time-lag covariance (right column) during trial of a given session (row). Blue bars indicate starting (S) and preparatory (P) periods used in the analysis (Figure 6 in the main text). First, second and fourth dashed lines indicate visual signals lighting up and the third dashed line indicates the removal of a visual cue and beginning of a waiting period.

Appendix 1—figure 4.. Rate and covariance stationarity during a reach-to-grasp trial: monkey N.

Analogous to Appendix 1—figure 3.

Appendix 1—figure 5.. Stationarity of single-unit activity during a reach-to-grasp trial (monkey N, session i140627-001).

Black lines indicate mean and gray areas +/- 1 standard deviation across trials of single unit activity in each panel (sliding window analysis with 5 ms step size and 100 ms window length). Blue bars indicate starting (S) and preparatory (P) periods used in the analysis (Figure 6 in the main text). First, second and fourth dashed lines (marked with stars) indicate visual signals lighting up and the third dashed line indicates the removal of a visual cue and beginning of a waiting period.

Appendix 1—figure 6.. Stationarity of single-unit activity during a reach-to-grasp trial (monkey N, session i140704-001).

Analogous to Appendix 1—figure 5.

Appendix 1—figure 7.. Comparison of simulation and theory.

(A) Variance of EE, EI, and II covariances as a function of distance. Darker dots are the results of the simulation. Lighter ones are the prediction of the discrete theory. (B) Variance of EE covariances as a function of distance (Equation (56) for variances). The lightest blue dots are the predictions of the discrete theory (${\mu }_a$ replaced by the discrete Fourier transform of Equation (52), taking into account Section 7), the medium blue line is the (0,2)-Padé prediction (${\mu }_a$ replaced by its Padé approximation Equation (53), taking into account Section 7), and the dark blue line is the higher order (1,2)-DLog-Padé prediction (${\mu }_a$ replaced by its Padé approximation Equation (53), using Equation (57), and taking into account Section 7). (C) Fitted slope of linear regions in panel B for different spectral bounds $R$ (light blue: discrete theory, medium blue: Padé approximation, dark blue: higher order Padé approximation).