Modeling Dynamic Neural Activity by combining Naturalistic Video Stimuli and Stimulus-independent Latent Factors

Thank you for your positive review, highlighting the novelty of our work and our “thorough experiments and analysis”.

Below we address your questions, and we will incorporate these replies into our manuscript:

Q1: You raise a valid point about baselines for the latent factors. The most relevant and commonly used "baseline" for internal brain states are observable behavioral variables, and we performed this exact comparison. However, to address your question, we compared it with CEBRA (Schneider et al., 2023). Please note, that CEBRA does not take visual stimuli as input or have held-out neurons. We trained CEBRA, which has the same latent dimension ($k=12$) as our model, on our data by fitting for each mouse one CEBRA model. Then, we performed the same canonical correlation analysis (CCA) with the CEBRA latents as we did for our latents (lines 322-327). “ Specifically, CCA “finds the linear combination of the latent variables $\mathbf{z_{(1,\dots,T)}}^{(1)},\dots,\mathbf{z_{(1,\dots,T)}}^{(k)}$ with maximal correlation to a chosen behavioral variable like pupil dilation or treadmill speed.” The CCA analysis was done with 5-fold cross-validation in an 80/20 split. Below, we can see the canonical correlation with pupil dilation and treadmill speed on the test set. The error bars indicate the standard error of the mean of cross-validation. Especially for treadmill speed, the latents of our model do correlate stronger with behavior compared to the CEBRA latents, probably because explicitly separating stimulus information actually helps.

Q2: Yes, you are right that mouse V1 has orientation-selective neurons (as well as phase-selective neurons). However, this information should be encoded in the model's core and readout weights and not in the latents [1] as it is a function of the neuron that does not change from trial to trial, changing the response based on the stimuli. However, the key idea of analysis in Table 3 is to show that latents are mainly encoding stimuli independent information. This is expected from the model and we wanted to test against the case that the model uses the latent state to explain variability that it cannot explain with the video-driven part. In theory, it could also do that for orientation, but we’d expect this to be more complex. In addition, orientation is not as easily extracted from natural video as local contrast. We thus focused on contrast because it's a fundamental, low-level stimulus property that can be directly computed from the video frames. This made it a straightforward control that latents do not encode contrast (e.g., stimulus) information - the low R² values in Table 3 confirm this. To further probe stimulus independence beyond contrast, we ran an additional test: we predicted the held‑out neurons for trial A while conditioning on the responses of trial B (same movie, different brain state). If latents mainly carry stimulus information, performance should be unaffected, as stimulus does not change; but if they carry brain-state‑specific information, accuracy should fall. We observed the decrease in performance (Fig. 3, lines 300‑315), confirming that the latents represent a video‑independent brain state.

Q3: Yes, your assumption is correct, and they are chosen randomly. We’ll clarify this in our final version. To answer your question, we ran additional experiments where we masked neurons based on their cortical coordinates. To test how sensitive the model is to that choice, we ran three experiments: we sorted the neurons by their cortical x, y, or z coordinate and then masked the first 50 % of each ordered list. Masking by x or y (cortical surface) caused a small drop in performance, whereas masking by z (cortical depth) left accuracy unchanged relative to random masking. This mirrors Fig. 5: latent patterns show clear topographic organisation across the cortical surface but remain largely constant with depth. Below are the results of the different masking strategies, we reported the standard error of the mean across the five mice.

We will incorporate this analysis into the appendix.

We hope that these additional comments and experiments will help you to understand our method better and revisit the evaluation.

References:

[1] Ustyuzhaninov, I., Burg, M. F., Cadena, S. A., Fu, J., Muhammad, T., Ponder, K., Froudarakis, E., Ding, Z., Bethge, M., Tolias, A. S., & Ecker, A. S. (2022, February 10). Digital twin reveals combinatorial code of non‑linear computations in the mouse primary visual cortex [Preprint]. bioRxiv. https://doi.org/10.1101/2022.02.10.479884

Thank you for your positive review, referring to our paper as “technically sound with a well-formulated probabilistic model” and “clearly written and well-organised”. We also appreciate your constructive feedback. We believe we can address all your concerns and questions:

Addressing raised weaknesses:

W1: Error bars on Table 4: Thanks for pointing this out. We ran the experiments again with models initialized on 2 more seeds and reported the standard error of mean along 3 seeds with the slightly changed results of Table 4 below:

W2 Introducing ZIG: Thanks for your suggestion. We would incorporate an intuitive introduction like this in the text: “In calcium imaging of neurons, the measured fluorescence signals often include many zero or near-zero values when neurons are inactive, along with positively skewed values when activity occurs. A zero-inflated gamma distribution models this by combining a point mass at zero (to capture silent periods) with a gamma distribution (to represent the continuous, positive-valued calcium activity during neural firing). This allows for more accurate modeling of the sparse and bursty nature of neural activity.”

