Data analysis

This protocol is extracted from research article:

The somatosensory cortex receives information about motor output

**
Sci Adv**,
Jul 10, 2019;
DOI:
10.1126/sciadv.aaw5388

The somatosensory cortex receives information about motor output

Procedure

We built models to reconstruct activity in S1, temporal changes in EMG, or joint kinematics using a partial dataset (training dataset) and tested them using the remainder of the same dataset (test dataset). For reconstruction of S1 activity and joint kinematics, we partitioned continuously recorded data of each session into 24 blocks (one block for 25 s of data). Among the 24 blocks, 20 randomly selected blocks were used for the training dataset and the remaining 4 blocks were used for the test dataset. For reconstruction of EMG signals, five-sixths of the full trials of each session were randomly selected as a training dataset and the remainder were selected as the test dataset. To assess the model, we calculated the correlation coefficient between observed data and their reconstruction in the test dataset. We also calculated variance accounted for (VAF) as follows$$\text{VAF}=1-\frac{\sum {(y(t)-f(t))}^{2}}{\sum {(y(t)-\overline{y(t)})}^{2}}$$(5)where *y*(*t*) is a vector of the actual activity in S1 at time index *t*, $\overline{y(t)}$ is the mean of *y*(*t*), and *f*(*t*) is the reconstructed activity at time index *t*. We performed sixfold cross-validation in the analysis of each session and used averaged values for the analysis. Then, we calculated averaged values of each electrode and each kinematic variable from data taken from 21 (Monkey T) and 7 (Monkey C) sessions in decoding S1 activity and joint kinematics and averaged values of each muscle from data taken from 24 (Monkey T) and 40 (Monkey C) sessions in decoding EMG activity. In control analyses of the model reconstruction, we created surrogate training datasets in which we shuffled temporal profiles of inputs independently across different blocks to generate a model and subsequently tested the model.

We built a model to predict the initial peak amplitude of EMG using a training dataset and tested it using a test dataset. Five-sixths of the full trials were randomly selected as a training dataset and the remainder were selected as the test dataset (Monkey T: 2847 trials for training, 570 trials for test; Monkey C: 4234 trials for training, 847 trials for test). To assess the model, we calculated the correlation coefficient between the observed and reconstructed EMG amplitude in the test dataset. We also calculated the largest possible variance (the first principal component) between the observed and reconstructed EMG amplitude to obtain the slope of the regression line. We drew regression lines that passed the centroid of the data in the plot (Fig. 6A). Among models of high-γ activity in M1 or S1 in sliding time windows in the premovement period (−500 to 0 ms around movement onset), we selected a model that reconstructed the initial peak amplitude of each muscle at the highest accuracy (Fig. 6, B and C). We performed sixfold cross-validation in the analysis. In control analyses of the model reconstruction, we created surrogate training datasets in which we randomized the trials of input array (trial shuffling) to generate a model and subsequently tested the model.

We built a model to predict the initial peak amplitude of EMG from M1, S1, or peripheral afferents activity using a training dataset and tested it using a test dataset. We first determined putatively the same units among different sessions according to the shape of waveforms and the distribution of interspike intervals. Thirty-one units were identified as putatively the same units among nine different sessions in Monkey T. Five-sixths of the full trials were randomly selected as a training dataset and the remainder were selected as the test dataset (1091 trials for the training, 219 trials for the test). To assess the model, we calculated the correlation coefficient and the slope of the regression line between the observed and reconstructed EMG amplitude in the test dataset. We performed sixfold cross-validation in the analysis. In control analyses of the model reconstruction, we created surrogate training datasets in which we randomized the trials of input arrays (trial shuffling) to generate a model and subsequently tested the model.

