LZeroSpikeInference: A package for estimating spike times from calcium imaging data using an L0 penalty ![Build Status](https://camo.githubusercontent.com/1d46a598efa346c1d6348f662db0b38615575665ce2df4a166b3c62c8867fb0f/68747470733a2f2f7472617669732d63692e636f6d2f6a6577656c6c7365616e2f4c5a65726f5370696b65496e666572656e63652e7376673f746f6b656e3d6f697856667462727132546b72534170524b6e32266272616e63683d6d6173746572)
This package implements an algorithm for deconvolving calcium imaging data for a single neuron in order to estimate the times at which the neuron spikes. See https://jewellsean.github.io/fast-spike-deconvolution/ for tutorials and additional information.
This algorithm solves the optimization problems
minimize_{c1,...,cT} 0.5 sum_{t=1}^T ( y_t - c_t )^2 + lambda sum_{t=2}^T 1_{c_t neq gamma c_{t-1} }
for the global optimum, where y_t is the observed fluorescence at the tth timepoint. We also solve the above problem with the constraint that c_t >= 0 (hardThreshold = T).
minimize_{c1,...,cT,b1,...,bT} 0.5 sum_{t=1}^T (y_t - c_t - b_t)^2 + lambda sum_{t=2}^T 1_{c_t neq gamma c_{t-1}, b_t neq b_{t-1} }
where the indicator variable 1_{(A,B)} equals 1 if the event A cup B holds, and equals zero otherwise.
In R, if devtools
is installed type
devtools::install_github("jewellsean/LZeroSpikeInference")
Once installed type
library(LZeroSpikeInference)
?LZeroSpikeInference
This package can be called from Python using the py2 package. To install LZeroSpikeInference and rpy2 for use in Python first
- Install R (for example
apt-get install r-base
)
and then from within R install this package (as above). Then pip install rpy2
- pip install --user rpy2
The following example illustrates use of the LZeroSpikeInference package from python
from numpy import array
import rpy2.robjects.packages
lzsi = rpy2.robjects.packages.importr("LZeroSpikeInference")
d = lzsi.simulateAR1(n = 500, gam = 0.998, poisMean = 0.009, sd = 0.15, seed = 8)
fit = lzsi.estimateSpikes(d[1], **{'gam':0.998, 'lambda':8, 'type':"ar1"})
spikes = array(fit[0])
fittedValues = array(fit[1])
Thanks to Luke Campagnola for suggesting this approach!
See Jewell and Witten, Exact Spike Train Inference Via L0 Optimization (2017)