W3 Comparisons with transformer cores: See Q1

W4 Limited impact: We would like to address your point regarding limited impact. First, advances in recording techniques have made large-scale datasets increasingly common, with several open-access resources now available, including MiCRONS [1], Sensorium 2023 [2], and the Allen Brain Observatory [3]. Deep learning based predictive models, shown to be superior in predictive performance, typically rely on such datasets [4], and hence we do not see data scale as a limiting factor for our approach—especially since the data we use is publicly available. Regarding the applicability to smaller-scale datasets, we believe our work can serve as a strong pretrained initialization and could adapt to smaller-scale datasets via fine-tuning. Second, our primary goal was not to outperform baselines on point-estimate metrics like correlation, but to more accurately model the full neuronal response distribution (log-likelihood), while additionally capturing the influence of unobserved latent states alongside stimulus-driven activity—a well-known strong signal in large-scale populations [5]. Lastly, there is evidence that the underlying latent state is effectively high-dimensional [5]. Therefore, it cannot easily be estimated using behavioural variables as proxy. Other existing methods ([3, 35, 40, 55, 56] from the original manuscript) either ignored the stimuli or were not probabilistic and thus ignored the noise correlations, both limiting the extraction and interpretability of the latent state vectors.

Since we show that our learned latents states correlate with behavior, this allows us to design realistic in-silico experiments, supporting in-vivo causal testing. We hence believe our model trained on a large-scale dataset has a direct experimental impact.

W4 Typos: Thanks for catching the typos! We will thoroughly go through the manuscript again.

Addressing raised questions:

Q1: Below, we show that the trends are qualitatively the same - e.g., the ZIG model with the ViT core is a bit worse than Poisson due to the optimization trade-off discussed in lines 199-204, and then the latent model actually improves the performance above the Poisson baseline. While we did not obtain the state-of-the-art performance using ViT during rebuttal---very likely attributable to limited hyperparameter tuning---we hope that this experiment is sufficient to show that the latent method we introduced is transferable towards other cores. Below, we rerun the Table 1 comparison, this time replacing the 3‑D CNN with a ViT backbone for every model.

Additionally, Conditioning on half the neurons further boosts correlation up to 0.194 (0.006), showing that encoding neuron responses works with a different core, too.

Q2: See W1

Q3: Yes, you are right, in the ‘Conditioned’ scenario, we give a subset of responses from the same time and movie frame as input for the model. We fully agree that it would be unfair to compare such a model to a baseline that lacks this extra information. We reported this experiment to show that our model additionally allows us to condition on other neurons and baselines that predict neural activity to stimuli currently cannot. Thus, our main experimental results in Table 1 report the performance of all models without conditioning on any response. See lines 179-180 for evaluation details of this “Non-Conditioned” scenario. In the ‘Non‑Conditioned’ scenario, we still report Pearson correlation, the field’s standard metric, and our video‑only ZIG model’s correlation remains close to the Poisson baseline, while also allowing exact log probability calculation (scoring the full-response distribution) and how including the learned latents improves the log probability.

Regarding existing works that allow conditioning on other neurons, to the best of our knowledge, there are currently no public baselines that support this. The closest is Neuroformer [6], but its original implementation does not support multisession training. Scaling it to a dataset like Sensorium 2023 would require substantial modifications—such as adjusting context window sizes and architectural changes—making it unsuitable as an out-of-the-box baseline.

However, as discussed in W4, our primary goal and the significance of the paper is not to outperform baselines on correlations but to create an accurate and biologically plausible model for the entire neuronal response distribution (log-likelihood) and the latent state. Such models could enable more previously unavailable analyses, such as the cortical topography analysis from Figure 5.

We hope these comments shed some light on our motivation, technical, and design choices. We also hope that it will help the reviewer reconsider some criticism.

References: [1] Zhou, P., Reimer, J., Zhou, D., Pasarkar, A., Kinsella, I., Froudarakis, E., Yatsenko, D. V., Fahey, P. G., Bodor, A., Buchanan, J., Bumbarger, D., Mahalingam, G., Torres, R., Dorkenwald, S., Ih, D., Lee, K., Lu, R., Macrina, T., Wu, J., da Costa, N., Reid, R. C., Tolias, A. S., & Paninski, L. (2020, March 25). EASE: EM‑assisted source extraction from calcium imaging data. bioRxiv. https://doi.org/10.1101/2020.03.25.007468

