Removing baseline trends with Increment-Shifted Averages
shows initial SpikeTAs of 6 neuron-muscle pairs selected to illustrate a variety of baseline trends. For each neuron-muscle pair, the left column shows the initial SpikeTA with the vertical scale adjusted to show the entire trace relative to zero average rectified EMG activity (
). The brief peak or trough following the trigger time (vertical line) in these examples represents the increase or decrease in
that results from short-latency synaptic input to the motoneuron pool time-locked to the spikes discharged by the trigger neuron. We refer to these peaks and troughs as SpikeTA effects. Each SpikeTA effect can be seen to ride on a baseline of
, which represents the average ongoing EMG activity around the time the trigger neuron discharged spikes. Although SpikeTAs used for standard analyses typically are 80 ms long, here the initial SpikeTAs are shown for 160 ms to emphasize the underlying baseline trends. Whereas in the uppermost SpikeTA this baseline is reasonably flat (A), the
baseline may show an approximately linear trend (ramp) that either rises (B) or falls (C), or may be overtly curvilinear with convexity either up (D) or down (E). In some cases, the baseline undergoes a sudden inflection at the time a SpikeTA effect might occur (F).
Figure 2 Adjustment of baseline trends with Increment Shifted Averages. A-F show SpikeTAs from 6 different neuron-muscle pairs with various baseline trends: A – flat; B – rising ramp; C – falling ramp; D – convex up; E – (more ...)
These baseline trends represent the average pattern of co-modulation of neuron firing rate and EMG activity. A rising baseline trend (B), for example, indicates that the neuron tended to discharge spikes at times of increasing EMG, such that on average more EMG activity occurred after a spike than before. Conversely, a falling baseline (C) indicates that the neuron tended to discharge spikes at times of decreasing EMG, such that on average more EMG activity occurred before a spike than after. Convex up curvilinear baseline trends (D) indicate that the neuron tended to discharge spikes during phasic bursts of EMG; whereas convex down trends (E) indicate that the neuron tended to discharge spikes between EMG bursts. Inflections in the baseline (F) indicate that the neuron tended to fire spikes when a steady level of EMG activity changed abruptly. The patterns of neuron-EMG co-modulation reflected by these baseline trends occur on a time scale an order of magnitude longer than the spike-locked synaptic inputs responsible for the ~10 – 20 ms peaks and troughs that follow the trigger time, i.e. the SpikeTA effects.
In order to quantify the features of a SpikeTA effect accurately, any baseline trend in the
reflecting long-term neuron-EMG co-modulation should be removed, for two reasons. First, in computing the peak percent increase (PPI) and mean percent increase (MPI), the baseline level of
is used to normalize the amplitude of the peak or trough. A rising or falling trend in the period used to estimate the baseline level of
would alter the values computed for PPI and MPI. For example, a rising trend preceding a SpikeTA effect would result in a lower estimate of the baseline level than a falling trend preceding the same SpikeTA effect. Second, the onset and offset times of the SpikeTA effect are identified when the peak or trough crosses some pre-defined level, typically set at ± 2 standard deviations from the baseline period mean. A rising or falling trend during the baseline sample period would increase the standard deviation of the baseline during that period. The larger standard deviation would result in onset and offset times closer to the maximum or minimum of the peak or trough, respectively, than for the same SpikeTA effect riding on a flat baseline.
Linear rising and falling trends customarily have been removed by calculating a best-fit line through a period of the SpikeTA that samples the baseline trend but not the SpikeTA effect per se (Lemon, Mantel et al., 1986
; McKiernan, Marcario et al., 2000
). This line then is subtracted from the initial SpikeTA, which flattens the baseline but also brings the average value down to zero. A constant (e.g. the original value of the SpikeTA at the trigger time) therefore is added back to the SpikeTA to restore the height of the baseline
above zero. This procedure removes a linear baseline trend while preserving the stochastic noise level of the SpikeTA waveform.
The ramp subtraction approach fails to remove curved baseline trends, however. We therefore developed a method of removing an arbitrary curvilinear baseline trend, while making no assumptions about the shape of that trend or the characteristics of any spike-locked effects. shows simultaneously recorded neural activity from a microelectrode in M1 and EMG activity from EDC. During this period the largest spike in the M1 recording discharged three times, generating three triggers that were used in forming a SpikeTA. These three triggers discriminated from the spike are shown in , along with the three corresponding segments of full-wave rectified EMG. These three segments of rectified EMG were aligned along with 6,011 additional segments around other spike triggers from the same neuron, and averaged to form the initial SpikeTA shown on the left in below.
Figure 3 Increment shifted averaging to estimate the baseline trend. A. A short segment of a microelectrode recording shows three large spikes discharged by a neuron. EMG recorded simultaneously from EDC is shown below. B. Triggers from each of the three neuron (more ...)
Our method of removing an arbitrary curvilinear baseline trend is illustrated in . Eighty-one artificial triggers were created at 1 ms intervals starting 40 ms before and ending 40 ms after each spike of the trigger neuron used to form the initial SpikeTA. The artificial triggers created around each of the three triggers shown in are illustrated in . Segments of rectified EMG beginning 30 ms before and ending 50 ms after each of these 81 artificial triggers then were aligned and averaged to create a single-spike increment shifted average (ISA), as illustrated in . The time sampled by a single-spike ISA thus extended from 30 ms before the first of the 81 artificial triggers to 50 ms after the last. The 1 kHz periodicity in the single-spike ISAs reflects the 1 ms increment used to generate these moving averages.
This approach distributes any spike-locked change in EMG evenly over the single-spike ISA, while capturing any baseline trend on the scale of 80 ms or more. The three single-spike ISAs shown in each capture the general trend in EDC EMG activity at the time each of the three spikes was discharged by the M1 trigger neuron. The first of these three spikes was discharged when EMG was relatively low but increasing; the second when EMG was higher and peaking, and the third when EMG was decreasing.
All 6,014 single-spike ISAs associated with the 6,014 spikes from the trigger neuron then were averaged together, resulting in an overall ISA estimate of the baseline trend, shown as the center trace in . This overall ISA was subtracted from the initial SpikeTA ( left), and the value of the initial SpikeTA at the trigger time was added back to restore the baseline level above zero (not shown here), resulting in the ISA adjusted SpikeTA shown on the right of . This procedure eliminated the curvilinear rising baseline trend from the adjusted SpikeTA.
Additional examples of overall ISA baseline estimates and adjusted SpikeTAs are shown in the middle and right columns, respectively, of . Note that whereas the initial SpikeTAs in (left) are shown for 160 ms to emphasize the baseline trends, the ISA baselines (center) and adjusted SpikeTAs (right) are 80 ms in duration, the duration of the SpikeTAs used for measurement of SpikeTA effects. Each ISA captures the baseline trend of each initial SpikeTA, while eliminating the brief peak or trough of the SpikeTA effect per se. In addition, because the ISA averages 81 EMG samples for each EMG sample used in the initial SpikeTA, the noise level of the ISA is markedly lower, and therefore the ISA-adjusted SpikeTA preserves the noise features of the initial SpikeTA.
For SpikeTAs in which the baseline trend is highly linear (e.g. ) the present ISA method has no advantage over linear ramp subtraction. As the baseline trend comes to include more curvature (e.g. Figure D, E and F), however, subtraction of a linear ramp fails to remove curvature from the baseline trend preceding the SpikeTA effect, as well as from the baseline trend underlying the SpikeTA effect (peak or trough) itself. Such curvature contributes extra variation to the baseline preceding the trigger time, producing inaccurate measurements of onset and offset of the SpikeTA effect. Furthermore, baseline trend curvature at the time of the SpikeTA effect itself results in over- or under-estimation of the amplitude of the SpikeTA effect. Using the ISA to adjust the baseline of the SpikeTA results in more accurate measurement of both the amplitude and temporal characteristics of the SpikeTA effect. ISA adjustment therefore has been applied to all SpikeTAs described below.
Identifying changes in SpikeTA effects with a Multi-Fragment Approach
shows two ISA-adjusted SpikeTAs compiled for the same neuron-muscle pair (e0030-FDPu) during the squeeze task (left) and during RPD of the neuron versus the muscle (right), both displayed on the same vertical scale. Both SpikeTAs show a clear post-spike facilitation (PSpF). Whereas the baseline levels of average EMG are similar, the PSpF on the right appears to be larger than that on the left, and appears longer in duration as well. Indeed, during the squeeze task, the PSpF had a measured PPI of 12.7% and a PWHM of 6.1 ms, whereas during RPD the PSpF had a PPI of 28.5% and a PWHM of 9.5 ms. Given that the squeeze SpikeTA was compiled from 16,109 EMG-filtered triggers, and the RPD SpikeTA from 24,470, one would expect these differences in PPI and PWHM to be genuine. But are these differences significant, or might they have arisen by chance alone?
Figure 4 A. SpikeTAs compiled for the same neuron-muscle pair during two different behavioral epochs are shown on the same vertical scale (Volts after amplification 10,000x). Left, Epoch 1 – squeeze task; Right, Epoch 5 – RPD of neuron e0030 versus (more ...)
We previously have used a multi-fragment approach to identify significant peaks or troughs in SpikeTAs (Poliakov and Schieber, 1998
). The spike train used to compile the SpikeTA is divided into multiple non-overlapping fragments, each including the same number of spikes, and a separate SpikeTA is compiled using the spikes in each fragment, resulting in multiple fragment spike triggered averages (fSpikeTAs
). From each fSpikeTA
, a test statistic, d
, is calculated as the mean value in a temporally fixed test window (chosen to encompass the typical peak or trough times of SpikeTA effects), minus the mean in two control windows that immediately precede and follow the test window. The test statistic, d
, thus assays whether the
in the test window is systematically different from that in the flanking control windows. If d
is statistically different from 0 across all the fragments, then we reject the null hypothesis that no peak or trough exists in the test window, and accept the peak or trough in the overall SpikeTA as significant.
shows 15 sequential ISA-adjusted fSpikeTAs, each compiled with 100 sequential EMG-filtered spikes during the squeeze task (left) and RPD (right). Each fSpikeTA is scaled to fill the same vertical height. Peaks in the fSpikeTAs corresponding in time to the peak in the overall SpikeTAs shown above in A are more evident in the RPD fSpikeTAs (right). Even though peaks are not obvious in the squeeze task fSpikeTAs (left), the same phenomena underlying the obvious peak in the overall SpikeTA can be assumed to be present, on average, in these fSpikeTAs.
(left) shows values of d (test window 6 to 16 ms after the trigger) for all 161 fSpikeTAs from the squeeze task (left). While not apparent upon visual inspection of the 15 fSpikeTAs shown above in shows that on average d was greater than 0. Furthermore, although the values of d from successive fSpikeTAs might be considered a time series, and therefore not independent of one another, shows that the variability of d was greater than any temporal trend. Indeed, the correlation between each value of d and its predecessor gave a correlation coefficient of r = 0.11 for the squeeze task and r = 0.12 for the RPD task, providing little suggestion that the successive d values depend on one another within either epoch.
Because the distribution of d often was significantly different from normal (Lilliefors test, P < 0.05), we used the Wilcoxon signed rank test to determine whether the median value of d was significantly different from 0. Applying this test to the 161 d values from squeeze task fSpikeTAs (, left), and to the 244 d values from the RPD task (, right) confirmed that the PSpF peaks in both overall SpikeTAs () were significant (squeeze task P = 2.3x10−23, RPD task P = 1.2x10−41).
We now extend this approach to compare measured amplitude and temporal features of SpikeTA effects obtained during different behaviors. The peak percent increase (PPI) of a SpikeTA effect is a standard measure of the amplitude of the peak (or trough) normalized as a percentage of the baseline EMG preceding the trigger time in the SpikeTA. As for d, even when no peak or trough is evident above the noise in a given fSpikeTA, we can measure a peak percent increase for the fragment. First, the mean value of the fSpikeTA in the baseline sample period (from 30 to 10 ms before the trigger time, see ) is subtracted from the mean value in the test window 6 to 16 ms after the trigger. If the result is positive (negative), a peak (trough) is assumed to be present. The fSpikeTA baseline mean is subtracted from the maximum (minimum) in the test window, providing a measure of the peak amplitude for the fragment. This fragment peak amplitude then is divided by the baseline mean and multiplied by 100, normalizing the peak amplitude as a percentage of the baseline mean, the fragment PPI, fPPI. This fPPI then can be used as a test statistic to compare PPIs of overall SpikeTAs compiled during different behavioral epochs.
The fPPIs of all fragment SpikeTAs from neuron-muscle pair e0030-FDPu are shown in as a function of time. Neuron-muscle pair e0030-FDPu was recorded during nine behavioral epochs. (Data from epochs 1 and 5 were shown in .) The overall ISA-adjusted SpikeTA compiled for each epoch is shown above in . To facilitate comparison of PPI in the nine epochs, here the vertical scale of each of the nine SpikeTAs has been adjusted separately to match the heights of the baselines above zero. Note that the SpikeTA effect of this neuron-muscle pair was particularly strong during epochs 3 and 5, PPIs for these two overall SpikeTAs being 30.0% and 28.5%, respectively.
Figure 5 A shows SpikeTAs for neuron-muscle pair e0030-FDPu compiled for each of nine successive behavioral epochs: 1) Squeeze task; 2) RPD of neuron e0030 versus FCR; 3) e0030 v FDPu; 4) e0030 v FCR; 5) e0030 v FDPu; 6) e0030 v APL; 7) e0030 v ECRB; 8) e0030 (more ...)
As for d, the variability of fPPI within a given epoch was substantially greater than any temporal trend within that epoch. This fragment-to-fragment variability primarily reflects variability in the behavior sampled during the relatively brief periods represented by individual fSpikeTAs, in the present example averaging 2.8 sec (range 1.5 to 33.4 sec). In spite of this variability, differences between certain epochs were evident. For example, fPPI values were systematically greater in epochs 3 and 5 than in epoch 1.
Because fPPIs in individual epochs were not normally distributed, we used the Kruskal-Wallis test of the null hypothesis that PPI did not vary among epochs. This hypothesis was rejected (P = 0). Post-hoc pairwise comparisons also confirmed that PPI in epoch 1 differed from that in epochs 3 and 5 (Wilcoxon rank sum tests, P = 1.4x10−30 and P = 1.1x10−31). We therefore can conclude that the values of PPI measured in the overall SpikeTAs of epoch 1 (PPI = 12.7%) indeed was significantly different from that measured in epoch 3 or 5 (PPI = 30.0% and 28.5%, respectively). In contrast, post-hoc testing showed no significant differences between the PPIs of epochs 3 and 5, or among PPIs of epochs 2, 4, 6, 7, 8 and 9 (P > 0.0015, i.e. P > 0.05 after Bonferroni correction for 36 pairwise comparisons).
A similar multi-fragment approach can be used to evaluate changes in a measure of the temporal width of the SpikeTA effect, the peak width at half maximum (PWHM). PWHM is a useful measure of width because its measurement is affected less by the residual noise level of the SpikeTA than is the total duration of the peak or trough from onset to offset. Again, although a peak or trough may not be evident in each fragment SpikeTA, the same algorithm used to measure the PWHM of an overall SpikeTA can be applied to each fSpikeTA (see ) The level halfway from the baseline period mean to the maximum (minimum) in the test window is determined (half the peak amplitude). The fSpikeTA then is followed both backward and forward from the maximum (minimum) to find the times at which the fSpikeTA crosses the half-maximum level before and after the maximum (minimum), and the difference between these two times provides a measure of the PWHM for the fragment, fPWHM, which can be used as a test statistic for comparing PWHMs. shows the fPWHMs for e0030-FDPu as a function of time. As for fPPI, categorizing these data according to the nine behavioral epochs and applying a Kruskal-Wallis test confirmed that PWHM varied across epochs (P = 0). Post-hoc testing confirmed that the PWHMs of 9.0 ms in epoch 3 and 9.5 ms in epoch 5 differed significantly from the PPI of 6.1 ms in epoch 1 (Wilcoxon rank sum tests, P = 3.2x10−26 and 6.9x10−25).
Additional SpikeTA effect measures—mean percent increase (MPI), onset time and offset time—also can be compared using the multi-fragment approach, subject to one caveat. Our algorithm for assessing these features requires that a fSpikeTA exceed 2 standard deviations (SDs) from its baseline mean during the test window. Then the fSpikeTA can be traced both backward and forward to its 2 SDs level to find the onset and offset times (fOnset and fOffset)of the fSpikeTA, respectively (see ). An MPI for the fragment (fMPI) then is computed between the fOnset and fOffset. In some fragments the fSpikeTA may not exceed 2 SDs during the test window. For such fragments, fOnset and fOffset, and consequently fMPI, cannot be determined. Nevertheless, the fragments in which the SpikeTA does exceed 2 SDs from the baseline mean still can be used to measure fOnset, fOffset, and fMPI. shows fOnset, fOffset, and fMPI for all measurable fragments from e0030-FDPu as functions of time. Each of these test statistics varied significantly across the nine behavioral epochs, permitting us to reject the null hypotheses that the Onset, Offset and MPI of the overall SpikeTA effects shown above in did not vary among epochs (Kruskal-Wallis tests: Onset P = 6.0x10−15; Offset P = 0; MPI, P = 0). Comparing epoch 1 with epoch 5 specifically, post-hoc testing confirmed significant differences in MPI (7.1% vs 12.7%, P = 2.8x10−21) and Offset time (17.4 ms vs 29.4 ms, P = 1.4x10−13), but not in Onset time (5.4 ms vs 4.6 ms, P > 0.0015)(Wilcoxon rank sum tests).
Values of fragment test statistics correlate with values measured in overall SpikeTAs
While fPPI, fPWHM, fMPI, fOnset, and fOffset are used as test statistics, the average values of these test statistics do not equal the values of PPI, PWHM, MPI, Onset and Offset, respectively, measured from the overall SpikeTAs. Differences result from the fact that the residual noise level of each fSpikeTA will be greater than that of the overall SpikeTA. As illustrated in , the higher noise level will result in a higher peak (lower trough) value of the fSpikeTA maximum (minimum), resulting in a tendency for fPPIs to be systematically greater than the corresponding PPI measured for the single overall SpikeTA. Furthermore, because the higher noise level of the fSpikeTAs moves the ± 2 standard deviation levels further from the baseline period mean, fOnsets will be found systematically later, and fOffsets systematically earlier, (both closer to the time of the maximum (minimum)) than the Onset and Offset times for the overall SpikeTA. And because the fOnsets and fOffsets will be closer to the time of the maximum (minimum), fMPIs will tend to be larger than the MPI measured in the overall SpikeTA. Similarly, the times at which the SpikeTA crosses the half maximum level before and after the maximum (minimum) both will be closer to the maximum (minimum), and consequently fPWHMs will tend to be shorter than PWHM measured in the overall SpikeTA.
Nevertheless, if the fragment test statistics reflect the same underlying process revealed in the overall SpikeTA, the test statistic values averaged across all the fragments of a given epoch should be correlated with the value measured from the overall SpikeTA. shows scatter plots of average fPPI plotted against the PPI of the overall SpikeTA effect, and of average fMPI plotted against the overall MPI, for neuron-muscle pair e0030-FDPu in the 9 epochs illustrated in . shows similar plots for PWHM, Onset and Offset. For all five of these parameters, the average of values measured from fSpikeTAs was correlated with the value measured from the overall SpikeTA (correlation coefficients: PPI, 0.97; MPI, 0.91; PWHM, 0.86; Onset, 0.80; Offset, 0.86), and linear regression showed a slope significantly different from 0 (P < 0.01, or P < 0.05 after Bonferonni correction for 5 tests). These observations support use of the fragment test statistics for comparison of measurements made on the overall SpikeTA effects.
Figure 6 Correlations between test statistics and measures of overall SpikeTA effects. Scatterplots illustrate the correlations between the average of fragment test statistics (ordinate: fPPI, fMPI, fPWHM, fOnset and fOffset) and measurements of amplitude (A) (more ...)
Selecting the number of spikes per fragment
The multi-fragment approach entails a trade-off between the number of spikes per fragment and the number of fragments. The more spikes allotted to each fragment, the better the peak-to-noise ratio of each fSpikeTA
, but the fewer the number of fragments. In the limit of having too few fragments (samples), statistical testing may not show significance, even in the presence of clear SpikeTA effects. As a reasonable compromise between these two extremes, for a recording that included a total number, N
, of EMG-sweep filtered spikes, we previously chose to use
fragments, each with
spikes to test an individual SpikeTA for the presence of a SpikeTA effect (Poliakov and Schieber, 1998
In the present study, however, we endeavor to compare SpikeTAs compiled for the same neuron-muscle pair during different epochs, a condition under which the number of EMG-filtered triggers available typically differs among epochs. Using
then would result in different numbers of spikes per fragment in different epochs. The peak-to-noise ratio of the fSpikeTAs
would vary among epochs simply because different numbers of spikes per fragment would be used in analyzing different epochs. We therefore chose here to use a fixed number of spikes per fragment. Because our epochs typically contained ~10,000 spikes, we chose to fix the number of spikes per fragment at 100 (close to
). In this way, although the number of fragments might vary among epochs, the peak-to-noise ratio of the fSpikeTAs
remains unaffected by differences in the number of spikes in various epochs.
To examine the effect of the number of spikes per fragment (hereafter abbreviated as n) and of the number of fragments (ν = N/n) on the outcome of multi-fragment analysis, we repeated the analysis using different n. The results are illustrated in for a neuron-muscle pair (e0030-Thenar) which showed clear SpikeTA effects in some behavioral epochs, but no SpikeTA effect in others (). shows how the P-value obtained for the presence of a significant peak (test-statistic d, Wilcoxon signed rank test) varied when we used n = 25, 50, 100, 200, 400, 800, 1600, or 3200 spikes per fragment. As expected, P-values generally were closer to 1 for larger n () and smaller ν (), and fell toward 0 as n was reduced increasing ν. Note, however, that for epochs 1 and 7 the SpikeTAs of which have no substantial peak, the P-value never fell below 0.1. Neither of these epochs would have been considered to have a significant SpikeTA effect using any value of n, even using a small n to produce a large ν. Epoch 9 provides an instructive borderline example. Visual inspection of the epoch 9 SpikeTA shown in suggested a questionable peak at the time a SpikeTA effect would be expected (*), but this peak was not substantially larger than features at other times in the same SpikeTA (^). As n was Page 27 reduced from 3200 to 25, increasing ν from 2 to 281, the P-value for epoch 9 fell from levels that could not be considered significant, through levels typically used as criteria for significance. With a criterion set at P < 0.005 (0.05 after Bonferonni correction for the 9 tests on 9 epochs), the SpikeTA of epoch 9 would not have been considered significant using n >100, but would have been considered significant using n≤ 100 (ν≥ 89). Using the same criterion (P < 0.005), the SpikeTA effects of the remaining epochs (2, 3, 4, 5, 6, 8) all would be considered significant using n ≤ 400 (ν≥ 41). Had these epochs with more obvious SpikeTA effects not been available, one would doubt the presence of an effect in epoch 9, but in view of the evident effect for this neuron-muscle pair in several epochs, the effect in epoch 9 might be considered marginally genuine.
Figure 7 Effect of the number of spikes per fragment and of the number of fragments on P-values. A. SpikeTAs from nine behavioral epochs are shown for a neuron-muscle pair selected for the presence of a clear SpikeTA effect in several epochs (2, 3, 4, 5, 6, and (more ...)
shows how the P-value obtained for variation among the nine epochs in measured parameters PPI, MPI, PWHM, Onset and Offset (test statistics fPPI, fMPI, PWHM, fOnset and fOffset; Kruskal-Wallis tests) changed when we used different n, and consequently different ν, respectively. With a criterion of P < 0.01 (P < 0.05 after Bonferonni correction for 5 tests), neither Onset nor Offset would have been considered to vary significantly among epochs using any n. With the same criterion, PWHM, MPI and PPI all would have been considered to vary significantly using any n < 3200 (average ν > 10).
Using a small n and a large ν thus will result in more SpikeTA effects, and more changes among epochs, being accepted as significant for any fixed criterion P-value. Using a large n and a small ν may result in some genuine SpikeTA effects and some real variation among epochs being rejected. For our present data, n = 100 appeared to be a reasonable compromise between these two extremes.