Decoding locomotion from population neural activity in moving C. elegans

  1. Kelsey M Hallinen
  2. Ross Dempsey
  3. Monika Scholz
  4. Xinwei Yu
  5. Ashley Linder
  6. Francesco Randi
  7. Anuj K Sharma
  8. Joshua W Shaevitz
  9. Andrew M Leifer  Is a corresponding author
  1. Department of Physics, Princeton University, United States
  2. Princeton Neuroscience Institute, Princeton University, United States
  3. Lewis-Sigler Institute of Integrative Genomics, Princeton University, United States
8 figures, 6 tables and 1 additional file

Figures

Figure 1 with 3 supplements
Population calcium activity and tuning of select neurons during spontaneous animal movement.

Recording AML310_A. (a) Calcium activity of 134 neurons is simultaneously recorded during locomotion. Activity is displayed as motion-corrected fluorescent intensity Fmc. Neurons are numbered according to agglomerative hierarchical clustering. White space indicates time-points where neural tracking failed. (b) Body bend velocity and body curvature derived from an eigenvalue decomposition, and (c) position on the plate during recording are shown. (d) Example neurons significantly tuned to velocity. Examples are those with the highest Pearson’s correlation coefficient in each category: activity (or its derivative) with positive (or negative) correlation to velocity. P-values are derived from a shuffling procedure that preserves per-neuron correlation structure. All tuning curves shown are significant at 0.05% after Bonferroni correction for multiple hypothesis testing (p<1.9×10-4). Boxplot shows median and interquartile range. Blue or orange shaded circles show neural activity at each time point during behavior. (e) Example neurons highly tuned to curvature were selected similarly. No neurons with negative dF/dt tuning to curvature passed our significance threshold.

Figure 1—figure supplement 1
Additional details and examples of velocity tuning.

Additional examples of velocity tuning curves for Fmc (top) and dFmc/dt (bottom) from recording AML310_A are shown. The correlation coefficient ρ captures the relation between each neuron’s activity and curvature and is plotted with a corresponding p-value (middle). The p-value is calculated from a time-lag shuffle and tests the null hypothesis that such a correlation would be found due to chance. Dashed line indicates a significance level of 0.05 after a Bonferonni correction for multiple hypothesis testing (p<1.9×10-4).

Figure 1—figure supplement 2
Additional details and examples of curvature tuning.

Additional examples of curvature tuning for Fmc (top) and dFmc/dt (bottom) from recording AML310_A are shown. The correlation coefficient ρ captures the relation between each neuron’s activity and velocity and is plotted with a corresponding p-value (middle). The p-value is calculated from a time-lag shuffle and tests the null hypothesis that such a correlation would be found due to chance. Dashed line indicates a significance level of 0.05 after a Bonferonni correction for multiple hypothesis test (p<1.9×10-4).

Figure 1—figure supplement 3
Number of significantly tuned neurons across recordings.

Pearson’s correlation coefficient ρ was calculated for each neuron in 11 GCaMP recordings and 11 GFP control recordings that lacked a calcium indicator. Neurons were counted as significantly tuned if their Pearson’s correlation coefficient passed a recording-specific multiple-hypothesis corrected significance test and exceeded an absolute value of 0.4. Bar shows median. Whiskers show inter-quartile range.

Figure 2 with 1 supplement
Neuron pair AVA is active during backward locomotion and exhibits expected tuning during moving population recordings.

(a) AVAR and AVAL are labeled by BFP under a rig-3 promoter in strain AML310. Two optical planes are shown from a single volume recorded during movement. Planes are near the top and bottom of the optical stack, corresponding to the animals’ extreme right and left. The recording is the same as in Figure 1. Top row shows BFP. Bottom row shows RFP in the nuclei of all neurons. Segmented neurons centered in the optical plane are labeled with ⊕, while neurons from nearby optical planes are labeled with ○. Arrow indicates AVAR or AVAL. Numbering corresponds to Figure 1a. (b) Calcium activity of AVAR and AVAL during locomotion in recording AML310_A, same as in Figure 1. (c) Aggregate tuning of AVA across four individuals (seven neurons). Boxplot shows median and interquartile range. Lightly shaded blue or orange circles show activity at each time point during behavior.

Figure 2—figure supplement 1
Sum of AVAL and AVAR activity.

Activity of AVAL and AVAR from AML310_A in Figure 2b are shown summed together. This permits comparison to recordings that do not resolve the two neurons separately.

Figure 3 with 5 supplements
Population neural activity decodes locomotion.