[2] Turishcheva, P., Fahey, P. G., Hansel, L., Froebe, R., Ponder, K., Vystrčilová, M., Willeke, K. F., Bashiri, M., Willeke, K. F., Bashiri, M. P., Willeke, K. F., Ecker, A. S., Sinz, F. H., Ding, Z., Tolias, A. S., & Ecker, A. S. (2024, July 12). The Dynamic Sensorium competition for predicting large‑scale mouse visual cortex activity from videos. arXiv. https://doi.org/10.48550/arXiv.2305.19654

[3] de Vries, S. E. J., Lecoq, J. A., Buice, M. A., Groblewski, P. A., Ocker, G. K., Oliver, M., … Koch, C. (2020). A large‑scale standardized physiological survey reveals functional organization of the mouse visual cortex. Nature Neuroscience, 23(1), 138–151. https://doi.org/10.1038/s41593-019-0550-5

[4] Lurz, K.‑K., Bashiri, M., Willeke, K. F., Jagadish, A., Wang, E., Walker, E. Y., Cadena, S. A., Muhammad, T., Cobos, E., Tolias, A. S., Ecker, A. S., & Sinz, F. H. (2021, January 12). Generalization in data‑driven models of primary visual cortex, ICLR 2021

[5] Stringer, C., Pachitariu, M., Steinmetz, N., Reddy, C. B., Carandini, M., & Harris, K. D. (2019). Spontaneous behaviors drive multidimensional, brain‑wide population activity. Science, 364(6437), 255. https://doi.org/10.1126/science.aav7893

[6] Antoniades, A., Yu, Y., Canzano, J., Wang, W., & Smith, S. L. (2023). Neuroformer: Multimodal and multitask generative pretraining for brain data. arXiv preprint arXiv:2311.00136.

We appreciate your review and constructive feedback. We believe we can address the raised weakness and question effectively. Specifically:

Addressing raised weakness:

Clarity of related works section: In the final manuscript, we will refine the related works section to more clearly position our contributions within the broader literature. Briefly, our method builds on core-readout models to predict dynamic neural activity to video stimuli, but introduces a probabilistic formulation that also includes learning of latent factors underlying neural activity. Existing probabilistic or latent-state models typically ignore the stimulus altogether, while stimulus-driven models are either deterministic or overlook latent states. We aim to clarify this landscape and clearly situate our work at the intersection of these approaches.

Addressing raised questions:

