Accurate Identification of Communication Between Multiple Interacting Neural Populations

We thank the reviewer for their thoughtful and constructive feedback. We’re encouraged by the reviewer’s comments highlighting MR-LFADS as a novel extension of LFADS, and by the positive assessments on our synthetic dataset. Below, we respond to the reviewer's concerns on the lack of comparison to experimental data. Rebuttal figures referenced throughout our response can be accessed at the following link: https://drive.google.com/drive/folders/1R94up1vl04bkkE12tpT80YeHo39vEbqw?usp=share_link

We agree that evaluation on experimental data is essential. We will update our paper to include demonstrations of MR-LFADS on multi-region electrophysiology from 5 simultaneously recorded Neuropixel probes in mice (Chen et al., Cell, 2024) performing a decision making task (Fig. 2a in our rebuttal figures PDF). We show that MR-LFADS predicts experimentally perturbed firing rates across multiple brain regions following photoinhibition of the anterior lateral motor (ALM) cortex (Fig. c-e). Critically, MR-LFADS was not trained on photoinhibition trials. Thus, MR-LFADS' ability to predict photoinhibition effects suggests identification of an accurate model of inter-region communication. Additionally, we show that MR-LFADS(R) yields more consistent solutions across random initializations (seeds) on these experimental data, as compared to the under-constrained MR-LFADS(G) (Fig. 2f,h).

• Photoinhibition study: We applied MR-LFADS-R to a $5$-region network spanning ALM, thalamus(A), thalamus(O), midbrain reticular nucleus (MRN), and superior colliculus (SC). Thalamus(A) includes the ventral medial (VM) and ventral anterior lateral (VAL) nuclei, which are strongly reciprocally connected with ALM (Guo et al., Nature, 2017), while Thalamus(O) includes other thalamic sub-areas, such as anterior ventral (AV) and lateral dorsal (LD) nuclei. In a subset of trials, ALM is briefly photoinhibited (Fig. 2a,c). We only fit MR-LFADS to unperturbed "control" trials and reserve the photoinhibition trials for post-training validation.

  We first confirm that MR-LFADS fits held-out control trials well (Fig. 2b). We then simulate photoinhibition in MR-LFADS by setting the MR-LFADS ALM-region rates to be the trial-average time-varying firing rates from the photoinhibition trials from the experiment (Fig. 2d), denoted $r^\text{inh}_*$.

  These perturbed firing rates then propagate through MR-LFADS via inferred messages, ultimately predicting photoinhibited firing rates $\hat{r}_*(\hat{m}^\text{inh})$. We find that the predicted ordering of ALM photoinhibition influence aligns with the experimental data (Fig. 2e).

• Consistency across seeds: We fit MR-LFADS-R and MR-LFADS-G to a $3$-region network (ALM, Thalamus, MRN). Both models produced consistent effectomes across different random seeds (Fig. 2f), but MR-LFADS-R exhibited significantly greater consistency in message content (Fig. 2g,h).

We thank the reviewer for the thoughtful and insightful comments. We are encouraged that our model was recognized as a valuable contribution and that the key claims were seen as well supported. Below, we respond to each point raised. Rebuttal figures referenced throughout our response can be accessed at the following link: https://drive.google.com/drive/folders/1R94up1vl04bkkE12tpT80YeHo39vEbqw?usp=share_link

1. On number of regions: We agree that evaluation on larger networks is important for assessing relevance to neuroscience applications. We will invest in improving MR-LFADS infrastructure to support larger networks, and we will provide these ongoing improvements to the community by releasing and continually updating an MR-LFADS code package. To provide immediate neuroscience relevance, see point 2 below.

2. On real-world applicability: We agree that evaluation on experimental data is essential. We will update our paper to include demonstrations of MR-LFADS on multi-region electrophysiology from 5 simultaneously recorded Neuropixel probes in mice (Chen et al., Cell, 2024) performing a decision making task (Fig. 2a in our rebuttal figures PDF). We show that MR-LFADS predicts experimentally perturbed firing rates across multiple brain regions following photoinhibition of the anterior lateral motor (ALM) cortex (Fig. c-e). Critically, MR-LFADS was not trained on photoinhibition trials. Thus, MR-LFADS' ability to predict photoinhibition effects suggests identification of an accurate model of inter-region communication. Additionally, we show that MR-LFADS(R) yields more consistent solutions across random initializations (seeds) on these experimental data, as compared to the under-constrained MR-LFADS(G) (Fig. 2f,h).

