|
nif:isString
|
-
Worm preparation and experiment design: All worms were grown and maintained under standard conditions [42], incubated with food at 20°C. Well fed worms were washed twice then gently spun down for 1 minute and the supernatant discarded by aspiration. We discovered empirically that ibuprofen affects the heat-induced escape response in our assay. For the drug application 100 μL of ibuprofen in M9 at 100 μM was added to the eppendorf tube. For the wild-type and mutant data set, M9 was used instead of the drug solution. Worms were then placed in an incubator for 30 minutes at 20°C. After that worms were poured onto a seeded agar plate and transferred to agar assay plates by a platinum wire pick. These assay plates were incubated at 20°C for 30 minutes, and then the experimental trials were done within the next 30 minutes. In total N = 201 worms for the control group, N = 441 worms for the ibuprofen group, and N = 100 worms for the mutant group (ocr-2(ak47) osm-9(ky10) IV; ocr-1(ak46)) group were tested. The mutant strain was obtained from the Caenorhabditis Genetics Center. The heat stimulation instrument has been described previously [29]. In summary, an infrared laser is directed to heat the head of a freely crawling worm (∼0.5mm FWHM) on an agar plate. The laser pulse is generated with a randomly chosen laser current between 0 to 200 mA, with a duration of 0.1 s. The heating of the worm is nearly instantaneous, and it is directly proportional to the current, between 0 and 2°C for the current range used in our experiments. The temperature change at 60 mA current is 0.4°C ± 0.03°C, 100 mA current is 0.89°C ± 0.05°C and 150mA current is 1.4°C ± 0.2°C. Worms were stimulated only once and not reused. The movements of the worms are imaged using a standard stereomicroscope with video capture and laser control software written in LabVIEW. For each stimulus trial, a random worm is selected on the plate and its motion is sampled at 60 Hz for 15 s, and the laser is engaged 1 s after the start of the video recording.
The recorded response videos were then processed with Matlab to calculate the time series of the worm centroid motion as described previously [29]. All the worms that were not stimulated near the head or were not moving forward in the beginning of the video were discarded. Numerical derivatives of the centroid times series were then taken and filtered with a custom 500 ms Gaussian filter, which was a one-sided Gaussian at the edges of the recorded time period, becoming a symmetric Gaussian away from the edges. This removed the noise due to numerical differentiation and also averaged out the spurious fluctuations in the forward velocity due to the imperfect sinusoidal shapes of the moving worm. We verified that different choices of the filter duration had little effect on the subsequent analysis pipeline. The direction of the velocity was determined by projecting the derivative of the centroid time series onto the head-to-tail vector for each worm, with the positive and negative velocity values denoting forward / backward motion, respectively. The filtered velocity profiles needed to be subsampled additionally. This was because the statistical model of the data, Eq (8), involved covariance matrices of the active and paused velocity profiles, Σp and Σa (note that velocity profiles are not temporally translationally invariant due to the presence of the stimulus, thus the full covariance matrix is needed, and not a simpler correlation function). To have a full rank covariance matrix, the number of trials must be larger than the number of time points. Alternatively, regularization is needed for covariance calculation. The autocorrelation function for all three worm types showed a natural correlation time scale of ≳ 0.2 s, whether the data was pre-filtered or not. Thus subsampling at a frequency > 5 Hz would not result in data loss. Therefore, instead of an arbitrary regularization, we chose to subsample the data at 12 Hz, leaving us with 37 data points to characterize the first 3 s of the worm velocity trace after the stimulus application, 1 ≤ t ≤ 4 s since the start of the trial. Eq (8) additionally needs knowledge of T, the number of effectively independent velocity measurements in the profile. This is obtained by dividing the duration of the profile by the velocity correlation time. An uncertainty of such procedure has a minimal effect on the model of the experiment since it simply changes log likelihoods of models by the same factor, not changing which model has the maximum likelihood. We then considered limiting the duration of the velocity profile used in model building: if velocities at certain time points do not contribute to the identification of I, they should be removed to decrease the number of unknowns in the model that must be determined from data (values of the templates at different time points). The first candidate for removal was the period of about 10 frames (0.16 s) after the laser stimulation since worms do not respond to the stimulus so quickly. However, removal of this time period had a negligible effect on the model performance, and we chose to leave it intact. In contrast, starting from 3.3 s (2.3 s after the stimulus) the fraction of explainable variance drops to nearly zero (Fig 5) since many worms already had turned by this time and resumed forward motion. Therefore, we eventually settled on the time in the 1.0…3.3 s range for building the model. Whenever needed, we estimated the variance of our predictions by bootstrapping the whole analysis pipeline [43]. For this, we created 1000 different datasets by resampling with replacement from the original control dataset and the mutant / ibuprofen datasets. Control statistical models (the scaling function f and the velocity templates) were estimated for each resampled control dataset. Standard deviations of these models were used as estimates of error bars in Fig 4. For Fig 7, we additionally needed to form the closest control / treatment worm pairs. These were formed between the resampled data sets for all worm types as well. Standard deviations of ΔItype evaluated by such resampling were then plotted in Fig 7 and used to estimate Z scores. Note that such resampling produces control / treatment paired worms that have slightly larger current differences than in the actual data; this leads to our error bars being overestimates. Model in Eq (1) requires knowing P(I). In principle, this is controlled by the experimentalist, and thus should be known. In our experiments, P(I) was set to be uniform. However, as described above, some of the worms were discarded in preprocessing, and this resulted in non-uniformly distributed current samples. To account for this, we used the empirical Pemp(I) in our analysis instead of P(I) = const. In turn, Pemp(I) was inferred using a well-established algorithm for estimation of one-dimensional continuous probability distributions from data [44]. All of this analysis was implemented using Matlab, and the code is available for download from a public GitHub repository https://github.com/EmoryUniversityTheoreticalBiophysics/C.-elegans.
Calculating the template velocities, the covariances, and the scaling function: The template for the paused state up is calculated by taking the average of all paused velocity profiles for each of the three worm datasets. The covariance Σp is then the covariance of the set of the paused velocity profiles. For active worms, we start with fixed putative parameter values I1 and I2. We then calculate the active template ua and the covariance matrix ∑a by maximizing the likelihood in Eq (7) ∂∑iNtype,alogP(vi|a,Ii)∂ua∝∑iNtype,avifI1,I2(Ii)-uafI1I22(Ii)=0,(14) ∂∑iNtype,alogP(vi|a,Ii)∂Σa∝∑iNtype,avi-uafI1,I2(Ii)2-(Σa)-1=0,(15) where Ntype,a is the number of active worms of the analyzed type. This gives: ua(I1,I2)=∑i=1Ntype,avifI1,I2(Ii)∑i=1Ntype,afI1,I22(Ii),(16) Σa=∑iNtype,avi-uafI1,I2(Ii)2. (17) Having thus estimated ua and Σa at fixed parameter values I1, I2, we maximize ∏i P(vi|a,Ii) over the parameters using standard optimization algorithms provided by MATLAB. We perform optimization from ten different initial conditions to increase the possibility that we find a global, rather than the local maximum.
|
![[This page as RDF]](http://data.gesis.org/somesci/static/rdf-icon.gif){kind=icon}