Nonlinear super-resolution signal processing allows intracellular tracking of calcium dynamics

Fluorescence imaging of calcium dynamics has revolutionized cellular biology, and especially neuroscience as it allows the study of neural activity across spatially well defined populations. Quantification of fluorescence signals is commonly performed using ratiometric measures. These measures, such as ΔF/F, are easy to implement, but they depend on the definition of a baseline, which is often not trivial outside of the measurement of evoked responses, and needs to be defined for every cell (or other cellular compartment) separately. Here, we present a new quantitative measure of fluorescence by taking advantage of the time dimension of the signal. The new method, which we named ARES (Autoregressive RESiduals), is based on the quantification of residuals after linear autoregression and does not require arbitrary baseline assignment. We detail the basic characteristics of the ARES signal, compare it against ΔF/F, and quantify the improvement in spatial and temporal resolution of recorded calcium data. We further exemplify its utility in the study of intracellular calcium dynamics by describing the propagation of a calcium wave inside a dendrite. As a robust and accurate method for quantifying fluorescence signals, ARES is particularly well-suited for the study of spontaneous as well as evoked calcium dynamics, subcellular localization of calcium signals, and spatiotemporal tracking of calcium dynamics. Highlights A baseline-free quantitative measure of fluorescence is introduced The method, called ARES, is based on autoregression ARES improves the spatiotemporal resolution of calcium imaging recordings Its utility for the localization and spatiotemporal tracking of calcium is shown