• Photoinhibition study: We applied MR-LFADS-R to a $5$-region network spanning ALM, thalamus(A), thalamus(O), midbrain reticular nucleus (MRN), and superior colliculus (SC). Thalamus(A) includes the ventral medial (VM) and ventral anterior lateral (VAL) nuclei, which are strongly reciprocally connected with ALM (Guo et al., Nature, 2017), while Thalamus(O) includes other thalamic sub-areas, such as anterior ventral (AV) and lateral dorsal (LD) nuclei. In a subset of trials, ALM is briefly photoinhibited (Fig. 2a,c). We only fit MR-LFADS to unperturbed "control" trials and reserve the photoinhibition trials for post-training validation.

  We first confirm that MR-LFADS fits held-out control trials well (Fig. 2b). We then simulate photoinhibition in MR-LFADS by setting the MR-LFADS ALM-region rates to be the trial-average time-varying firing rates from the photoinhibition trials from the experiment (Fig. 2d), denoted $r^\text{inh}_*$.

  These perturbed firing rates then propagate through MR-LFADS via inferred messages, ultimately predicting photoinhibited firing rates $\hat{r}_*(\hat{m}^\text{inh})$. We find that the predicted ordering of ALM photoinhibition influence aligns with the experimental data (Fig. 2e).

• Consistency across seeds: We fit MR-LFADS-R and MR-LFADS-G to a $3$-region network (ALM, Thalamus, MRN). Both models produced consistent effectomes across different random seeds (Fig. 2f), but MR-LFADS-R exhibited significantly greater consistency in message content (Fig. 2g,h).

3. On inclusion of Granger causality (GC) and Dynamic Causal Modeling (DCM): We agree that both GC and DCM offer valuable approaches for modeling inter-region interactions. While we were very interested in including these methods, a careful and fair implementation -- particularly of the many nonlinear variants of GC and the more involved generative modeling in DCM -- was beyond the time frame of this rebuttal. Rather than risk an incomplete comparison, we will update our Related Work section to highlight these approaches, and we will pursue them in our ongoing work in this space.

4. On spurious connections: We agree that it would be helpful to further analyze MR-LFADS’s behavior -- specifically, whether the model tends to miss true connections or infer spurious ones, as the reviewer suggests. We have added the following to our paper: "It would also be interesting, from an interpretability perspective, to assess whether MR-LFADS tends to under- or over-estimate connectivity -- that is, whether it is more likely to miss true connections or to infer spurious ones -- and if these mismatches are driven by particular statistical properties of the signals or the dynamics of the downstream areas.”

   Regarding the reviewer’s question about Exp.1: yes, the model does infer weak connections even when no ground truth connectivity is present. However, these inferred message norms are two orders of magnitude smaller than those associated with true connections. To illustrate this, we have added heatmap visualizations showing the full distribution of inferred message norms (Fig. 3).

We thank the reviewer for their thorough and constructive feedback. We appreciate the positive assessment of the manuscript’s contribution, clarity, and experimental design. Below, we respond to each of the points raised. Rebuttal figures referenced throughout our response can be accessed at the following link: https://drive.google.com/drive/folders/1R94up1vl04bkkE12tpT80YeHo39vEbqw?usp=share_link

1. On Poisson outputs: We thank the reviewer for this helpful suggestion. We agree that it is important to demonstrate that the model can operate on spike train data, not just continuous activity. To this end, we applied MR-LFADS to large-scale electrophysiology recordings from Chen et al. (Cell, 2024) (Fig. 2a in our rebuttal figures PDF), which consist of spiking data and were modeled using a Poisson output distribution (Fig. 2b). We elaborate on these results below.

2. On real-world applicability: We agree that evaluation on experimental data is essential. We will update our paper to include demonstrations of MR-LFADS on multi-region electrophysiology from 5 simultaneously recorded Neuropixel probes in mice (Chen et al., Cell, 2024) performing a decision making task (Fig. 2a). We show that MR-LFADS predicts experimentally perturbed firing rates across multiple brain regions following photoinhibition of the anterior lateral motor (ALM) cortex (Fig. c-e). Critically, MR-LFADS was not trained on photoinhibition trials. Thus, MR-LFADS' ability to predict photoinhibition effects suggests identification of an accurate model of inter-region communication. Additionally, we show that MR-LFADS(R) yields more consistent solutions across random initializations (seeds) on these experimental data, as compared to the under-constrained MR-LFADS(G) (Fig. 2f,h).

• Photoinhibition study: We applied MR-LFADS-R to a $5$-region network spanning ALM, thalamus(A), thalamus(O), midbrain reticular nucleus (MRN), and superior colliculus (SC). Thalamus(A) includes the ventral medial (VM) and ventral anterior lateral (VAL) nuclei, which are strongly reciprocally connected with ALM (Guo et al., Nature, 2017), while Thalamus(O) includes other thalamic sub-areas, such as anterior ventral (AV) and lateral dorsal (LD) nuclei. In a subset of trials, ALM is briefly photoinhibited (Fig. 2a,c). We only fit MR-LFADS to unperturbed "control" trials and reserve the photoinhibition trials for post-training validation.

  We first confirm that MR-LFADS fits held-out control trials well (Fig. 2b). We then simulate photoinhibition in MR-LFADS by setting the MR-LFADS ALM-region rates to be the trial-average time-varying firing rates from the photoinhibition trials from the experiment (Fig. 2d), denoted $r^\text{inh}_*$.

  These perturbed firing rates then propagate through MR-LFADS via inferred messages, ultimately predicting photoinhibited firing rates $\hat{r}_*(\hat{m}^\text{inh})$. We find that the predicted ordering of ALM photoinhibition influence aligns with the experimental data (Fig. 2e).

