Improved inference of latent neural states from calcium imaging data

Improved inference of latent neural states from calcium imaging data | bioRxiv /* */ /* */ <!-- <!-- /*! * yepnope1.5.4 * (c) WTFPL, GPLv2 */ (function(a,b,c){function d(a){return"[object Function]"==o.call(a)}function e(a){return"string"==typeof a}function f(){}function g(a){return!a||"loaded"==a||"complete"==a||"uninitialized"==a}function h(){var a=p.shift();q=1,a?a.t?m(function(){("c"==a.t?B.injectCss:B.injectJs)(a.s,0,a.a,a.x,a.e,1)},0):(a(),h()):q=0}function i(a,c,d,e,f,i,j){function k(b){if(!o&&g(l.readyState)&&(u.r=o=1,!q&&h(),l.onload=l.onreadystatechange=null,b)){"img"!=a&&m(function(){t.removeChild(l)},50);for(var d in y[c])y[c].hasOwnProperty(d)&&y[c][d].onload()}}var j=j||B.errorTimeout,l=b.createElement(a),o=0,r=0,u={t:d,s:c,e:f,a:i,x:j};1===y[c]&&(r=1,y[c]=[]),"object"==a?l.data=c:(l.src=c,l.type=a),l.width=l.height="0",l.onerror=l.onload=l.onreadystatechange=function(){k.call(this,r)},p.splice(e,0,u),"img"!=a&&(r||2===y[c]?(t.insertBefore(l,s?null:n),m(k,j)):y[c].push(l))}function j(a,b,c,d,f){return q=0,b=b||"j",e(a)?i("c"==b?v:u,a,b,this.i++,c,d,f):(p.splice(this.i++,0,a),1==p.length&&h()),this}function k(){var a=B;return a.loader={load:j,i:0},a}var l=b.documentElement,m=a.setTimeout,n=b.getElementsByTagName("script")[0],o={}.toString,p=[],q=0,r="MozAppearance"in l.style,s=r&&!!b.createRange().compareNode,t=s?l:n.parentNode,l=a.opera&&"[object Opera]"==o.call(a.opera),l=!!b.attachEvent&&!l,u=r?"object":l?"script":"img",v=l?"script":u,w=Array.isArray||function(a){return"[object Array]"==o.call(a)},x=[],y={},z={timeout:function(a,b){return b.length&&(a.timeout=b[0]),a}},A,B;B=function(a){function b(a){var a=a.split("!"),b=x.length,c=a.pop(),d=a.length,c={url:c,origUrl:c,prefixes:a},e,f,g;for(f=0;f<d;f++)g=a[f].split("="),(e=z[g.shift()])&&(c=e(c,g));for(f=0;f<b;f++)c=x[f](c);return c}function g(a,e,f,g,h){var i=b(a),j=i.autoCallback;i.url.split(".").pop().split("?").shift(),i.bypass||(e&&(e=d(e)?e:e[a]||e[g]||e[a.split("/").pop().split("?")[0]]),i.instead?i.instead(a,e,f,g,h):(y[i.url]?i.noexec=!0:y[i.url]=1,f.load(i.url,i.forceCSS||!i.forceJS&&"css"==i.url.split(".").pop().split("?").shift()?"c":c,i.noexec,i.attrs,i.timeout),(d(e)||d(j))&&f.load(function(){k(),e&&e(i.origUrl,h,g),j&&j(i.origUrl,h,g),y[i.url]=2})))}function h(a,b){function c(a,c){if(a){if(e(a))c||(j=function(){var a=[].slice.call(arguments);k.apply(this,a),l()}),g(a,j,b,0,h);else if(Object(a)===a)for(n in m=function(){var b=0,c;for(c in a)a.hasOwnProperty(c)&&b++;return b}(),a)a.hasOwnProperty(n)&&(!c&&!--m&&(d(j)?j=function(){var a=[].slice.call(arguments);k.apply(this,a),l()}:j[n]=function(a){return function(){var b=[].slice.call(arguments);a&&a.apply(this,b),l()}}(k[n])),g(a[n],j,b,n,h))}else!c&&l()}var h=!!a.test,i=a.load||a.both,j=a.callback||f,k=j,l=a.complete||f,m,n;c(h?a.yep:a.nope,!!i),i&&c(i)}var i,j,l=this.yepnope.loader;if(e(a))g(a,0,l,0);else if(w(a))for(i=0;i (function(w,d,s,l,i){w[l]=w[l]||[];w[l].push({'gtm.start':new Date().getTime(),event:'gtm.js'});var f=d.getElementsByTagName(s)[0];var j=d.createElement(s);var dl=l!='dataLayer'?'&l='+l:'';j.src='//www.googletagmanager.com/gtm.js?id='+i+dl;j.type='text/javascript';j.async=true;f.parentNode.insertBefore(j,f);})(window,document,'script','dataLayer','GTM-M677548'); Skip to main content Home About Submit ALERTS / RSS Search for this keyword Advanced Search New Results Improved inference of latent neural states from calcium imaging data View ORCID Profile Stephen Keeley , David Zoltowski , View ORCID Profile Adam Charles , View ORCID Profile Jonathan Pillow doi: https://doi.org/10.1101/2025.10.17.682993 Stephen Keeley 1 Department of Natural Sciences, Fordham University Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for Stephen Keeley For correspondence: skeeley1{at}Fordham.edu David Zoltowski 2 Department of StatiWu Tsai Neurosciences Institute, Stanford University Find this author on Google Scholar Find this author on PubMed Search for this author on this site Adam Charles 3 Department of Biomedical Engineering, Johns Hopkins University Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for Adam Charles Jonathan Pillow 4 Princeton Neuroscience Institute, Princeton University Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for Jonathan Pillow Abstract Full Text Info/History Metrics Preview PDF Abstract Calcium imaging (CI) is a standard method for recording neural population activity, as it enables simultaneous recording of hundreds-to-thousands of individual somatic signals. Accordingly, CI recordings are prime candidates for population-level latent variable analyses, for example using models such as Gaussian Process Factor Analysis (GPFA), hidden Markov models (HMMs), and latent dynamical systems. However, these models have been primarily developed and fine-tuned for electrophysiological measurements of spiking activity. To adapt these models for use with the calcium signals recorded with CI, per-neuron fluorescence time-traces are typically either de-convolved to approximate spiking events or analyzed directly under Gaussian observation assumptions. The former approach, while enabling the direct application of latent variable methods developed for spiking data, suffers from the imprecise nature of spike estimation from CI. Moreover, isolated spikes can be undetectable in the fluorescence signal, creating additional uncertainty. A more direct model linking observed fluorescence to latent variables would account for these sources of uncertainty. Here, we develop accurate and tractable models for characterizing the latent structure of neural population activity from CI data. We propose to augment HMM, GPFA, and dynamical systems models with a CI observation model that consists of latent Poisson spiking and autoregressive calcium dynamics. Importantly, this model is both more flexible and directly compatible with standard methods for fitting latent models of neural dynamics. We demonstrate that using this more accurate CI observation model improves latent variable inference and model fitting on both CI observatons generated using state-of-the-art biophysical simulations as well as imaging data recorded in an experimental setting. We expect the developed methods to be widely applicable to many different analysis of population CI data. Introduction Electrophysiological recordings have historically been the de-facto approach for observing the activity of single neurons. Consequently, many statistical models for neuronal data at the single-cell level are formulated for spike train observations. That is, they describe a mapping from stimuli or unobserved latent variables to a probability distribution over discrete spike events. This class of models includes (but is not limited to) Poisson regression models ( Truccolo et al.,2005 ; Pillow et al.,2008 ; McFarland et al.,2013 ; Park et al.,2013 ; Zoltowski and Pillow, 2018 ), non-Poisson spike count regression models ( Pillow and Scott, 2012 ; Goris et al.,2014 ; Williamson et al.,2015 ; Gao et al.,2015 ; Linderman et al.,2016a ; Stevenson, 2016 ; Charles et al.,2018 ), linear dynamical system (LDS) ( Macke et al.,2011 ), nonlinear dynamics models ( Linderman et al.,2017 ; Pandarinath et al.,2018 ; Duncker et al.,2019 ; Zhao and Park, 2019 ), Gaussian process factor analysis (GPFA) models ( Yu et al.,2009a ; Zhao and Park, 2017 ; Keeley et al.,2019 ), and nonlinear tuning curve models ( Zhang et al.,1998 ; Zemel et al.,1998 ; Cronin et al.,2010 ; Rad and Paninski, 2011 ; Calabrese et al.,2011 ; Park et al.,2014 ; Savin and Tkacik, 2016 ; Rad et al.,2017 ; Wu et al.,2017 ). Calcium imaging (CI) is a popular approach for recording large neural populations due to its ability to image large areas ( > 0.5 μm 2 ) at micron level resolution, enabling simultaneously recording many neurons (100-1000 typical, 10 9 in the most advanced systems ( Demas et al.,2021 )) and tracking the same population of cells for days at a time. CI, however, measures neural spiking activity indirectly. Specifically, fluorescence changes in recorded CI time series represent fluctuations in intra-cellular calcium concentrations that result from the biophysics of action potentials; each spike results in a rapid rise in the calcium concentration. The fluorescence then jumps as the calcium is bound to the fluorescent proteins followed by a slower decay as the calcium unbinds ( Song et al.,2021 ; Helmchen and Tank, 2015 ). While theoretically the relationship between neural spiking and calcium-based fluorescence is well characterized, practically variability in concentrations and noise considerations complicates the ability to discern single spikes (or even small bursts of 2-3 spikes) from calcium traces ( Song et al.,2021 ; Ledochowitsch et al.,2019 ). Due to the indirect relationship between neural firing and CI traces, point-process models are no longer directly applicable to CI data. Instead methods to statistically analyze CI datasets commonly take two approaches: (1) ignore spiking and resort to Gaussian noise models operating directly on the calcium traces; or (2) apply Poisson models to estimates of the spike train obtained from calcium inference methods ( Smith and Häusser, 2010 ; Ko et al.,2013 ; Pnevmatikakis et al.,2016 ). The former approach is suboptimal because the statistics of CI data are asymmetric with high-skew and long tails, and thus widely deviate from the statistics of Gaussian distributions ( Wei et al.,2020 ). Additionally, i . i . d . noise observations are not appropriate for the long time-scale autocorrelations observed in CI data. The latter approach, spike-time estimation, is suboptimal because it does not take into account the uncertainty of the unobserved spikes. Specifically, single spikes have been shown to only be visible in the CI traces a fraction of the time ( Huang et al.,2021 ), and the nonlinearities in calcium buffering in bursting can make spike estimation highly unreliable in many settings. Our approach instead uses an observation model that is more faithful to the data-generating process, which we demonstrate provides more accurate inference and scientific insight in use. Specifically, here we extend a calcium observation likelihood first presented in Ganmor et al. (2016 ) to latent variable models where the firing rate is a function of the latent variables. The three primary latent variable models we consider are hidden Markov models (HMMs) ( Smith and Brown, 2003 ; Escola et al.,2011 ; Krause and Drugowitsch, 2022 ), Gaussian Process Factor Analysis (GPFA) ( Yu et al.,2009a ), and Latent Factor Analysis via Dynamical Systems (LFADS) ( Pandarinath et al.,2018 ), though our approach also applies to switching dynamical systems ( Linderman et al.,2017 ) and more general nonlinear dynamics models ( Mudrik et al.,2024 , 2025 ; Yezerets et al.,2025 ; Chen et al.,2024 ). The models each are fit using different methods, including maximum likelihood (ML) for the HMM, and variational inference for GPFA and LFADS, which demonstrates the range of methods that can be adapted with our approach. We evaluate the models using simulated datasets, a state-of-the-art biophysical simulator NAOMi ( Song et al.,2021 ), and in-vivo 2-photo calcium imaging recordings. In each of these cases, we show that the calcium likelihood can be used in conjunction with various LVMs is able to better capture the underlying latent structure in the simulated or experimental neural population traces than competing observation likelihoods. Results Calcium LVM Framework We propose the following framework for adapting LVMs developed for spiking data to the setting of calcium imaging recordings y t ∈ ℝ N of N neurons for t ∈{1, …, T }. We consider models that prescribe a generative distribution over Poisson spike counts s t ∈ ℕ 0 N via latent variables x where the latent variables x may be continuous or discrete and the firing rate λ θ depends on the value of the latent variables, potentially through a mapping parameterized by θ . This formulation includes many common LVMs used in neuroscience that we will describe in detail in the following sections. Previously, Ganmor et al. (2016 ) proposed a model for fitting firing rates from calcium traces via approximate marginalization over unobserved Poisson spike counts. Here we propose to augment the generative models over spike counts with this model such that for a given neuron n where we have generalized the approach of Ganmor et al. (2016 ) to the setting of higher-order AR processes. This model has three sets of parameters per-neuron. The AR coefficients α i account for autoregressive calcium dynamics. If the model is AR (1),then α 1 determines the exponential decay of fluorescence. Next, c describes the fluorescence increase due to a single spike and σ 2 is the variance of additive Gaussian noise. Download figure Open in new tab Figure 1. Schematic for the calcium LVMs. A low-dimensional latent variable maps to Poisson firing rates for time-series neural population activity. These Poisson spike probabilities are marginalized over and fed through an auto-regressive process to describe the evolution of the calcium traces. Importantly, in this model the unobserved spike count variables can be approximately marginalcross-validation (60 trainingi zed out via numerical integration such that it is tractable to evaluate . Accordingly, with automatic differentiation approaches and current approximate inference methods this observation model can be generally applied across a variety of neural population models. In the following sections we will detail three specific applications of this framework. Calcium Hidden Markov Models Hidden Markov Models (HMMs) are latent sequence models for identifying discrete structure over time. Here we describe the calcium HMM, an HMM with a more appropriate observation distribution for identifying discrete sequential structure in calcium imaging recordings. The generative model for a population of neurons is where z ∈{1, …, K } is a discrete variable taking on one of K values, ⊙ indicates element wise multiplication, and Ψ is a diagonal matrix of per-neuron noise variances. The parameters π 0 and π k correspond to the initial discrete state distribution and transition distribution conditioned on state k , respectively. Here, we have rewritten the per-neuron observation model in Equation (3) in a vectorized form for an AR(1) process. In the model, each discrete state prescribes a separate Poisson spike rate for every neuron. The spike counts generated from the model are then input to the AR observation process, whose parameters are shared across states. The model is fit via the Expectation-Maximization (EM) algorithm Dempster et al. (1977 ). The spike counts are numerically marginalized to get state-dependent likelihoods p ( y t ∣ y t −1 , z t = k ) used in the E-step. In the M-step, the expected joint log likelihood is optimized with respect to the model parameters using automatic differentiation. Simulated Data We demonstrate the potential benefits of the Calcium HMM in a simulated experiment where a population of neurons exhibits sequential firing activity ( Fig. 2a-b ). The data is simulated from an HMM with 5 states and a separate cluster of neurons has high firing rates in each state. The transition matrix enforces sequential transitions between clusters. We first simulated the spike counts and then simulated the calcium observations given the spike counts. The calcium observations were generated from Equation (21) with additional independent measurement noise per time point. Importantly, this additional noise is biophysically motivated, meaning that the data do not fully arise from any of the models we will consider and are more faithful to real data collected in neuroscience experiments. Download figure Open in new tab Figure 2. Simulated calcium HMM. (a) Schematic for the overall calcium-HMM model. (b) The simulated data and (c) the inferred and true discrete states of the underlying latent states using different calcium models. (d) the inferred test log-likeilhoods as we vary the number of discrete latent states for different models. (e) true neural rates for each state of the underlying data alongside the inferred AR parameter mean (middle) and expected calcium fluorescence value (right). We fit four different models to the simulated data. First, we fit a Poisson HMM with five states to the underlying spike counts. This comparison point quantifies how informative the underlying spikes are for recovering the model parameters and latent states, providing a reasonable upper bound on performance. We next fit three HMMs with different observation models to the generated calcium imaging data. The three observation models are independent Gaussian, autoregressive (AR) Gaussian, and the calcium observation model from Equation (21). The AR Gaussian model can be thought of as a special case of the calcium model that ablates the spiking component to test the relative utility of including the discrete spiking. We found that modeling the calcium imaging data as an autoregressive model with a latent discrete spiking component best explained the simulated data. The recovered latent states from the calcium HMM on test data (correlation between true discrete states and inferred discrete states ρ =0.89) more closely matched the true latent states than those inferred from the Gaussian HMM ( ρ =0.53) or autoregressive HMM ( ρ =0.73), and approached the performance of the Poisson HMM ( ρ =0.91) fit to the spiking data ( Figure 2c ). The calcium HMM achieved the test highest log-likelihood out of the models fit to the simulated calcium data and correctly identified the number of latent discrete states ( Fig. 2d ). Finally, for each discrete state we compared the true neuron firing rates with the “calcium influx” c ⊙ λ z , which are the firing rates inferred from the Calcium HMM model scaled by the influx parameters. We find that the estimated calcium influx approximately matches the true simulated firing rates and correctly identifies the different clusters of neurons, which is in contrast to the AR model that learns a less clear clustering ( Fig. 2e ). Modeling Piriform Cortex Recordings During Odor Presentation After validating calcium HMM in simulation, we next demonstrate fitting the model to calcium imaging responses recorded in mouse piriform cortex during passive odor presentation ( Daste and Pierré, 2022 ; Srinivasan et al.,2023 ). In this dataset, ten different odors were presented eight times over the course of a session containing 80 trials. The duration of each trial was 30 seconds and the odor was passively presented for 1 second starting 10 s after the trial onset. We fit HMMs with calcium or Gaussian observation models to the responses of 284 neurons in piriform cortex recorded via calcium imaging ( Fig. 3a ). We first selected the number of discrete states K for each model via cross-validation (60 training trials, 19 test trials; K =13 for calcium and K =12 for Gaussian). The test log likelihood for the calcium HMM is substantially higher than the Gaussian observation model due to the more expressive observation model ( Fig. 3b ). Download figure Open in new tab Figure 3. HMM comparison on odor response data. The calcium HMM identifies an odor onset state that is more tightly coupled with the actual odor onset. (a) Example Δ F / F traces for the population of recorded neurons on one trial. (b) Test log likelihoods for calcium and Gaussian HMMs as a function of the number of discrete states. Arrows indicate the number of discrete states with the highest test log likelihood. (c) Inferred most likely states on both training and test trials for each model. Each model identifies a consistent “odor onset” state linked to the time of odor presentation at 10s. (d) The fraction of trials in the odor onset state at each time point for each model. The calcium HMM odor onset state peaks more closely to the odor presentation window and has a shorter width (black arrows denote calculation of width). We next analyzed the resulting calcium and Gaussian HMMs given the optimal number of states for each model. In each model, we found a discrete state that appears locked to odor onset across most trials ( Fig. 3c ). This “odor onset” state appears to mark a transition in the population response that occurs across most odor presentations. Notably, we found that the temporal extent of the odor onset state for the calcium HMM was more closely locked to the odor delivery window than the odor onset state for the Gaussian HMM ( Fig. 3d ). For the calcium HMM, trials generally transitioned into the odor onset state shortly after the odor presentation started and transitioned out of the odor onset state after a couple of seconds (peak onset state occupancy at 10.63s and onset state occupancy width of 2.65s). However, for the Gaussian HMM the transition into the odor onset state was generally delayed relative to the odor presentation and had a longer duration (peak onset state occupancy at 12.85s and onset state occupancy width of 6.64s). The differences in the inferred odor onset state highlight the potential utility of the calcium LVMs described in this paper. In this dataset, the calcium HMM identified a discrete state that was time-locked to a behavioral variable of interest (odor presentation) in terms of both onset and duration of the discrete state. The calcium HMM odor onset state is more consistent with a population state that is time locked to the odor delivery window than the Gaussian HMM odor onset state. These results match our intuitions from the simulated calcium HMM example, where the inferred discrete states using the calcium HMM closely match the true discrete states but the inferred discrete states using the Gaussian HMM are generally delayed relative to the true discrete states. Calcium Gaussian Process Factor Analysis Gaussian Process Factor Analysis (GPFA) is a standard tool for identifying smooth continuous latent structure underlying neural population data ( Yu et al.,2009b ). Many modern variants of GPFA include count-observation likelihoods and are typically evaluated on spiking data Zhao and Park ( 2016 ); Duncker and Sahani ( 2018 ); Keeley et al. (2020b ); Wu et al. (2017 ). Here, we apply the calcium observation model to GPFA. In GPFA the prior on each latent dimension over time is a Gaussian Process. For a discrete set of time points, the prior is a multivariate Gaussian where the covariance K j is a T × T matrix whose entries are defined by a covariance function k ( t, t ′ ). Here, we use the function k ( t, t ′ ) = exp(−( t − t ′ ) 2 /(2 ℓj 2 )) with a length-scale parameter ℓj .However, the proposed approach applies to other covariance functions as well. The dimensions are concatenated at each time point to produce a vector . The neural rates are generated via a linear map from the latent space followed by a non-linearity f which generates provides the rate parameter of the Poisson distribution. The observations are then modeled with using the Calcium AR process as before, however, for this example we include an additional AR-2 comparison, where the CI dynamics evolve dependent on the previous two time-steps. where p is set to either 1 or 2 for our GPFA comparisons. The linear transformation contains, weights or “loadings” C and offsets d . Here, Ψ is a diagonal covariance matrix ( Fig. 4(b) ). The model is fit using a variational inference scheme which samples over the expectation term in the objective function (a so-called ‘black-box’ approach, see Keeley et al. (2020a )). Download figure Open in new tab Figure 4. Calcium GPFA simulated experiment using biophysical calcium imaging simulator. (a) Graphical depiction of the biophysical calcium imaging simulator. (b) The Calcium GPFA model. (c) The temporal evolution of the three true underlying latent variables and the inferred latents from the population data using different observation likelihoods (left) and the overall estimation error of the latent variables under each model (right). To evaluate the model, we used a biophysical calcium imaging simulator NAOMi ( Song et al.,2021 ) to simulate a population of 30 neurons whose spiking activity was generated from a Poisson GPFA model with three continuous latent dimensions ( Fig. 4(a) ). The simulator was critical for obtaining a realistic recording of CI data with known population ground truth spiking. Specifically, NAOMi simulates data via the biological and optical processes underlying CI, including variability in expression, optical aberrations, and calcium dynamics, providing a fair comparison that is not simply data sampled from the same statistical model used to fit the data. We fit GPFA models with CI and Gaussian observation models to the data. The model with CI observations more accurately recovers both the continuous latent states ( Fig. 4c,d ), though in this example with the NaOMI-simulated traces, we find that the order of the AR process in the likelihood plays an important role in identifying the true underlying latent structure, with the CI AR2 likelihood outperforming AR1 as well as GPFA run on the deconvolved traces. Modeling nonlinear dynamics with Calcium LFADS Latent Factor Analysis via Dynamical Systems (LFADs Sussillo et al. (2016 ); Pandarinath et al. (2018 )) is a model of neural population spiking data with nonlinear, recurrent neural network (RNN) dynamics. The model is general and provides accurate fits to neural population activity in multiple real-world settings Pandarinath et al. (2018 ); Pei et al. (2021 ); Keshtkaran et al. (2022 ). As in the previous models, the standard LFADS formulation prescribes a distribution over Poisson firing rates via a Poisson GLM readout from the RNN. For calcium imaging data, previous work has proposed to incorporate an additional set of continuous latent variables corresponding to approximate spiking ( Prince et al.,2021 ) or to fit the model to deconvolved spiking events ( Zhu et al.,2022 ). Here we proposed an alternative adaptation of LFADS to calcium imaging data by adding the calcium observation model on top of the firing rates. The schematic of the model extension can be seen in Figure 5a . Importantly, the standard amortized variational inference fitting procedure for LFADS Kingma and Welling ( 2014 ); Kingma et al. (2015 ); Sussillo et al. (2016 ) does not need to be modified with this change since the latent Poisson spike counts are marginalized out, in contrast to the model in Prince et al. (2021 ). Additionally, the model is fit directly to fluorescence measurements without deconvolution, in contrast to Zhu et al. (2022 ). Download figure Open in new tab Figure 5. (a) Calcium LFADS model. (b) Generated latents variables as well as Poisson firing rates, spiking activity, and observed calcium traces. (c) inferred latent dynamics under an (AR1) calcium model and Gaussian likelihood. (c) latent state prediction peroformance as well as inferred and true calcium trace hyperparameters. To demonstrate the approach, we synthetically generated latent time series from a 3D Lorenz attractor. We mapped the latent time series to Poisson spiking rates via a Poisson GLM. Then, we simulated Poisson spiking followed by autoregressive calcium dynamics ( Fig. 5b ). We compare the model with LFADS fit using a Gaussian observation model. In Figure 5c,d we show that the model with calcium observations infers more accurate latent variables than the model with Gaussian observations (mean latent state reconstruction R 2 =0.77 for calcium model compared to R 2 =0.37 for the Gaussian model, with 78/80 trials better reconstructed by the calcium model). Additionally, the inferred parameters of the calcium autoregressive dynamics are qualitatively similar to the ground truth parameters across the simulated neurons ( Fig. 5d ). Discussion Here we have demonstrated that a tractable likelihood for calcium imaging data can be used to adapt a variety of latent variable models in neuroscience to the setting of calcium imaging recordings. The proposed models can be fit with similar computational and inference requirements to the equivalent spiking versions of the models and do not require deconvolution methods. Importantly, our work complements a variety of recent efforts for improved estimation of neural activity from calcium imaging data. Many of these approaches work with the implicit goal of best inferring spikes from calcium traces ( Pachitariu et al.,2018 ; Wei et al.,2020 ). Although such efforts are useful, they do not take into account the uncertainty underlying spiking in estimating neural rates. Our work depends on extracting the calcium traces from fluorescence videos, which is itself an active area of research with a number of available methods ( Pnevmatikakis et al.,2016 ; Charles et al.,2022 ; Pachitariu et al.,2016 ; Dinç et al.,2021 ). As such, errors in calcium imaging source extraction, e.g., due to false transients ( Gauthier et al.,2022 ), can impact the outputs of our model in much the same way that multi-unit activity in electrophysiology. Care, therefore, should be taken to curate or validate the calcium traces, or use a robust estimation method ( Gauthier et al.,2022 ; Inan et al.,2017 ). In contrast to these spike-estimation approaches, there has been recent work that also uses generative models to describe calcium dynamics themselves, avoiding the two-step approach of first de-convolving, and then inferring latent structure from spike trains. In particular, Prince et al. (2021 ) uses a variational auto-encoder-style model to directly infer latent dynamics from raw calcium traces using latent Poisson rates and observed Gaussian likelihoods on the resultant calcium traces. Additionally, Koh et al. (2022 ) derives a generative autoregressive model of calcium dynamics from an underlying latent variables which evolve via linear dynamics. Here, there is no spiking explicitly modeled, and the underlying latent dynamics directly prescribe correlations in the observed calcium traces. Similarly, Rupasinghe et al. (2021 ) uses signal and noise correlations in the calcium population activity to extract latent time-series which generate latent (Bernoulli) spike rates, which, similar to our approach, leverage uncertainty in latent spike counts. This again avoids the two-step procedure and directly infers latent structure from population covariance. Our work here as well as in Ganmor et al. (2016 ) complements these approaches, but in contrast to them, does not enforce a specific model of latent structure. Instead, we have shown how a versatile calcium likelihood can be used in conjunction with a wide variety of latent variable models used in neuroscience. We do not assume a specific latent variable structure and instead leave that choice to the practitioner. While we demonstrate the effectiveness of this approach for GPFA, LFADS, and HMMs, there are many other existing latent variable models that could be adapted for our approach including switching dynamical systems ( Linderman et al.,2016b ; Zoltowski et al.,2020 ; Karniol-Tambour et al.,2022 ), other nonstationary dynamical systems ( Mudrik et al.,2024 , 2025 ), extensions to GPFA ( Keeley et al.,2020a ; Gokcen, 2023 ), or other auto-encoder style latent variable models ( Zhou and Wei, 2020 ; Schneider et al.,2023 ). Additionally, our proposed approach may also be relevant for large-scale models of neural recordings ( Azabou et al.,2023 ; Ye et al.,2023 ; Zhang et al.,2024 ). Such models often incorporate Poisson spiking rates and neural reconstruction losses often uses Poisson or cross-entropy observation models. The CI observation model incorporated in this paper could be applied to directly adapt such models for CI data and even presents the opportunity to jointly train large-scale models on electrophyisological and CI recordings. Overall, we demonstrate that models developed for direct use on neural spike-trains can be adapted to calcium imaging data using a simple plug-and-play approach. Our public repository 1 is available in both JAX and PyTorch implementations, with tutorials to demonstrate how they can be integrated into existing models, and it is our hope that this method will accelerate the application of powerful latent variable models to calcium imaging data. Methods Calcium observation model A likelihood originally proposed in 2016 ( Ganmor et al.,2016 ) defines a conditional probability of measured calcium fluorescence y (the so-called “dF over F”) given a Poisson firing rate λ .This model is an autoregressive (AR) model whose output depends linearly on its previous value one timestep in the past (so-called AR(1) model). This basic AR(1) version additionally has the calcium level depend on Poisson spiking with additive independent Gaussian noise on each time step: This model has three parameters: α , the AR coefficient, which determines exponential decay of fluorescence; c , the fluorescence increase due to a single spike; σ 2 , the variance of the additive Gaussian noise. Practically, this model can be interpreted as a process where an individual spike causes an instantaneous rise in calcium florescence, followed by an exponential decay due to the AR coefficient. Here, additive Gaussian noise is fed through the AR process. Importantly, the model allows us to marginalize over spike counts n , so we can consider the probability of the fluorescence given the rate, λ , and we need not consider individual spikes. We elaborate on the details of the model below. Likelihood Evaluation For model fitting and inference, we must be able to efficiently evaluate the likelihood of observed calcium responses marginalized over the unobserved spike count vector n . The independence across time bins in ( ?? ) means that these marginals can be computed independently for each time bin. That is, we can compute likelihood by summing over spike count in each bin from 0 to some maximum possible spike count R . Thus, numerical evaluation of the likelihood can be achieved practically as: where is the Poisson probability over spike counts in each time bin. In practice we will of course compute the log-likelihood, given by: Importantly, this marginalization is amenable to automatic differentiation packages, making this truly a ‘plug and play’ observation distribution for various models. Synthetic HMM Dataset and Experimental Details The synthetic HMM dataset had K =5 states and D =25 observed neurons. We generated two sequences of length T = 2000; the first was observed for training and the second was held-out for testing. The state transition matrix was designed to generate a repeating chain structure in the latent states such that the transition matrix P was The synthetic neurons were grouped into 5 groups of 5 neurons corresponding to each of the latent states. Each neuron’s firing rate was 0.2 spikes per bin. The calcium observation model parameters were selected as α ∼ 𝒩 (0.8, 0.1), c =1.0,and σ 2 =1 e −2. The AR coefficients were clipped to be within the range α ∈[0.6, 0.95]. Importantly, we added simulated measurement noise to the simulated calcium traces via a Gaussian noise model with zero mean and standard deviation of 0.2. Therefore, the simulated dataset is generated from a model that is mismatched to all models considered. We initialized the calcium AR parameters via the following procedure. The AR coefficients were initialized via a linear regression predicting the next calcium observation from the previous observation for each neuron. The initial variance was set to the squared residual error of this linear regression. The initial fluorescence increases c were set randomly from the distribution 𝒩 (1.0, 0.2). Piriform Cortex Recodings During Odor Presentation We applied the Calcium HMM model to a publicly available dataset of piriform cortex calcium imaging recordings during passive odor presentation ( Daste and Pierré, 2022 ). This dataset consists of 8 repeated presentations of 10 different odors, yielding 80 total trials. We used publicly available code to process the code into Δ F / F ( https://gitlab.com/fleischmann-lab/datasets/daste-odor-set2021-11 ). We additionally filtered out one anomalous trial and one anomalous neuron. After processing, the dataset consisted of recordings of 284 neurons across 79 trials. Each trial lasted 30 seconds and the imaging sampling frequency was 4.53 Hz (see ( Srinivasan et al.,2023 ), ‘Two-photon microscopy’ methods for additional details). We used the first 30 and last 30 trials as training trials and tested on the middle 19 trials. The observation model parameters were initialized with a two-step procedure. First, initial discrete state sequences were set using the known odor sequence or using K-means clustering. Then, we optimized the observation model parameters (both Poisson rates and calcium observation parameters) by maximizing the likelihood of the calcium observations with the discrete state sequence fixed. After this initialization, we then optimized all parameters with respective to the HMM marginal likelihood using SGD for 5000 iterations. The calcium HMM was implemented in JAX using HMM tools from the Dynamax repository 2 . This was important for computationally efficiency, as the optimization step was automatically batched across trials and compiled for speed. GPFA Synthetic Dataset with NAOMi Simulator To use the NAOMi simulator, we generated a neural volume using the anatomy module and simulated the light propagation using the optics module. Rather than use the built in Hawkes process to simulate arbitrary dynamics, we imposed the spike times from the GPFA spike generation using calcium dynamics module’s feature that enables user-defined spikes. The calcium module then generated the ground-truth fluorescence traces from the provided spikes and we generated the CI simulated videos using the scanning module. To recover single fluorescence traces from the video we used profile-assisted least-squares (PALS) as in the original NAOMi paper ( Song et al.,2021 ). In this process we scan the volume under no-noise conditions with only one neuron “on” at a time and no neuropil or other contamination. These scans provide the ground-truth spatial profiles that can be used to identify the temporal profiles by solving a per-frame least-squares optimization. To reduce noise, we regularize the time-traces with a sparsity-promoting ℓj 1 norm, also as in ( Song et al.,2021 ). ThePALStracesremove the confounding factor of cell detection, which can vary significantly between approaches and can result in bleed-through errors that can effect the inferred coding properties of the cells ( Gauthier et al.,2022 ). To generate the simulated data for Calcium GPFA inference, generate Poisson spiking activity from 30 neurons for 4000 timespoints derived from a 3 dimensional latent space, where each latent is governed by a Gaussian Process with a different temporal length scale ℓj = {250, 450, 500}. Spike-times from these 30 neurons were then used as the ground-truth spikes in the NAOMi simulator. However, after simulation only 22 calcium traces had nonzero calcium dynamics. Therefore, inference for calcium GPFA was used on 22 calcium traces. The competing models either used Gaussian observations or Poisson observations from spikes that were determined by deconvolution via SpikeML Deneux et al. (2016 ). GP length scales were all initialized to 350, calcium AR1 and AR2 parameters were initialized to.51 and {1.81, −.81} for all neurons, respectively. The noise parameter per neuron was initialized by setting it to the variance of the calcium values determined across consecutive timepoints, and the amplitude of the calcium influx due to a spike was initialized to 1 for all neurons. The model was learned via black-box variational inference for 20000 iterations. All inferred latents were regressed to the true latent before calculating the mean squared error. Calcium LFADS Model and Experimental Details To evaluate the Calcium LFADS model, we simulated synthetic calcium observations from the Lorenz dynamical system across 400 trials each of length 100 with Δ t = 0.025 using the Runge-Kutta method (RK4). We simulated observations from 30 neurons. The average number of spikes per bin was 0.42. The calcium observation model parameters for each neuron were randomly sampled with α ∼ Unif(0.8, 0.95) and c ∼ Unif(0.8, 1.2). The autoregressive dynamics noise was set to 1 e −3. After simulating the calcium traces, we added zero-mean Gaussian measurement noise with a standard deviation of 0.2 independently to each time step. Importantly, this measurement noise is not present in the generative model, as in the HMM example. The full generative model of the LFADS model with calcium observations is where F is a recurrent neural network (GRU) and the random inputs u 1: T are generated from an autoregressive process with . The model is fit as a sequential variational autoencoder as in ( Pandarinath et al.,2018 ) with encoder networks inferring approximate posterior distributions over the initial state x 0 and sequence of random inputs u 1: T . Before fitting the model, we initialized the calcium observation model hyperameters α and c . For each neuron, the AR parameter α was estimated via a linear regression predicting y t −1 from y t for timepoints where y t −1 <y t . The influx parameter c was estimated by sweeping over a range of possible values from 0.8 to 2.2 and identifying the value the best aligned the differences y t −1 − y t with quantized values 0, c , 2 c , 3 c . Acknowledgments SK was supposed by the NIH BRAIN Initiative (F32MH115445-03). DMZ was funded by the Wu Tsai Interdisciplinary Postdoctoral Research Fellowship. ASC was supported by the NSF under Grant No. 2340338 (Faculty Early Career Development Program -CAREER). JWP was supported by the Simons Collaboration on the Global Brain (SCGB AWD543027), the NIH BRAIN Initiative (9R01DA056404), and a U19 NIH-NINDS BRAIN Initiative Award (U19NS104648). Footnotes ↵ 1 https://github.com/skeeley/Calcium_likelihoods ↵ 2 https://github.com/probml/dynamax References ↵ Azabou , M. , Arora , V. , Ganesh , V. , Mao , X. , Nachimuthu , S. , Mendelson , M. , Richards , B. , Perich , M. , Lajoie , G. , and Dyer , E. ( 2023 ). A unified, scalable framework for neural population decoding . Advances in Neural Information Processing Systems , 36 : 44937 – 44956 . OpenUrl ↵ Calabrese , A. , Schumacher , J. W. , Schneider , D. M. , Paninski , L. , and Woolley , S. M. N. ( 2011 ). A generalized linear model for estimating spectrotemporal receptive fields from responses to natural sounds . PLoS One , 6 ( 1 ): e16104 . OpenUrl CrossRef PubMed ↵ Charles , A. S. , Cermak , N. , Affan , R. O. , Scott , B. B. , Schiller , J. , and Mishne , G. ( 2022 ). Graft: graph filtered temporal dictionary learning for functional neural imaging . IEEE Transactions on Image Processing , 31 : 3509 – 3524 . OpenUrl CrossRef PubMed ↵ Charles , A. S. , Park , M. , Weller , J. P. , Horwitz , G. D. , and Pillow , J. W. ( 2018 ). Dethroning the fano factor: A flexible, model-based approach to partitioning neural variability . Neural computation , 30 ( 4 ): 1012 – 1045 . OpenUrl CrossRef PubMed ↵ Chen , Y. , Mudrik , N. , Johnsen , K. A. , Alagapan , S. , Charles , A. S. , and Rozell , C. ( 2024 ). Probabilistic decomposed linear dynamical systems for robust discovery of latent neural dynamics . Advances in Neural Information Processing Systems , 37 : 104443 – 104470 . OpenUrl ↵ Cronin , B. , Stevenson , I. H. , Sur , M. , and Körding , K. P. ( 2010 ). Hierarchical bayesian modeling and markov chain monte carlo sampling for tuning-curve analysis . J Neurophysiol , 103 ( 1 ): 591 – 602 . kording. OpenUrl CrossRef PubMed Web of Science ↵ Daste , S. and Pierré , A. ( 2022 ). Two photon calcium imaging of mice piriform cortex under passive odor presentation. (Version 0.220928.1306) [Data set] . DANDI archive . doi: 10.48324/dandi.000167/0.220928.1306 . OpenUrl CrossRef ↵ Demas , J. , Manley , J. , Tejera , F. , Barber , K. , Kim , H. , Traub , F. M. , Chen , B. , and Vaziri , A. ( 2021 ). High-speed, cortex-wide volumetric recording of neuroactivity at cellular resolution using light beads microscopy . Nature Methods , 18 ( 9 ): 1103 – 1111 . OpenUrl PubMed ↵ Dempster , A. , Laird , N. , and Rubin , R. ( 1977 ). Maximum likelihood from incomplete data via the EM algorithm . J. Royal Statistical Society, B , 39 ( 1 ): 1 – 38 . OpenUrl ↵ Deneux , T. , Kaszas , A. , Szalay , G. , Katona , G. , Lakner , T. , Grinvald , A. , Rózsa , B. , and Vanzetta , I. ( 2016 ). Accurate spike estimation from noisy calcium signals for ultrafast three-dimensional imaging of large neuronal populations in vivo . Nature communications , 7 ( 1 ): 12190 . OpenUrl PubMed ↵ Dinç , F. , Inan , H. , Hernandez , O. , Schmuckermair , C. , Hazon , O. , Tasci , T. , Ahanonu , B. O. , Zhang , Y. , Lecoq , J. , Haziza , S. , et al. ( 2021 ). Fast, scalable, and statistically robust cell extraction from large-scale neural calcium imaging datasets . BioRxiv , pages 2021 – 03 . ↵ Duncker , L. , Bohner , G. , Boussard , J. , and Sahani , M. ( 2019 ). Learning interpretable continuous-time models of latent stochastic dynamical systems . arXiv preprint arXiv: 1902.04420 . Duncker , L. and Sahani , M. ( 2018 ). Temporal alignment and latent gaussian process factor inference in population spike trains . In Advances in Neural Information Processing Systems , pages 10445 – 10455 . ↵ Escola , S. , Fontanini , A. , Katz , D. , and Paninski , L. ( 2011 ). Hidden markov models for the stimulus-response relationships of multistate neural systems . Neural Computation , 23 ( 5 ): 1071 – 1132 . OpenUrl CrossRef PubMed ↵ Ganmor , E. , Krumin , M. , Rossi , L. F. , Carandini , M. , and Simoncelli , E. P. ( 2016 ). Direct estimation of firing rates from calcium imaging data . arXiv preprint arXiv: 1601.00364 . ↵ Gao , Y. , Busing , L. , Shenoy , K. V. , and Cunningham , J. P. ( 2015 ). High-dimensional neural spike train analysis with generalized count linear dynamical systems . In Advances in neural information processing systems ,pages 2044 – 2052 . ↵ Gauthier , J. L. , Koay , S. A. , Nieh , E. H. , Tank , D. W. , Pillow , J. W. , and Charles , A. S. ( 2022 ). Detecting and correcting false transients in calcium imaging . Nature Methods , 19 ( 4 ): 470 – 478 . OpenUrl PubMed ↵ Gokcen , E. ( 2023 ). Disentangling communication across populations of neurons . Carnegie Mellon University . ↵ Goris , R. L. T. , Movshon , J. A. , and Simoncelli , E. P. ( 2014 ). Partitioning neuronal variability . Nat Neurosci , 17 ( 6 ): 858 – 865 . OpenUrl CrossRef PubMed ↵ Helmchen , F. and Tank , D. W. ( 2015 ). A single-compartment model of calcium dynamics in nerve terminals and dendrites . Cold Spring Harbor Protocols , 2015 ( 2 ): pdb – top085910 . OpenUrl ↵ Huang , L. , Ledochowitsch , P. , Knoblich , U. , Lecoq , J. , Murphy , G. J. , Reid , R. C. , de Vries , S. E. , Koch , C. , Zeng , H. , Buice , M. A. , et al. ( 2021 ). Relationship between simultaneously recorded spiking activity and fluorescence signal in gcamp6 transgenic mice . Elife , 10 : e51675 . OpenUrl CrossRef PubMed ↵ Inan , H. , Erdogdu , M. A. , and Schnitzer , M. ( 2017 ). Robust estimation of neural signals in calcium imaging . Advances in neural information processing systems , 30 . ↵ Karniol-Tambour , O. , Zoltowski , D. M. , Diamanti , E. M. , Pinto , L. , Tank , D. W. , Brody , C. D. , and Pillow , J. W. ( 2022 ). Modeling communication and switching nonlinear dynamics in multi-region neural activity . bioRxiv ,pages 2022 – 09 . ↵ Keeley , S. , Aoi , M. , Yu , Y. , Smith , S. , and Pillow , J. W. ( 2020a ). Identifying signal and noise structure in neural population activity with gaussian process factor models . Advances in Neural Information Processing Systems , 33 : 13795 – 13805 . OpenUrl ↵ Keeley , S. , Zoltowski , D. , Yu , Y. , Smith , S. , and Pillow , J. ( 2020b ). Efficient non-conjugate gaussian process factor models for spike count data using polynomial approximations . In International Conference on Machine Learning , pages 5177 – 5186 . PMLR . ↵ Keeley , S. L. , Zoltowski , D. M. , Yu , Y. , Yates , J. L. , Smith , S. L. , and Pillow , J. W. ( 2019 ). Efficient non-conjugate gaussian process factor models for spike count data using polynomial approximations . arXiv preprint arXiv: 1906.03318 . ↵ Keshtkaran , M. R. , Sedler , A. R. , Chowdhury , R. H. , Tandon , R. , Basrai , D. , Nguyen , S. L. , Sohn , H. , Jazayeri , M. , Miller , L. E. , and Pandarinath , C. ( 2022 ). A large-scale neural network training framework for generalized estimation of single-trial population dynamics . Nature Methods , 19 ( 12 ): 1572 – 1577 . OpenUrl PubMed ↵ Kingma , D. P. , Salimans , T. , and Welling , M. ( 2015 ). Variational dropout and the local reparameterization trick . In Advances in Neural Information Processing Systems , pages 2575 – 2583 . Kingma , D. P. and Welling , M. ( 2014 ). Auto-encoding variational bayes . arXiv : arxiv: 1312.6114 . ↵ Ko , H. , Cossell , L. , Baragli , C. , Antolik , J. , Clopath , C. , Hofer , S. B. , and Mrsic-Flogel , T. D. ( 2013 ). The emergence of functional microcircuits in visual cortex . Nature , 496 ( 7443 ): 96 – 100 . OpenUrl CrossRef PubMed Web of Science ↵ Koh , T. H. , Bishop , W. E. , Kawashima , T. , Jeon , B. B. , Srinivasan , R. , Kuhlman , S. J. , Ahrens , M. B. , Chase , S. M. , and Byron , M. Y. ( 2022 ). Dimensionality reduction of calcium-imaged neuronal population activity . bioRxiv . ↵ Krause , E. L. and Drugowitsch , J. ( 2022 ). A large majority of awake hippocampal sharp-wave ripples feature spatial trajectories with momentum . Neuron , 110 ( 4 ): 722 – 733 . OpenUrl PubMed ↵ Ledochowitsch , P. , Huang , L. , Knoblich , U. , Oliver , M. , Lecoq , J. , Reid , C. , Li , L. , Zeng , H. , Koch , C. , Waters , J. , et al. ( 2019 ). On the correspondence of electrical and optical physiology in in vivo population-scale two-photon calcium imaging . BioRxiv , page 800102 . ↵ Linderman , S. , Adams , R. P. , and Pillow , J. W. ( 2016a ). Bayesian latent structure discovery from multi-neuron recordings . In Advances in neural information processing systems , pages 2002 – 2010 . ↵ Singh , A. and Zhu , J. Linderman , S. , Johnson , M. , Miller , A. , Adams , R. , Blei , D. , and Paninski , L. ( 2017 ). Bayesian Learning and Inference in Recurrent Switching Linear Dynamical Systems . In Singh , A. and Zhu , J. , editors, Proceedings of the 20th International Conference on Artificial Intelligence and Statistics, volume 54 of Proceedings of Machine Learning Research , pages 914 – 922 , Fort Lauderdale, FL, USA . PMLR . ↵ Linderman , S. , Miller , A. C. , Adams , R. P. , Blei , D. M. , Paninski , L. , and Johnson , M. J. ( 2016b ). Recurrent switching linear dynamical systems . Arxiv preprints . ↵ Macke , J. H. , Buesing , L. , Cunningham , J. P. , Byron , M. Y. , Shenoy , K. V. , and Sahani , M. ( 2011 ). Empirical models of spiking in neural populations . In Advances in neural information processing systems , volume 24 , pages 1350 – 1358 . OpenUrl ↵ McFarland , J. M. , Cui , Y. , and Butts , D. A. ( 2013 ). Inferring nonlinear neuronal computation based on physiologically plausible inputs . PLoS Comput Biol , 9 ( 7 ): e1003143. #x002B;. OpenUrl CrossRef PubMed ↵ Mudrik , N. , Chen , Y. , Yezerets , E. , Rozell , C. J. , and Charles , A. S. ( 2024 ). Decomposed linear dynamical systems (dlds) for learning the latent components of neural dynamics . Journal of Machine Learning Research , 25 ( 59 ): 1 – 44 . OpenUrl CrossRef ↵ Mudrik , N. , Ly , R. , Ruebel , O. , and Charles , A. S. ( 2025 ). Creimbo: Cross-regional ensemble interactions in multi-view brain observations . In The Thirteenth International Conference on Learning Representations . ↵ Pachitariu , M. , Stringer , C. , and Harris , K. D. ( 2018 ). Robustness of spike deconvolution for neuronal calcium imaging . Journal of Neuroscience , 38 ( 37 ): 7976 – 7985 . OpenUrl Abstract / FREE Full Text ↵ Pachitariu , M. , Stringer , C. , Schröder , S. , Dipoppa , M. , Rossi , L. F. , Carandini , M. , and Harris , K. D. ( 2016 ). Suite2p: beyond 10,000 neurons with standard two-photon microscopy . BioRxiv , page 061507 . ↵ Pandarinath , C. , O’Shea , D. J. , Collins , J. , Jozefowicz , R. , Stavisky , S. D. , Kao , J. C. , Trautmann , E. M. , Kaufman , M. T. , Ryu , S. I. , Hochberg , L. R. , Henderson , J. M. , Shenoy , K. V. , Abbott , L. F. , and Sussillo , D. ( 2018 ). Inferring single-trial neural population dynamics using sequential auto-encoders . Nature Methods , 15 ( 10 ): 805 – 815 . OpenUrl PubMed ↵ Burges , C. J. C. , Bottou , L. , Welling , M. , Ghahramani , Z. , and Weinberger , K. Q. Park , I. M. , Archer , E. W. , Priebe , N. , and Pillow , J. W. ( 2013 ). Spectral methods for neural characterization using generalized quadratic models . In Burges , C. J. C. , Bottou , L. , Welling , M. , Ghahramani , Z. , and Weinberger , K. Q. , editors, Advances in Neural Information Processing Systems 26 , pages 2454 – 2462 . Curran Associates, Inc . OpenUrl ↵ Park , M. , Weller , J. P. , Horwitz , G. D. , and Pillow , J. W. ( 2014 ). Bayesian active learning of neural firing rate maps with transformed gaussian process priors . Neural Computation , 26 ( 8 ): 1519 – 1541 . OpenUrl CrossRef PubMed ↵ Pei , F. , Ye , J. , Zoltowski , D. , Wu , A. , Chowdhury , R. H. , Sohn , H. , O’Doherty , J. E. , Shenoy , K. V. , Kaufman , M. T. , Churchland , M. , et al. ( 2021 ). Neural latents benchmark’21: evaluating latent variable models of neural population activity . arXiv preprint arXiv: 2109.04463 . ↵ Bartlett , P. , Pereira , F. , Burges , C. , Bottou , L. , and Weinberger , K. Pillow , J. and Scott , J. ( 2012 ). Fully bayesian inference for neural models with negative-binomial spiking . In Bartlett , P. , Pereira , F. , Burges , C. , Bottou , L. , and Weinberger , K. , editors, Advances in Neural Information Processing Systems 25 , pages 1907 – 1915 . OpenUrl ↵ Pillow , J. W. , Shlens , J. , Paninski , L. , Sher , A. , Litke , A. M. , and Chichilnisky , E. J. Simoncelli , E. P. ( 2008 ). Spatiotemporal correlations and visual signaling in a complete neuronal population . Nature , 454 : 995 – 999 . OpenUrl CrossRef PubMed Web of Science ↵ Pnevmatikakis , E. A. , Soudry , D. , Gao , Y. , Machado , T. A. , Merel , J. , Pfau , D. , Reardon , T. , Mu , Y. , Lacefield , C. , Yang , W. , et al. ( 2016 ). Simultaneous denoising, deconvolution, and demixing of calcium imaging data . Neuron , 89 ( 2 ): 285 – 299 . OpenUrl CrossRef PubMed ↵ Prince , L. Y. , Bakhtiari , S. , Gillon , C. J. , and Richards , B. A. ( 2021 ). Parallel inference of hierarchical latent dynamics in two-photon calcium imaging of neuronal populations . bioRxiv . ↵ Rad , K. R. , Machado , T. A. , Paninski , L. , et al. ( 2017 ). Robust and scalable bayesian analysis of spatial neural tuning function data . The Annals of Applied Statistics , 11 ( 2 ): 598 – 637 . OpenUrl ↵ Shawe-Taylor , J. , Zemel , R. , Bartlett , P. , Pereira , F. , and Weinberger , K. Rad , K. R. and Paninski , L. ( 2011 ). Information rates and optimal decoding in large neural populations . In Shawe-Taylor , J. , Zemel , R. , Bartlett , P. , Pereira , F. , and Weinberger , K. , editors, Advances in Neural Information Processing Systems 24 , pages 846 – 854 . OpenUrl ↵ Rupasinghe , A. , Francis , N. , Liu , J. , Bowen , Z. , Kanold , P. O. , and Babadi , B. ( 2021 ). Direct extraction of signal and noise correlations from two-photon calcium imaging of ensemble neuronal activity . Elife , 10 : e68046 . OpenUrl CrossRef PubMed ↵ Lee , D. D. , Sugiyama , M. , Luxburg , U. V. , Guyon , I. , and Garnett , R. Savin , C. and Tkacik , G. ( 2016 ). Estimating nonlinear neural response functions using gp priors and kronecker methods . In Lee , D. D. , Sugiyama , M. , Luxburg , U. V. , Guyon , I. , and Garnett , R. , editors, Advances in Neural Information Processing Systems 29 , pages 3603 – 3611 . Curran Associates, Inc . OpenUrl ↵ Schneider , S. , Lee , J. H. , and Mathis , M. W. ( 2023 ). Learnable latent embeddings for joint behavioural and neural analysis . Nature , 617 ( 7960 ): 360 – 368 . OpenUrl CrossRef PubMed ↵ Smith , A. C. and Brown , E. N. ( 2003 ). Estimating a state-space model from point process observations . Neural computation , 15 ( 5 ): 965 – 991 . OpenUrl CrossRef PubMed Web of Science ↵ Smith , S. L. and Häusser , M. ( 2010 ). Parallel processing of visual space by neighboring neurons in mouse visual cortex . Nature neuroscience , 13 ( 9 ): 1144 – 1149 . OpenUrl CrossRef PubMed ↵ Song , A. , Gauthier , J. L. , Pillow , J. W. , Tank , D. W. , and Charles , A. S. ( 2021 ). Neural anatomy and optical microscopy (naomi) simulation for evaluating calcium imaging methods . Journal of Neuroscience Methods , 358 : 109173 . OpenUrl CrossRef PubMed ↵ Srinivasan , S. , Daste , S. , Modi , M. N. , Turner , G. C. , Fleischmann , A. , and Navlakha , S. ( 2023 ). Effects of stochastic coding on olfactory discrimination in flies and mice . PLoS Biology , 21 ( 10 ): e3002206 . OpenUrl CrossRef PubMed ↵ Stevenson , I. H. ( 2016 ). Flexible models for spike count data with both over-and under-dispersion . Journal of computational neuroscience , 41 ( 1 ): 29 – 43 . OpenUrl CrossRef PubMed ↵ Sussillo , D. , Jozefowicz , R. , Abbott , L. , and Pandarinath , C. ( 2016 ). Lfads-latent factor analysis via dynamical systems . arXiv preprint arXiv: 1608.06315 . ↵ Truccolo , W. , Eden , U. T. , Fellows , M. R. , Donoghue , J. P. , and Brown , E. N. ( 2005 ). A point process framework for relating neural spiking activity to spiking history, neural ensemble and extrinsic covariate effects . J. Neurophysiol , 93 ( 2 ): 1074 – 1089 . OpenUrl CrossRef PubMed Web of Science ↵ Wei , X.-X. , Zhou , D. , Grosmark , A. , Ajabi , Z. , Sparks , F. , Zhou , P. , Brandon , M. , Losonczy , A. , and Paninski , L. ( 2020 ). A zero-inflated gamma model for deconvolved calcium imaging traces . arXiv preprint arXiv: 2006.03737 . ↵ Williamson , R. S. , Sahani , M. , and Pillow , J. W. ( 2015 ). The equivalence of information-theoretic and likelihood-based methods for neural dimensionality reduction . PLoS Comput Biol , 11 ( 4 ): e1004141 . OpenUrl CrossRef PubMed ↵ Guyon , I. , Luxburg , U. V. , Bengio , S. , Wallach , H. , Fergus , R. , Vishwanathan , S. , and Garnett , R. Wu , A. , Roy , N. G. , Keeley , S. , and Pillow , J. W. ( 2017 ). Gaussian process based nonlinear latent structure discovery in multivariate spike train data . In Guyon , I. , Luxburg , U. V. , Bengio , S. , Wallach , H. , Fergus , R. , Vishwanathan , S. , and Garnett , R. , editors, Advances in Neural Information Processing Systems 30 , pages 3496 – 3505 . Curran Associates, Inc . OpenUrl PubMed ↵ Ye , J. , Collinger , J. , Wehbe , L. , and Gaunt , R. ( 2023 ). Neural data transformer 2: multi-context pretraining for neural spiking activity . Advances in Neural Information Processing Systems , 36 : 80352 – 80374 . OpenUrl ↵ Yezerets , E. , Mudrik , N. , and Charles , A. S. ( 2025 ). Decomposed linear dynamical systems (dlds) models reveal instantaneous, context-dependent dynamic connectivity in c . elegans. Communications Biology , 8 ( 1 ): 1218 . OpenUrl PubMed ↵ Yu , B. M. , Cunningham , J. P. , Santhanam , G. , Ryu , S. I. , Shenoy , K. V. , and Sahani , M. ( 2009a ). Gaussian-process factor analysis for low-dimensional single-trial analysis of neural population activity . Journal of Neurophysiology , 102 ( 1 ): 614 . OpenUrl CrossRef PubMed Web of Science ↵ Yu , B. M. , Cunningham , J. P. , Santhanam , G. , Ryu , S. I. , Shenoy , K. V. , and Sahani , M. ( 2009b ). Gaussian-process factor analysis for low-dimensional single-trial analysis of neural population activity . In Advances in Neural Information Processing Systems , pages 1881 – 1888 . ↵ Zemel , R. S. , Dayan , P. , and Pouget , A. ( 1998 ). Probabilistic interpretation of population codes . Neural Comput , 10 ( 2 ): 403 – 430 . OpenUrl CrossRef PubMed Web of Science ↵ Zhang , K. , Ginzburg , I. , McNaughton , B. , and Sejnowski , T. ( 1998 ). Interpreting neuronal population activity by reconstruction: Unified framework with application to hippocampal place cells . Journal of Neurophysiology , 79 : 1017 – 1044 . OpenUrl CrossRef PubMed Web of Science ↵ Zhang , Y. , Wang , Y. , Jiménez-Benetó , D. , Wang , Z. , Azabou , M. , Richards , B. , Tung , R. , Winter , O. , Dyer , E. , Paninski , L. , et al. ( 2024 ). Towards a” universal translator” for neural dynamics at single-cell, single-spike resolution . Advances in Neural Information Processing Systems , 37 : 80495 – 80521 . OpenUrl PubMed Zhao , Y. and Park , I. M. ( 2016 ). Variational latent gaussian process for recovering single-trial dynamics from population spike trains . arXiv preprint arXiv: 1604.03053 . ↵ Zhao , Y. and Park , I. M. ( 2017 ). Variational latent gaussian process for recovering single-trial dynamics from population spike trains . Neural Computation , 29 ( 5 ): 1293 – 1316 . GP prior over latent variable, LNP observation model. OpenUrl CrossRef PubMed ↵ Zhao , Y. and Park , I. M. ( 2019 ). Variational joint filtering . arXiv preprint arXiv: 1707.09049 . Memming’s latest on nonlinear latent dynamics w/Poisson noise. ↵ Zhou , D. and Wei , X.-X. ( 2020 ). Learning identifiable and interpretable latent models of high-dimensional neural activity using pi-vae . Advances in Neural Information Processing Systems , 33 : 7234 – 7247 . OpenUrl ↵ Zhu , F. , Grier , H. A. , Tandon , R. , Cai , C. , Agarwal , A. , Giovannucci , A. , Kaufman , M. T. , and 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 . OpenUrl CrossRef PubMed ↵ Zoltowski , D. , Pillow , J. , and Linderman , S. ( 2020 ). A general recurrent state space framework for modeling neural dynamics during decision-making . In International Conference on Machine Learning , pages 11680 – 11691 . PMLR . ↵ Bengio , S. , Wallach , H. , Larochelle , H. , Grauman , K. , Cesa-Bianchi , N. , and Garnett , R. Zoltowski , D. and Pillow , J. W. ( 2018 ). Scaling the poisson glm to massive neural datasets through polynomial approximations . In Bengio , S. , Wallach , H. , Larochelle , H. , Grauman , K. , Cesa-Bianchi , N. , and Garnett , R. , editors, Advances in Neural Information Processing Systems 31 , pages 3521 – 3531 . Curran Associates, Inc . OpenUrl View the discussion thread. Back to top Previous Next Posted October 17, 2025. Download PDF Email Thank you for your interest in spreading the word about bioRxiv. NOTE: Your email address is requested solely to identify you as the sender of this article. Your Email * Your Name * Send To * Enter multiple addresses on separate lines or separate them with commas. You are going to email the following Improved inference of latent neural states from calcium imaging data Message Subject (Your Name) has forwarded a page to you from bioRxiv Message Body (Your Name) thought you would like to see this page from the bioRxiv website. Your Personal Message CAPTCHA This question is for testing whether or not you are a human visitor and to prevent automated spam submissions. Share Improved inference of latent neural states from calcium imaging data Stephen Keeley , David Zoltowski , Adam Charles , Jonathan Pillow bioRxiv 2025.10.17.682993; doi: https://doi.org/10.1101/2025.10.17.682993 Share This Article: Copy Citation Tools Improved inference of latent neural states from calcium imaging data Stephen Keeley , David Zoltowski , Adam Charles , Jonathan Pillow bioRxiv 2025.10.17.682993; doi: https://doi.org/10.1101/2025.10.17.682993 Citation Manager Formats BibTeX Bookends EasyBib EndNote (tagged) EndNote 8 (xml) Medlars Mendeley Papers RefWorks Tagged Ref Manager RIS Zotero Tweet Widget Facebook Like Google Plus One Subject Area Neuroscience Subject Areas All Articles Animal Behavior and Cognition (7635) Biochemistry (17691) Bioengineering (13892) Bioinformatics (41936) Biophysics (21452) Cancer Biology (18588) Cell Biology (25504) Clinical Trials (138) Developmental Biology (13378) Ecology (19899) Epidemiology (2067) Evolutionary Biology (24320) Genetics (15609) Genomics (22506) Immunology (17736) Microbiology (40394) Molecular Biology (17181) Neuroscience (88605) Paleontology (666) Pathology (2832) Pharmacology and Toxicology (4824) Physiology (7641) Plant Biology (15153) Scientific Communication and Education (2045) Synthetic Biology (4294) Systems Biology (9825) Zoology (2271)