Introduction
Neuronal activity is traditionally studied using either electrophysiological methods or fluorescence, most typically, calcium imaging.Electrical recording techniques provide the best temporal resolution with sampling rates of > 20 kHz.However, they rely on inference methods (e.g.spike sorting) to determine cell identities.[1,2].Calcium imaging, on the other hand, follows the slower time course of calcium dynamics, at the point of having to rely on an independent set of inference methods to attempt to reconstruct the underlying spiking activity of neurons [3].However, thanks to the direct photonic observation of networks, calcium imaging enables spatial reconstruction of functional networks [2,4], including in subcellular resolution [5].
Notwithstanding its importance, calcium imaging has several limitations.First, the nonlinear and slow dynamics of cellular calcium concentrations and indicators hamper the efforts to reconstruct the spiking activity of neurons from calcium traces [3,6,7].Second, chromophores undergo "photobleaching", where the structure of the molecules themselves is modified during the fluorescence process, making them permanently lose their light emission properties [8].Photobleaching is particularly prominent during prolonged recording sessions, which is common in functional studies and live recordings of neurons [9].Light scattering also contributes to the decrease of the spatial signal-to-noise ratio and is often difficult to control [2].
The current standard analysis methods for calcium imaging recordings are ratiometric measures, methods that are based on computing the difference between the signal of interest and a baseline [10], or between two different wavelengths as applied during Fluorescence Resonance Energy Transfer (FRET) imaging [11,12].∆F/F, which is the difference in fluorescence between any frame and a reference (∆F) in a region of interest (ROI), normalized by the reference value (/F), is by far the most common analysis method used for functional calcium imaging.The simplicity and speed of implementation, and the ease of interpretation of this ratiometric method ensured its widespread implementation.However, ∆F/F has several limita-tions.A ratiometric measure of this kind has the most efficient application when the baseline remains largely static, otherwise the determination of a correct baseline is non-trivial, requiring for instance kernel density estimates [13], and might greatly affect the outcome of a study.∆F/F also fails to correct for photobleaching effects in case of long time intervals between the current signal and the baseline value (i.e.prolonged recordings), or unless the baseline is continuously estimated, which is often sub-optimal.Furthermore, ∆F/F does not take advantage of the time dimension present in many applications of the calcium imagine technique which limits its ability to reconstruct the non-linear calcium signals.In an attempt to overcome some of these shortcomings, and to take advantage of the time dimension present in prolonged recordings, we developed ARES (Autoregression RESiduals; https://github.com/DepartmentofNeurophysiology/ARES).ARES is based on the analysis of a sliding time window of fixed length, to which we apply a pixel-wise linear autoregression (predicting a pixel value on the basis of its past values).The time series of the residuals obtained by the comparison between the predicted signal and the measured signal is averaged over a time window to reduce the noise and isolate the transient (novel) events in the signal.Application of ARES to the analysis of calcium dynamics during single-cell stimulation showed that ARES improves the spatial and temporal resolution in functional signals with respect to the raw signal and standard ∆F/F methods.Because ARES provides subcellular information about calcium dynamics, it opens new avenues in high-resolution imaging using calcium indicators.

Experimental methods
All animal procedures were approved by the Radboud University Animal Experiment Committee.Cultures of dissociated cortical neurons were prepared from postnatal day 0 mice on C57BL/6J background.The mother was sacrificed by cervical dislocation, pups by decapitation before their brains were rapidly removed.

Calcium Imaging.
Neurons were visualized using LED illumination (Coolled, pE-100) and an exi blue fluorescence microscopy camera (Q-imaging, Model number: EXI-BLU-R-F-M-14-C, CA) coupled to a light microscope (Nikon, FN1) placed on an active vibration isolation table (Table Stable;  Cell cultures were used 9-15 days after plating, i.e. 7-13 days after viral transfection.Intracellular whole-cell current-clamp recordings were performed as described before [17,18,19].In short, pyramidal neurons were visually selected according to their somatic shape and dendritic morphology.

Method Description
The main algorithm uses a sliding time window over the frames composing the film in single-pixel resolution.In the following, the time series of a single pixel will be denoted with P , and the value of P at time (frame) t will be denoted with P t .The total number of frames (and the total duration of the time series) is T .We compute a linear autoregression of order k on the single pixel time series P .Calling the k parameters of the linear model U i , and the residuals ϵ, we have T equations and k parameters (one equation with k parameters for each time point): For k < T , there are more equations than parameters: this is an overdetermined system and it can not be solved exactly.To find the set of parameters U = U 1 , ..., U k that minimizes the residuals, we use a simple linear least-squares method [20].Letting P be a (T × k) matrix, having the i-lagged time series P t−i as the i-th column, U the vector of the k model coefficients (k × 1) and E the vector of the residuals (T × 1), we can rewrite the system in matrix multiplication form: P * = P U + E. To find the set Ū that satisfies with the assumption that the columns of P are linearly independent, Ū can be computed as: Ū = ( P T P ) −1 P T P * which is solvable by QR decomposition.Having found the coefficients of the linear regression Ū , we construct the time series of the residuals E(T × 1) by comparing the measured pixel time series with the ones predicted by the model.Finally, we compute the average over time of the residuals time series E, which constitutes the value of a single pixel in an ARES frame: The residuals can be interpreted as noise or unpredictable (novel) events.
Thanks to averaging, the noise tends to vanish, especially with longer time windows, while any transient increase or decrease in the signal is still detected.Another way to interpret the average residual time series is to consider that we are using a certain time window (of size T ) of the signal to predict its future: the less the future time points can be predicted, the more the residuals will differ from zero.

Performance of ARES relative to ∆F/F
We introduce ARES, a new method for the quantitative analysis of functional calcium imaging data; it does not require any arbitrary assignment of a baseline and deals with the non-linearity in calcium signals.The method is based on a linear prediction of the signal by using its past (autoregression): the prediction error indicating the presence of significant variations, and therefore of an informative signal.This is especially true for calcium signals, as most biological diffusion phenomena, are essentially exponential increases, followed by exponential decays.To avoid considering the error in predicting noise as part of the resulting signal, the prediction errors are averaged over a short temporal window, which contributes to the 'lag' between ARES and the raw signal or ∆F/F.This lag decreases with increasing the intensity of the calcium response of a cell (Fig 1 A, C).By using the past of a pixel's time series to predict its future, we take advantage of the time dimension of the data, including information on the known properties of calcium dynamics (exponential rise and decay), resulting in a better signal to noise ratio (SNR), an improved granularity, and an increased sharpness (amplitude/duration of the signal) with respect to the standard ratiometric normalization ∆F/F methods (Fig 1B , C).

Effects of process order and time window
ARES has two main hyper-parameters that need to be determined before the analysis: the order of the autoregressive process used in the autoregression (process order), and the length of the time window on which the analysis is performed (see Method description for details).The first (process order) influences the sensitivity of ARES to shorter or faster dynamics: a

Discussion
Ratiometric methods such as ∆F/F are intrinsically limited by the necessity of a comparison with a single baseline value.We herein introduced a method, called ARES, that extrapolates the temporal dynamics of the signal in a pixel-wise fashion, and uses this knowledge to detect transient events.
We showed that the method improves the temporal and spatial resolution of calcium imaging signals, while maintaining the same linearity properties as ∆F/F (Fig 1C ).As an example application, we showed that ARES allows for the study of intracellular calcium dynamics, in the form of wave propagation from spines into dendrites.
The ARES signal is characterized by small positive increments in the signal, representing sudden rises in the calcium signal, and by larger negative troughs, representing the beginning of the signal decay.Once the decay in calcium concentration becomes linearly predictable, ARES takes a value of zero as there is no longer any difference between the predicted and actual signal trace.Linear regressions are not able to correctly capture the dynamics of calcium, which are defined by exponential rises or decays, but as the linear prediction fails whenever there is a novel calcium event, the non-linear dynamics are captured in the residuals.In considering the residuals themselves, instead of the regression coefficients, anything that is not fitted by the linear regression (including the calcium transients) is now represented in the processed signal.
The improved signal resolution and transient detection of ARES, together with its baseline independence and time locality, make the method ideal for the detection of small and non-stimulus locked events, and for studying their dynamic properties.Furthermore, the toolbox accompanying the method provides an intuitive user interface for the use with in-vitro recordings, and an explicit example of the application of the method.The super-resolution imaging of signal transients might offer other interesting advantages such as the differentiation between intra-and extra-cellular calcium sources, and a better reconstruction of the spiking activity, which remains a major challenge in the field.[6,7].

Supplementary Material
Supplemental Video 1. ARES allows the tracking of a dendritic calcium wave front The video shows imaging of a dendrite taken from primary neuronal cultures (see Material and Methods) at 10 Hz.Four videos are presented side by side for comparison: the raw signal, a relative heatmap, and the heatmaps of ∆F/F and the positive part of the ARES analyzed signal, with tracking of the wavefront displayed as a small red dot.2) The amplitude of each trace, normalize to the maximum signal reached for comparison.
3) The duration of the response to the stimuli provided.4) The rising slope of the signals.5) The signal to noise ratio (noise computed by taking a sample of 100 pixels outside the cell soma and dendrites).6) The coefficient of variation (CoV, measure of granularity of the signal, more negative is a more granular signal).Data shown are the averages and standard errors over 6 cells recordings.Each cell value is, in turn, the average over 3 repetitions of the same stimulation protocol.
TS-150).The data was acquired at full field at 10 fps in MicroManager (https://micro-manager.org/) at 14-bit with a readout time of 30 MHz.Data acquisition was triggered using a custom Arduino interface, coupling the whole-cell recording (see below) software Patchmaster (HEKA) with MicroManager, the program controlling the camera and the excitation light source.Electrical recordings, electrical stimulation and calcium imaging were time aligned with clock signals generated in the Patchmaster.2.1.4.Whole-cell intracellular recording and stimulation.
smaller process order corresponds to an increased sensibility to shorter transients (Fig 2).The latter, that is the number of consecutive frames taken into consideration in each autoregression (size of the time window), has a similar effect as the process order (Fig 3).Varying the window size has a linear effect on the signal: shortening the time window results in sharper signal and improved SNR, but choosing time windows that are too short (<10 frames recorded at 10 Hz) does not allow capturing the entire duration of the calcium transient, hence ultimately reducing the SNR.The parameters can be easily adjusted with the open-source toolbox provided herein (https://github.com/DepartmentofNeurophysiology/ARES).3.3.Tracking wavefrontsSpatial projections of ARES (see Fig1B) suggest that ARES also provides high-resolution spatial information about the source of the signal.Therefore we tested the algorithm for tracking calcium wavefronts.In-vitro bright field calcium imaging data were collected in primary neuron cultures as described (see Materials and Methods).ARES presents an increased SNR and at the same time a signal that maintains the same amplitude of the calcium signal and the ∆F/F response, but without a very long tail decay (Fig1).This allows us to observe and measure the intracellular calcium influx in the form of a propagating wave coming from a dendritic spine (Fig4and Supplementary Video 1).Unlike ARES, ∆F/F is not able to resolve how the propagation of the wavefront along the dendrite.

Figure 1 :
Figure1: ARES signal and comparison with ∆F/F and raw signals.(A) A current pulse (duration: 500ms; intensity: increasing with increments of +40pA) was injected every 7 seconds.The traces displayed are the average projection on the time axis of each calcium response over the whole region of interest (ROI).Each trace is normalized to the global max in the time series for direct comparisons.(B) Average projection of the absolute value of each signal over the whole time lapse in the ROI.(C) Comparison between ARES, ∆F/F and the raw signal.From top-left to bottom right: 1) the lag between a detectable response in the signal and the stimulus onset (1 frame = 0.1 seconds).2) The amplitude of each trace, normalize to the maximum signal reached for comparison.3) The duration of the response to the stimuli provided.4) The rising slope of the signals.5) The signal to noise ratio (noise computed by taking a sample of 100 pixels outside the cell soma and dendrites).6) The coefficient of variation (CoV, measure of granularity of the signal, more negative is a more granular signal).Data shown are the averages and standard errors over 6 cells recordings.Each cell value is, in turn, the average over 3 repetitions of the same stimulation protocol.