(a–d) Performance of the best single neuron (BSN) is compared to a linear population model in decoding velocity and body curvature for the exemplar recording AML310_A shown in Figure 1. (a) Predictions on the held-out test set are compared to measured velocity. Light green shaded region indicates held-out test set. Red arrows indicate examples of features that the population captures better than the BSN. (b) Performance is reported as a coefficient of determination RMS2 evaluated on the mean-subtracted held-out test data (green points). (c,d) Model predictions are compared to measured curvature. (e) Performance of velocity decoding is shown for recordings of n=11 individuals (strain AML310 and AML32) and for recordings of n=11 GFP control animals lacking a calcium indicator (strain AML18). Two-sided Wilcoxon rank test is used to test significance of population performance compared to BSN, p=3.9×103. Welch’s unequal variance t-test is used to test significance of population performance compared to GFP control, p=3.2×10-2. (f) Performance of curvature decoding is shown for all recordings. Each recording is colored the same as in e. p=3.2×10-2 and p=1.8×10-3 for comparisons of population performance to that of BSN, and GFP control, respectively.

Figure 3—figure supplement 1
Performance correlates with maximal GCaMP Fano Factor, a metric of signal.

Decoding performance is plotted against maximal GCaMP Fano Factor for each recording for velocity and curvature. Maximal GCaMP Fano Factor is the Fano Factor of the raw GCaMP activity for the neuron in each recording with the highest Fano Factor, maxi(σ2[Fi,GCaMP]μ[Fi,GCaMP]). Labels for each recording are shown. Dashed red line is the line of best fit (correlation coefficient between fit and data is ρ=0.46 for velocity and ρ=0.49 for curvature).

Figure 3—figure supplement 2
Neural activity and behavior for all moving calcium imaging recordings (AML310 and AML32).
Figure 3—figure supplement 3
Neural activity and behavior for all moving GFP control recordings (AML18).

Neural activity and behavior for all moving GFP control recordings (AML18).

Figure 3—figure supplement 4
Alternative population models.

Performance of alternative population models for decoding velocity. Traces are shown for exemplar recording AML310_A. Mean across all moving GCaMP recordings is also listed. Gray shading shows held-out test set. (a) The population model used throughout the paper. This model uses ridge regression with fluorescence signals and their temporal derivatives as features. (b) A linear model using ridge regression, with only fluorescence (not temporal derivative) signals as features. (c) A linear model using fluorescence signals and their temporal derivatives as features, regularized with a combination of a ridge penalty and the squared error of the temporal derivative of behavior. (d) The model in c, but using only fluorescence signals as features. (e) A linear model using fluorescence signals and their temporal derivatives as features, regularized with an ElasticNet penalty with an L1 ratio of 10-2. (f) The model in e, but using only fluorescence signals as features. (g) The multivariate adaptive regression splines (MARS) model, using fluorescence signals and their temporal derivatives as features. (h) A linear model together with a shallow decision tree, using fluorescence signals and their temporal derivative as features.

Figure 3—figure supplement 5
Nonlinear fits using best single neuron.

Performance of polynomial regression models for decoding velocity on 11 GCaMP recordings using the best single neuron. The best single neuron is defined as the one with the best decoding performance using a linear model on the training data. We compare the performance of these regression models to that of the population model using a two-sided Wilcoxon rank test. The population model significantly (p < 0.05) outperforms polynomial regression models up to fourth order.

Example where population decoded a fuller range of animal behavior.

(a) The decoding from the best single neuron and the population model are compared to the measured velocity for example recording AML32_A. (b) Predictions from the best single neuron saturate at a velocity of approximately 0.1 mm s−1.

Figure 5 with 3 supplements
Weights assigned to neurons by the population model in the exemplar recording, and their respective tuning.

(a) The weight W assigned to each neuron’s activity (Fmc) or its temporal derivative (dFmc/dt) by the velocity population decoder is plotted against its Pearson’s Correlation coefficient ρ which characterizes its tuning to velocity. Recording AML310_A is shown, same as in Figure 1. Dashed red line shows line of best fit. Right panel shows the observed distribution of weights. A zero-mean Gaussian with standard deviation set to the empirically observed standard deviation is also shown. (b) Same as in a, but for curvature. (c) Tuning and activity of the top highest amplitude weighted neurons for velocity is shown. Activity of each neuron is time aligned to the observed behavior (top row). Neurons are labeled corresponding to their number in the heatmap in Figure 1. Their rank and weight W in the decoder is listed. Red arrows highlight peaks in the temporal derivative of activity of neuron #24 and #110, while cyan arrows highlight peaks of neuron #44. Y- and X-axes labels and scales are preserved within individual rows and columns, respectively. Light green shading indicates the held-out portion of the recording. (d) Same as c but for curvature. Red and cyan arrows show two sets of deep ventral bends that are captured by different neurons. Green arrows show dorsal bends.

Figure 5—figure supplement 1
Comparison of weights assigned to a neuron’s activity versus its temporal derivative.