• Consistency across seeds: We fit MR-LFADS-R and MR-LFADS-G to a $3$-region network (ALM, Thalamus, MRN). Both models produced consistent effectomes across different random seeds (Fig. 2f), but MR-LFADS-R exhibited significantly greater consistency in message content (Fig. 2g,h).

3. On notational issues and typos: We greatly appreciate the reviewer’s attention to detail. We have corrected the identified issues and revised the notation for the priors over inferred inputs and initial conditions in the updated manuscript (Fig. 4). Typographical errors on lines 28 and 633 have also been fixed. Additionally, the reviewer was correct to point out the inconsistency in the KL divergence description on line 211 -- we have updated this to reflect the correct direction of the penalty: "Therefore, we carefully set the KL penalty weight for inferred inputs ($\lambda_\text{KL}^{\hat{u}}$) to be higher than that for messages ($\lambda_\text{KL}^{\hat{m}}$)."

4. On clarifications for incorrect inferences due to incomplete recordings: We appreciate this thoughtful question from the reviewer. The interpretation offered — “I suppose you mean that if you were to add more complete recordings, the inferred connectivity would change” — captures our intended meaning precisely. In practice, we often cannot guarantee that all relevant areas are recorded. As a result, inputs from unobserved regions must be accounted for (as latent variable inferred inputs) to minimize the deleterious effects of unobserved nuisance variables. We will update our paper to clarify this point.

We thank the reviewer for their thoughtful and constructive feedback. We also appreciate the reviewer’s positive assessment of our effort to address limitations in common architectural choices for communication models. Below, we respond to each of the points raised. Rebuttal figures referenced throughout our response can be accessed at the following link: https://drive.google.com/drive/folders/1R94up1vl04bkkE12tpT80YeHo39vEbqw?usp=share_link

1. On the use of Jensen-Shannon divergence (JSD): Indeed, JSD is classically defined over probability distributions. We used it heuristically by treating the normalized message-norm vectors as categorical distributions and used JSD to assess how closely they match in shape. We agree that this is not the most principled application of JSD. We have now re-evaluated our approach using cosine similarity, $S_{\cos}$, as an alternative metric for comparing inferred message-norm vectors to ground truth. We find that all model comparisons and trends remain consistent when using cosine similarity. We will update our paper, exchanging JSD for cosine similarity throughout. We have provided the corresponding comparisons between JSD and cosine similarity in Figure 1 of our rebuttal figures PDF.

2. On inferred latent inputs in MR-SDS: We appreciate the reviewer's insight into inferred latent noise in MR-SDS (Karniol-Tambour el al., ICLR, 2024). Eq 4 of the MR-SDS paper states: $x^k_t = f_{kk}^{z_t}(x^k_{t-1}) + \ldots + f_{ku}^{z_t}(u_t) + \epsilon^k_t,\quad \epsilon^k_t \sim \mathcal{N}(0, Q_k^{z_t}).$ In the absence of provided external inputs, $u_t$, the components that could account for unobserved inputs are the dynamics, $f_{kk}^{z_t}$, or the covariance of the noise term, $Q_k^{z_t}$. Since the dynamics function is agnostic to external inputs, as the reviewer points out, $Q_k^{z_t}$ is the only parameter that can reflect the influence of unobserved inputs, potentially through a low-rank parameterization to reflect low-dimensional inputs. However, this term captures only the variance of the inferred inputs, ignoring the time-varying mean component. That is, the inputs in LFADS are parameterized as $\hat{u}_t \sim \mathcal{N}(u_t | \mu_t^u, \sigma_t^u)$ (Eq 23 in the LFADS paper), with both $\mu_t^u$ and $\sigma_t^u$ being inferred. Therefore, it is not entirely clear to us whether the noise in MR-SDS serves an identical role to the inferred inputs in LFADS. We would greatly appreciate further clarification from the reviewer in case we are misunderstanding important details here.

3. On generative model notation (Eq. 5–7):
   We appreciate the reviewer’s attention to this detail and we appologize for the confusion. We do have a proper generative model, and we will update our paper to clarify as much. In our submission, Equation (6) describes the prior over the messages (and other latent variables), and Equation (7) provides the approximate posterior over the messages (which is used for inference and not if generating data from the model). The latter was not clear in our initial submission, and we will clarify this in the revised paper. In addition, we have updated Equation (5) in our paper: $\hat{x}_t^i \sim \mathcal{N}(\mu_t^{x_i}, \mathrm{diag}(\sigma_t^{x_i}))$, to indicate the diagonal covariance, as suggested.

4. On effectome visualizations: We thank the reviewer for this helpful suggestion. We agree that the thresholded diagrams serve primarily a qualitative purpose. In response, we have added unthresholded heatmap visualizations (as in Fig. 2e of our paper) to provide a more rigorous presentation of the inferred effectome (Figure 3 of our rebuttal figures PDF).