Figure 2 :
Figure2: Effects of process order (p).The order of the autoregressive process used to predict the signal is a fundamental parameter for the ARES computation.Here we show how its values influence notable characteristics of the output, depending on the stimulation (and therefore response) level, and compare it with the raw signal and ∆F/F.In general, with the exception of the lower possible process order of 1, the results are quite stable and comparable.Most notably, the amplitude of the response (C) from ARES and ∆F/F are the same across all stimulus levels.(A) The number of frames between the stimulus and a detectable response in the signal.(B) The duration in frames of the calcium response.(C) The amplitude of the signal, normalized to the absolute maximum calcium level.(D) The amplitude of the signal normalized to the very first detectable response of the stimuli series.(E) The slope of the rising phase of the signal.(F) The slope of the decay phase.(G) The signal to noise ratio (SNR).(H) The granularity measure of the image.The process order we used by default in other figures was 3.

Figure 3 :
Figure 3: Effects of window size (w).The length of the time window considered for the signal prediction is a fundamental parameter for the ARES computation.Here we showed how its length influence notable characteristics of the output, depending on the stimulation level.In general, the effect of varying the window size is mostly linear, with the exception of the amplitude of the response (C) which remains stable across all stimulus levels.(A) The number of frames between the stimulus and a detectable response in the signal.(B) The duration in frames of the calcium response.(C) The amplitude of the signal, normalized to the absolute maximum calcium level.(D) The amplitude of the signal normalized to the very first detectable response of the stimuli series.(E) The slope of the rising phase of the signal.(F) The slope of the decay phase.(G) The signal to noise ratio (SNR).(H) The granularity measure of the image.The window length we used in our other examples by default in our figures was 25 frames.

Figure 4 :
Figure 4: Spatial analysis of subcellular calcium dynamics with ARES.(A) The ROI in a dendritic spine.(B) The propagation of a calcium wave in consecutive frames (see Supplemental Film 2).Top row: the color code represents ARES values.The red circles indicate the position of the center of mass of the activity (the wavefront).Bottom row: corresponding images analyzed with ∆F/F.(C) The 2-D movement tracking of the wave front, tracked as the position of their center of mass, and (D) its trajectory together with time axis.Different shading of the same color represent different events recorded from the same spine.In comparison to ARES, ∆F/F (red) is unable to pick up any dynamic of the wave propagation.
Centers of Mass of activity -Movement.ARES CM Trajectory -1st event ARES CM Trajectory -2nd event ARES CM Trajectory -3rd event F/F CM Trajectory -1st event F/F CM Trajectory -2nd event F/F CM Trajectory -