Comparison of weights assigned to a neuron’s activity versus its temporal derivative for velocity (left) or curvature (right) decoders. Comparison of weights assigned to a neuron’s activity |WF| by the population decoder, versus the weights assigned to its temporal derivative |WdF/dt| for each neuron in recording AML310_A.

Figure 5—figure supplement 2
Comparison of weights assigned for decoding velocity vs decoding curvature.

Comparison of weights assigned for decoding velocity vs decoding curvature. (a) The magnitude of the weight assigned to each neuron in recording AML310_A for velocity |Wvel| is compared to the magnitude of its assigned weight for curvature |Wcurv|. Each neuron is plotted twice, once for the weight assigned to its activity and once for the weight assigned to the temporal derivative of its activity. (b) Same as in (a), except here the higher weight of either activity or its temporal derivative is plotted.

Figure 5—figure supplement 3
Example traces of highly weighted neurons used to decode curvature in AML32_A.

Traces of top five highest weighted neurons used to decode curvature in AML32_A. Same recording as in Figure 4. Arrows indicate activity peaks corresponding to ventral (blue shades, top) or dorsal turns (red shades, bottom). Different neurons contribute activity peaks to different sets of turns.

Figure 6 with 1 supplement
Number of neurons needed by the model to decode velocity and curvature.

(a) The minimum number of neurons needed for a restricted model to first achieve a given performance is plotted from recording AML310_A in Figure 1. Performance, RMS,all2 is reported separately for velocity (blue) and curvature (green) and is calculated on the entire recording (test and train). Intersect refers to the intersection of the set of neurons included in both partial models (velocity and curvature) for a given performance. Red dashed line, N90, indicates number of neurons needed to achieve 90% of full model performance. (b) N90 is computed for velocity and curvature for all recordings. The number of neurons present in both populations at 90% performance level (intersection) is shown. Box shows median and interquartile range. (c) N90 for all recordings is shown plotted versus the performance of the full population velocity or curvature decoder, respectively. Number of intersection neurons (red ’x’) is plotted at the higher of either the velocitys or curvature’s performance.

Figure 6—video 1
Animation showing partial model performance as neurons are added, corresponding to Figure 6a.

Top panel shows performance RMS,all2 evaluated on both test and training set. Bottom left shows measured velocity (black) and decoded velocity (blue). Gray shading indicates test set. Bottom right shows measured velocity compared to decoded velocity for training (blue) and test set (green).

Figure 7 with 2 supplements
Immobilization alters the correlation structure of neural activity.

(a) Calcium activity is recorded from an animal as it moves and then is immobilized with a paralytic drug, recording AML310_E. (b) Activity of AVAL and AVAR from (a). (c) Population activity (or its temporal derivative) from (a) is shown projected onto its first three PCs, as determined by only the immobilized portion of the recording. (d) Neural state space trajectories from (c) are shown split into moving and immobile portions, color coded by time. Scale, axes and PCs are the same in both plots. (e) Pairwise correlations of neural activity ρi,j are shown as heatmaps for all neurons during movement and immobilization, sorted via a clustering algorithm. Top row is sorted to movement, bottom row is sorted to immobilization. (f) Dissimilarity between correlation matrices for moving and immobile portions of a recording are shown compared to the dissimilarity observed between correlation matrices taken at similar time windows within moving-only recordings. Dissimilarity is (ρi,j-ρi,j)2. Dissimilarity was measured in three moving-immobile recordings with paralytic and 11 moving-only recordings. p=1.2×102, Welch’s unequal variance t-test. Boxes show median and interquartile range.

Figure 7—figure supplement 1
Example from additional moving-to-immobile recording.

Calcium activity is recorded from an animal as it moves and then is immobilized with a paralytic drug, recording AML32_H. Activity and behavior. (b) Population activity (or its derivative) from (a) is shown projected onto its first three PCs, as determined by only the immobilized portion of the recording. (c) Neural state space trajectories from (b) are plotted in 3D and shown split into moving and immobile portions. (d) Pairwise correlations of neural activity ρi,j are shown as heatmaps for all neurons during movement and immobilization, sorted via clustering algorithm. Top and bottom rows are sorted to movement or immobilization, respectively.

Figure 7—figure supplement 2
Immobile-only recording.

Calcium activity is recorded from an animal immobilized with nano-beads, recording AML310_G. (a) Calcium activity. (b) Activity of neurons AVAL and AVAR. (c) Population activity (or its temporal derivative) from (a) is shown projected onto its first three PCs, as determined by only the immobilized portion of the recording. (d) Neural state space trajectories from (b) are plotted in 3D.

Correlations with respect to AVAL and AVAR during movement and immobilization.