To obtain the onset time of the activity or the reconstruction, we first calculated the average of the aligned waveform in a test dataset [Monkey T: 21.4 ± 0.7 trials (mean ± SD, *n* = 21 sessions); Monkey C: 19.8 ± 1.0 trials (*n* = 7 sessions)]. Then, we calculated a threshold of the averaged aligned waveform by adding an average of activity during the baseline period (−1250 to −250 ms around movement onset) to one-fifth of the amplitude of the activity from 250 ms before to 100 ms after movement onset. If the activity was over the threshold in five consecutive bins, then the first of these bins was set as the onset of the activity. The calculated onset corresponded well with that based on visual inspection. We calculated the average values of the onset from six test datasets in one session and finally obtained their average from the whole sessions. We also used the average activity during the baseline period to calculate the area above the baseline. To obtain the end point of the activity in S1, we first calculated a threshold of the aligned waveform by adding an average of activity during the baseline period (−1250 to −250 ms around movement onset) to ^{3}/_{10} of the amplitude of the activity from 250 ms before to 100 ms after movement onset. If the activity 1000 ms after movement onset was below the threshold in five consecutive bins, the first of these bins was set as the end point of the activity in S1. The end point of the activity in S1 was 1307 ± 27 ms (mean ± SD, *n* = 16 signals; eight electrodes, two bands) after movement onset for Monkey T and 1450 ± 51 ms (mean, *n* = 12 signals; six electrodes, two bands) for Monkey C.

To calculate the point when M1 or S1 started to encode EMG burst amplitude (Fig. 6), we used a common threshold in the encoding of M1 and S1. For both the correlation coefficient and the slope of the regression line, we calculated the encoding magnitude by subtracting the baseline activity during the baseline period (−845 to −255 ms around movement onset) from the magnitude from 245 ms before to 255 ms after movement onset. Then, we calculated the ^{1}/_{5} (correlation coefficient) or ^{1}/_{20} (slope) of the encoding magnitude and used a larger one of the values for M1 and S1 as a common threshold. Last, we added the common threshold to the baseline activity for M1 and S1 and used this value as respective thresholds. If the indices were over the threshold in three consecutive bins, then the first of these bins was set as the point when M1 or S1 started to encode EMG burst amplitude. Reconstruction accuracy (slope) of the flexor digitorum profundus muscle of Monkey C was only 0.02 in the decoding from M1 and 0.04 from S1, so we did not calculate the point when M1 or S1 started to encode burst amplitude for this muscle.

To calculate the point when M1, S1, or peripheral afferents started to encode EMG burst amplitude (fig. S8), we used a common threshold for their encoding. For both the correlation coefficient and the slope of the regression line, we calculated the encoding magnitude by subtracting the baseline activity during the baseline period (−845 to −255 ms around movement onset) from the magnitude from 245 ms before to 255 ms after the movement onset. Then, we calculated the ^{1}/_{5} (correlation coefficient) or ^{1}/_{20} (slope) of the encoding magnitude and used a larger one of the values for M1 and S1 as a common threshold. Last, we added the common threshold to the baseline activity for M1, S1, and peripheral afferents and used this value as respective thresholds. If the indices were over the threshold in three consecutive bins, then the first of these bins was set as the point when M1, S1, or afferents started to encode EMG burst amplitude.

To calculate the variability of the averaged profile of the MCx and Afferent components, a principal component analysis (PCA) was performed on the MCx and Afferent components throughout the premovement and movement periods (−500 to 2000 ms around movement onset) of eight (Monkey T) or six (Monkey C) electrodes. Values for high-γ 1 and 2 were averaged before PCA was conducted. More than 90% of the total variance was captured by the first principal component for the MCx component and by the first two principal components for the Afferent component.

To obtain a weight value given to each electrode, weight vector *w*_{j,k,l} in Eq. 1 or 3 in the manuscript was averaged across time points. Values for high-γ 1 and 2 were averaged. The one-way ANOVA was used to determine whether there are any statistically significant differences between the means of weight values of different electrodes.

Note: The content above has been extracted from a research article, so it may not display correctly.

Q&A

Your question will be posted on the Bio-101 website. We will send your questions to the authors of this protocol and Bio-protocol community members who are experienced with this method. you will be informed using the email address associated with your Bio-protocol account.