Uncovering statistical structure in large-scale neural activity with Restricted Boltzmann Machines

TL;DR

Restricted Boltzmann Machines are used to model the activity of simultaneously recorded neurons from the Allen Institute Visual Behavior Neuropixels dataset, spanning multiple cortical and subcortical regions of the mouse brain, and reproduce the global relaxation dynamics of neural activity.

Abstract

Large-scale electrophysiological recordings now allow simultaneous monitoring of thousands of neurons across multiple brain regions, revealing structured variability in neural population activity. Understanding how these collective patterns emerge from microscopic neural interactions requires models that are scalable, predictive, and interpretable. Statistical physics provides principled frameworks to address this complexity, including maximum-entropy models that offer transparent descriptions of collective neural activity but remain largely limited to pairwise interactions and modest system sizes. Here, we use Restricted Boltzmann Machines (RBMs) to model the activity of $\sim1500$-$2000$ simultaneously recorded neurons from the Allen Institute Visual Behavior Neuropixels dataset, spanning multiple cortical and subcortical regions of the mouse brain. RBMs extend the maximum-entropy framework through latent variables, enabling the capture of higher-order dependencies while allowing explicit extraction of effective interaction networks. Recent advances in efficient Markov Chain sampling and training enable accurate learning of these models at this scale. RBMs reproduce the complex statistics of neural recordings with high accuracy. Generated samples match empirical pairwise and higher-order correlations, as well as global statistics such as the distribution of population activity. The inferred parameters provide direct access to effective neuronal interactions, revealing coordination patterns in population activity. These couplings display clear anatomical structure: neurons within visual cortical areas show stronger interactions, consistent with visually driven behavior, while cross-area couplings are weaker. Despite being trained on temporally shuffled data, Markov Chain Monte Carlo simulations also reproduce the global relaxation dynamics of neural activity.

Figures (8)

• Figure 1: Schematic representation of a Restricted Boltzmann Machine.
• Figure 2: Model selection. We plot the evolution of the train and test log-likelihood (resp. in dashed and solid lines) along the training trajectory, for different trainings performed by varying the number of hidden nodes ${N_\mathrm{h}}$. For each curve the star symbol indicates the update at which the maximum test log-likelihood is reached.
• Figure 3: Generated versus real data statistics. Quality of generated samples assessed by comparing empirical statistics of $N=1560$ neurons with those of the training and test sets, using $10{,}000$ randomly selected sequences from each of the training, test, and generated datasets. Panel (a): projection of training data (in black) and generated samples (in red) along the first two principal components (PCA directions) of the training set of data; histograms at the top and at the right highlight the one-dimensional marginal distributions over the first and second direction, respectively. Panel (b): single-neuron firing rate; (c): neuron-neuron covariance matrix; (d): three-body correlations. Red points denote training–test comparisons; green points denote training–generated comparisons. The slope of a linear fit and the Pearson correlation coefficient are reported in the legend. Panel (e): histogram of the number of simultaneously active neurons; (f): histogram of sample energies assigned by the RBM; (g): probability that a neuron is active conditioned on the network state: points are obtained by averaging over bins with similar values of $h_\mathrm{eff}$ (defined in Eq. \ref{['eq:heff']}); error bars correspond to the standard error of the mean.
• Figure 4: Extraction of effective model. Panels (a)-(b)-(c) show the histograms of the effective fields, $2-$body and $3-$body interactions of the the trained RBM, in the Ising-like parametrization. (e) shows the evolution of the norms of each coupling set computed according to Eq. \ref{['eq:frobenius']}. Finally, (d) shows the coarse-grained set of effective 2-body interactions within each brain area and between pairs of brain areas with at least $n=30$ neurons, computed according to Eq. \ref{['eq:Jij_ba_formula']}.
• Figure 5: Interacting network across brain areas. Histogram of effective pairwise couplings $J_{ij}^{(2)}$ ((a)) and of $R_{ij}^{(2)}$ values indicating the involvement of neuron pairs in higher-than-two-body interactions ((b)), computed from the trained RBM and grouped by brain area. Each panel shows the histogram of the a pairwise quantity, $J^{(A)}_{ij}$ or $R^{(A)}_{ij}$, for $i,j \in A$ where $A$ represents one of the brain areas of interest (shown in the legend of each panel together with the number of neurons belonging to it, i.e. $n_A = \sum_i \mathbb{I}\left[ {i\in A} \right]$). For the $J$s, we show in the index the value of the skewness $\gamma$ and the excitation-to-inhibition ratio (E/I) ratio. We only show brain areas having at least $30$ neurons. The RBM is the same as in Fig. \ref{['fig:genstats']}, where its generative performance is evaluated, and Fig. \ref{['fig:effective_model']}, where histograms for the entire brain (i.e., without separating regions) are shown.