(a) The Pearson’s correlation of each neuron’s activity to AVAR and AVAL is shown during movement and immobilization. Selected neurons are numbered as in Figure 7 (same recording, AML310_E). Neurons are sorted according to their correlation during movement. (b) Scatter plot shows relation between a neuron’s correlation to AVA during movement and its correlation during immobilization. Gray squares and blue circles indicate correlation to AVAL and AVAR, respectively. (c) On average, neurons become more positively correlated to AVA upon immobilization, p = 0.019 Wilcoxon ranked test. Box shows median and interquartile range. (d) Activity traces of selected neurons are shown time aligned to AVA. Green and purple shading indicate positive or negative correlation to AVA, respectively.

Tables

Table 1
Number of neurons needed to achieve 90% of full model performance, N90, reported as (median ± standard deviation), across all 11 recordings.
Velocity N90Curvature N90Intersection N90Total recorded
27 ± 1631 ± 187 ± 5121 ± 12
Key resources table
Reagent type (species) or resourceDesignationSource or referenceIdentifiersAdditional information
Strain, strain background (C. elegans)AML310this workDetails in Table 2
Strain, strain background (C. elegans)AML32Nguyen et al., 2017RRID:WBI-STRAIN:WBStrain00000192
Strain, strain background (C. elegans)AML18Nguyen et al., 2016RRID:WBI-STRAIN:WBStrain00000191
Table 2
Strains used.

Associated Research Resource Identifiers are listed in Key Resources.

StrainGenotypeExpressionRoleReference
AML310wtfIs5[Prab-3::NLS::GCaMP6s; Prab-3::NLS::tagRFP]; wtfEx258 [Prig-3::tagBFP::unc-54]tag-RFP and GCaMP6s in neuronal nuclei; BFP in cytoplasm of AVA and some pharyngeal neurons (likely I1, I4, M4 and NSM)Calcium imaging with AVA labelThis Study
AML32wtfIs5[Prab-3::NLS::GCaMP6s; Prab-3::NLS::tagRFP]tag-RFP and GCaMP6s in neuronal nucleiCalcium imagingNguyen et al., 2017
AML18wtfIs3[Prab-3::NLS::GFP, Prab-3::NLS::tagRFP]tag-RFP and GFP in neuronal nucleiControlNguyen et al., 2016
Table 3
Recordings used in this study.
Unique IDStrainDuration (mins)Notes
AML310_AAML3104Ca2+ imaging w/ AVA label, moving
AML310_B4
AML310_C4
AML310_D4
AML310_EAML3108Ca2+ imaging w/ AVA label, moving-to-immobile
AML310_F8
AML310_GAML31015Ca2+ imaging w/ AVA label, immobile
AML32_AAML3211Ca2+ imaging, moving
AML32_B11
AML32_C10
AML32_D11
AML32_E4
AML32_F5
AML32_G4
AML32_HAML3213Ca2+ imaging, moving-to-immobile
AML18_AAML1810GFP control, moving
AML18_B10
AML18_C7
AML18_D5
AML18_E5
AML18_F6
AML18_G9
AML18_H6
AML18_I7
AML18_J6
AML18_K6
Table 5
Alternative models explored.

Most are linear models, using either the Ridge or ElasticNet regularization. In some cases, we add an additional term to the cost function which penalizes errors in the temporal derivative of model output (which, for velocity models, corresponds to the error in the predicted acceleration). For features, we use either the neural activities alone, or the neural activities together with their temporal derivatives. We also explore two nonlinear models: MARS Friedman, 1991, and a shallow decision tree which chooses between two linear models.

ModelPenaltyFeaturesNumber of parameters
LinearRidgeF and dF/dt2Nn+1
LinearRidgeFNn+1
LinearRidge + Acceleration PenaltyF and dF/dt2Nn+1
LinearRidge + Acceleration PenaltyFNn+1
LinearElasticNetF and dF/dt2Nn+1
LinearElasticNetFNn+1
MARS (nonlinear)MARSF and dF/dtvariable
Linear with Decision Tree (nonlinear)RidgeF and dF/dt4Nn+9

Additional files

Download links

A two-part list of links to download the article, or parts of the article, in various formats.

Downloads (link to download the article as PDF)

Open citations (links to open the citations from this article in various online reference manager services)

Cite this article (links to download the citations from this article in formats compatible with various reference manager tools)

  1. Kelsey M Hallinen
  2. Ross Dempsey
  3. Monika Scholz
  4. Xinwei Yu
  5. Ashley Linder
  6. Francesco Randi
  7. Anuj K Sharma
  8. Joshua W Shaevitz
  9. Andrew M Leifer
(2021)
Decoding locomotion from population neural activity in moving C. elegans
eLife 10:e66135.
https://doi.org/10.7554/eLife.66135