Code accompanying the paper "Dimensionality reduction of calcium-imaged neuronal population activity"[1]. Runs on MATLAB2019a. Note that this code is currently still in preliminary form.
CILDS generative model, see methods for variable definitions and dimensions.
CILDS allows estimates of shared activity among neurons (i.e., the latent variables) to influence the estimates of deconvolved spiking activity, and vice versa. In other words CILDS performs deconvolution for all neurons and dimensionality reduction jointly, in a unified framework. This is in contrast to deconv-LDS, which deconvolves the activity of each neuron independently. With low-dimensional latent variables that are jointly estimated with the model of calcium decay, CILDS is better able to peer through the calcium decay to more clearly identify the shared activity among neurons, as compared to deconv-LDS and LDS applied directly on fluorescence. The differences between models are illustrated in the figure below.
add CILDS and necessary folders to the search path of MATLAB.
File location: main
>> cilds_setup
If using deconvolution or 'ldsInit' for CILDS, you must download OASIS[2] from https://github.com/zhoupc/OASIS_matlab. Add OASIS to the search path of MATLAB. The demo does not need this.
File location: https://github.com/zhoupc/OASIS_matlab
>> oasis_setup
File location: core_cilds
>> script_cilds_demo
An example of CILDS being used is
[EstParam, Result] = cilds(data, RunParam);-
data - Structure containing recorded data (fluorescence traces)
- Dimensions: 1 x N_TRIAL
- Fields:
- y (N_NEURON x T) -- neural data
-
RunParam - Structure containing dimension of observations (no. neurons) and dimension of latent variables, training indices and testing indices
- Dimensions: 1 x 1
- Fields:
- N_LATENT -- desired latent dimension
- TRAININD (only if splitting training and testing) -- trial indices of training data
- TESTIND (only if splitting training and testing) -- trial indices of testing data
Note that to use 'ldsInit' for initializing CILDS, OASIS must be downloaded.
In the given demo, a few tests are given for running CILDS. For instance, running test 5 checks that without the true parameters, the posteriors approach the ground truth latent variables with enough EM iterations (given a reasonable SNR)
-
EstParam - Structure containing estimated model parameters from maximization step
- Dimensions: 1 x 1
- Fields:
- A (N_NEURON x N_LATENT),
- B (N_NEURON x N_NEURON),
- G (N_NEURON x N_NEURON) - Gamma,
- D (N_LATENT x N_LATENT),
- Q (N_NEURON x N_NEURON),
- R (N_NEURON x N_NEURON),
- P (N_NEURON x N_LATENT),
- b (N_NEURON x 1),
- mu_1 (N_NEURON x 1),
- cov_1 (N_NEURON x N_NEURON)
- h_2 (N_LATENT x 1),
- G_2 (N_LATENT x N_LATENT)
-
Result - Structure containing estimated posteriors (latent variables) from expectation step
- Dimensions - 1 x N_TRIAL
- Fields (If leaveoneout toggled, flProj and frProj saved (predicted fluorescence and firing rate):
- z (N_LATENT x T) - latent variables
- c (N_NEURON x T) - calcium
Sample figure produced, in a manner similar to Fig. 2 of paper (also shown in this document), right panel with estimated latent variables plotted over ground truth latent variables over time:
Additional expected outputs are included in the figures folder under the header cilds_check*.
File location: main
>> script_simdata
Run CILDS/CIFA/LDS/deconv-LDS on simulated data (note that this won't run if you haven't run script_simdata and if you haven't downloaded OASIS and added it to the path)
If you do want to run it without installing OASIS, add 'initType','randInit' to cilds_crossvalidate.
File location: main
>> script_simdimred
This performs a 2-fold cross-validation using the chosen dimensionality reduction method. For example, if using CILDS, the function called would be
cilds_crossvalidate(data,RunParam(iData),'zDimList',RunParam(iData).N_LATENT,...
'numFold',2,'maxIter',maxIter,'fileHeader',resultFile,...
'InitParam',InitParam);The example performs 2-fold cross-validation using N_LATENT number of latent dimensions, a maximum number of EM iterations as specified, saving the results in resultFile and initializing certain parameters as specified in InitParam.
-
data - Structure containing recorded data (fluorescence traces)
- Dimensions: 1 x N_TRIAL
- Fields:
- y (N_NEURON x T) -- neural data
-
RunParam - Structure containing dimension of observations (no. neurons) and dimension of latent variables, training indices and testing indices
- Dimensions: 1 x 1
- Fields:
- N_LATENT (not necessary if 'zDimList' is specified) -- desired latent dimension
- TRAININD (necessary in cross-validation) -- trial indices of training data
- TESTIND (necessary in cross-validation)-- trial indices of testing data
-
'numFold' - Scalar indicating number of crossvalidation folds
-
'zDimList - Scalar or vector containing latent dimensions to be used in dimensionality reduction
-
'maxIter' - Scalar indicating number of expectation-maximization iterations to run (default is 500)
-
'InitParam' - Structure containing user-defined initialization parameters
- Dimensions: 1 x 1
- Possible fields: A,B,G,D,Q,R,P,b,mu_1,cov_1,h_2,G_2(see above EstParam description for dimensions)
-
'fileHeader' - String containing start of file name for saving (default: nan)
-
'initType' - String to choose initialization type. (default (won't run if OASIS is not installed):'ldsInit').
Compute R2 between estimated latent variables and ground truth latent variables in given sample simulated data
Note that the process of data generation and dimensionality reduction takes awhile so the results from the run are provided in sim_stat just for plotting purposes. This step produces something like the figure below. File location: main
>> script_simstat
Matlab (2019a) using Intel(R) Xeon(R) CPU processors (Gold 6230, 2.1 GHz) with 250 GB of RAM.
[1] Koh, T.H., Bishop, W.E., Kawashima, T. et al. Dimensionality reduction of calcium-imaged neuronal population activity. Nat Comput Sci 3, 71–85 (2023). https://doi.org/10.1038/s43588-022-00390-2
[2] Vogelstein, J. T. et al. Fast nonnegative deconvolution for spike train inference from population calcium imaging. J. Neurophysiol. 104, 3691–3704, 10.1152/jn.01073.2009 (2010). PMID: 20554834, https://doi.org/10.1152/jn.01073.2009.