Q1: Thanks for raising this point. In fact, we had the same question. Therefore, this experiment is already in Appendix G. We keep the 3D CNN and Gaussian read‑out unchanged and train the latent model with a Poisson loss. As we cannot directly compute likelihood on a discrete Poisson distribution for continuous signals, we insert the means of the predicted ZIG distributions into the Poisson loss. For each neuron $i$ and timepoint $t$, our ZIG model predicts the ZIG-parameters $q_{it}(\mathbf{x};\psi)$ and $\theta_{it}(\mathbf{x};\psi)$ (as in Figure 2A) based on the video input $\mathbf{x}$ and its parameters $\psi$. With these, we analytically compute the mean of the ZIG distribution $\hat{r_{it}}\bigl((q_{it}(\mathbf{x};\psi),\theta_{it}(\mathbf{x};\psi)\bigr)$. During training, we then use Poisson Loss between those predicted responses $\hat{r_{it}}$ and the observed neural responses $r_{it}$ as an objective. Basically, we parametrize a mean firing rate (needed for Poisson) via a ZIG distribution. In addition to Appendix G below, we report both correlation and likelihood, showing that a latent model trained with Poisson loss has a worse likelihood performance and a similar correlation performance compared to the ZIG latent model. Moreover, adding a latent improves the likelihood compared to the video-only baselines, regardless of the used loss function.

We hope that these comments helped to understand our paper a bit better, and we will improve the writing based on your feedback.

Thank you for your review, and for viewing our work to have “clear contribution” and “sufficient experimental validation” and providing valuable feedback.

Here we address the weaknesses and questions raised during the review.

Addressing raised weaknesses:

Performance on key metrics: We acknowledge that the ZIG model slightly underperforms the Poisson baseline in lines 240-245 of our main text (L240--245) due to the optimization trade-off, since ZIG is not an exponential family distribution (see [1]), and maximizing its likelihood does not automatically maximize its correlation since mean is not a sufficient statistic of the ZIG distribution [1].
However, our primary goal was not to outperform baselines on point-estimate predictions (correlation) but to more accurately model the entire neuronal response distribution (log-likelihood) along with latent factors, which are known to be a strong signal in neuronal populations [2]. Models that are solely based on sensory input cannot capture these noise correlations. We chose ZIG as a more biologically fitting distribution for deconvolved calcium trace signals, which often include many zero or near-zero values when neurons are inactive, along with positively non-integer skewed values when activity occurs, and our choice also agrees with the literature [3, 4]. Compared to Poisson, which is a discrete distribution, for ZIG, we can also compute the exact likelihood. We chose ZIG for these reasons. While we report the correlation to facilitate comparison with other models, we want to highlight that likelihood is a complete metric as it captures the full stimulus-conditioned distribution rather than correlation, which just captures the conditional mean. Thus, we believe likelihood aligns better with our aim to model both variability and higher moments of neural activity.

Model complexity: We agree that our approach involves a two step procedure—training a standard video-to-neuron model before training the latent module--instead of training everything from scratch, as we show that we get better results this way (reported in Appendix L). However, we do not think that this is necessarily a weakness since using a pretrained backbone is now standard practice in many modern machine learning pipelines, models consisting of multiple separate components are seldom trained fully from scratch. Additionally, we view our two step approach to also provide us with modularity: one can simply use an existing pre-trained video-neuron encoder and fine-tune with the latent module on top, further promoting practical adoption.

Purity of latent states: Our claim about the latent separation stems from the independence structure of the latent variables that follows a factor analysis model (Fig. 2C), where latents are a priori independent of the stimulus – giving a bias toward stimulus independence. Our objective function (ELBO), $E_z [\log p(y \mid x, z)] - \mathrm{KL}\bigl[q(z \mid y) | \mathcal{N}(0,1)\bigr]$, encourages this via the KL term that tells $z$ to remain isotropic unless it helps predict the responses. In the analysis you are referring to, we sample from the approximate posterior $q(z \mid y)$. Due to the collider structure $x \to y \leftarrow z$, the posterior $q(z \mid y)$ can still carry stimulus information. Our control experiment checks whether these samples encode stimulus signals, and we will make this subtlety clearer.

Addressing raised questions:

Q1: Thank you for this suggestion—we had very similar thoughts and in fact ran a series of follow‑up experiments in which we replaced the Gamma component of the ZIG distribution by a more flexible normalizing‑flow-based distribution, while keeping everything else fixed (latent-state and stimulus dependence). The normalizing flow consists of a Gaussian base distribution, followed by a series of invertible transforms, $T$. To capture the peak at zero, the responses below a threshold $\rho$ are still modeled by a uniform distribution as in our original ZIG setup. For each neuron, the flow has learnable parameters. Thus, for responses $y_i>\rho$ of neuron $i$, we apply the flow $T_i$: $p(T_i(y_i)\mid z,x)= \mathcal{N} \left( \mu(x, z), 1 \right)$, where $z$ denotes the latent variable, $x$ the video and $\mu(x, z)$ is the predicted mean of the Gaussian. The actual response likelihood follows from the change--of--variables rule $p(y_i \mid z,x)= p\bigl(T_i(y_i)\mid z,x\bigr) \cdot \bigl|\det\nabla_{y_i}T_i(y_i)\bigr|.$

We also tested replacing the Gaussian base distribution with a Gamma positive‑tail, and a different number of flow layers, ranging from 2-9. The flow layers alternated affine with non‑linear functions, including {log⁡,exp⁡,tanh⁡,softplus,ELU}. Our best combination, a two‑layer $\text{affine}\rightarrow \text{softplus}^{-1}$ flow with a Gaussian base distribution, could not improve upon our ZIG model (see below).

Thus, we concluded that the ZIG distribution is already a pretty good fit to the neuron responses. Further, we concluded that it is tricky to close the approximation gap- the trade-off between correlation and likelihood- that arises for non-exponential distributions [1].

We will integrate these experimental results as part of the appendix.

Q2: Thank you for the interesting suggestion. We agree that a model producing discrete spikes and followed by convolutional calcium dynamics would mimic the biological process and data collection better and is scientifically interesting. We tested this briefly in the past. However, because we could not make it work well, we abandoned this idea.

We hope these comments make our goals and choices clearer and help address some of the concerns raised.

References:

[1] Konstantin-Klemens Lurz, Mohammad Bashiri, Edgar Y. Walker, and Fabian H. Sinz. (2023) Bayesian oracle for bounding information gain in neural encoding models. In The Eleventh International Conference on Learning Representations.

[2] Stringer, C., Pachitariu, M., Steinmetz, N., Reddy, C. B., Carandini, M., & Harris, K. D. (2019). Spontaneous behaviors drive multidimensional, brain‑wide population activity. Science, 364(6437), 255. https://doi.org/10.1126/science.aav7893

[3] Zhou, Ding, and Xue-Xin Wei. "Learning identifiable and interpretable latent models of high-dimensional neural activity using pi-VAE." Advances in neural information processing systems 33 (2020): 7234-7247.

[4] Zhu, F., Grier, H. A., Tandon, R., Cai, C., Agarwal, A., Giovannucci, A., ... & Pandarinath, C. (2022). A deep learning framework for inference of single-trial neural population dynamics from calcium imaging with subframe temporal resolution. Nature neuroscience, 25(12), 1724-1734.