Computer Simulation of Post-Spike Facilitation in Spike-Triggered Averages of Rectified EMG

S. N. Baker and R. N. Lemon

Sobell Department of Neurophysiology, Institute of Neurology, Queen Square, London WC1N 3BG, United Kingdom

    ABSTRACT
Abstract
Introduction
Methods
Results
Discussion
References

Baker, S. N. and R. N. Lemon. Computer simulation of post-spike facilitation in spike-triggered averages of rectified EMG. J. Neurophysiol. 80: 1391-1406, 1998. When the spikes of a motor cortical cell are used to compile a spike-triggered average (STA) of rectified electromyographic (EMG) activity, a post-spike facilitation (PSF) is sometimes seen. This is generally thought to be indicative of direct corticomotoneuronal (CM) connections. However, it has been claimed that a PSF could be caused by synchronization between CM and non-CM cells. This study investigates the generation of PSF using a computer model. A population of cortical cells was simulated, some of which made CM connections to a pool of 103 motoneurons. Motoneurons were simulated using a biophysically realistic model. A subpopulation of the cortical cells was synchronized together. After a motoneuron discharge, a motor unit action potential was generated; these were summed to produce an EMG output. Realistic values were used for the corticospinal and peripheral nerve conduction velocity distribution, for slowing of impulse conduction in CM terminal axons, and for the amount of cortical synchrony. STA of the rectified EMG from all cortical neurons showed PSF; however, these were qualitatively different for CM versus non-CM cells. Using an epoch analysis to determine reliability in a quantitative manner, it was shown that the onset latency of PSF did not distinguish the two classes of cells after 10,000 spikes because of high noise in the averages. The time of the PSF peak and the peak width at half-maximum (PWHM) could separate CM from synchrony effects. However, only PWHM was robust against changes in motor unit action-potential shape and duration and against changes in the width of cortical synchrony. The amplitude of PSF from a CM cell could be doubled by the presence of synchrony. It is proposed that, if a PSF has PWHM <7 ms, this reliably indicates that the trigger is a CM cell projecting to the muscle whose EMG is averaged. In an analysis of experimental data where macaque motor cortical cells facilitated hand and forearm muscle EMG, 74% of PSFs fulfilled this criterion. The PWHM criterion could be applied to other STA studies in which it is important to exclude the effects of synchrony.

    INTRODUCTION
Abstract
Introduction
Methods
Results
Discussion
References

Experiments in which spontaneous neural activity is recorded from awake, behaving animals have added considerably to our understanding of the function of the nervous system. In such studies, it is always advantageous to identify the inputs and outputs of recorded cells so that the functional role of their firing pattern can be interpreted. Identification of the output is of particular importance because it reveals the target structure of information generated by a particular neuron. In the motor cortex, for example, neurons whose axons descend as far as the pyramidal tract can be antidromically activated from electrodes chronically implanted in the tract (e.g., Evarts 1964).

The technique of spike-triggered averaging (STA) was proposed to allow further identification of pyramidal tract neurons (PTNs) (Fetz and Cheney 1980). This approach reveals any facilitation or suppression of muscle activity that is exerted by the PTN. An electromyogram (EMG) is recorded from a number of muscles simultaneously with the cell discharge. The rectified EMG is then averaged with respect to the neuron spikes. A brief, short-latency post-spike facilitation (PSF) in the average is conventionally taken as indicating that the trigger cell makes monosynaptic connections to the motoneurons innervating the muscle whose EMG was averaged. STA thus allows both identification of the output neuron and its target muscles.

However, whereas antidromic activation is an uncontroversial technique, the use of STA in this way is not without its difficulties. It has been repeatedly pointed out (Fetz and Cheney 1980; Kirkwood 1994) that a PSF in a STA does not necessarily indicate that a cell connects to motoneurons but rather could occur because the trigger cell is synchronized in firing with corticomotoneuronal (CM) cells. The presence of such synchrony is well established (Baker et al. 1997; Smith and Fetz 1989). This issue is addressed here using a computational model of the production of PSF. It is shown that a non-CM cell can produce a significant PSF but that this is qualitatively and quantitatively different from PSF produced from a CM cell. Means of distinguishing between synchrony PSFs and those caused by direct connections from trigger cell to motoneurons are investigated and compared.

    METHODS
Abstract
Introduction
Methods
Results
Discussion
References

A diagram of the cortical and spinal network modeled in this study is given in Fig. 1. A variety of cortical cell types had input to the simulated motoneuron pool. From the discharge pattern of the motoneuron pool, a simulated EMG was produced. This permitted study of the transform from cortical input to EMG output and provided a direct and realistic comparison with experimental data. The following describes the details of the model.


View larger version (26K):
[in this window]
[in a new window]
 
FIG. 1. Overview of the model. Cortical cells (black-square) are modeled after Ashby and Zilm (1982). A common input cell connects to a population of target cells, synchronizing them. Most of these cells are corticomotoneuronal (SCM) and make monosynaptic connections to the motoneurons (black-square). The S cell is synchronized with the SCM cells, but has no connections to the motoneurons. The CM cell is corticomotoneuronal but is not synchronized with the other cortical cells. Its central conduction time from cortex to spinal cord is 1.4 ms. The motoneurons innervate motor units, which generate realistic motor unit action potentials (MUAPs). These are summed to produce the simulated electromyogram (EMG). The peripheral conduction time varies linearly with motor unit recruitment number.

Neurons

The neurons were simulated by using two complementary models published previously. The cortical neurons (Fig. 1, black-square) were based on the description by Ashby and Zilm (1982). The membrane potential rose linearly at a constant rate of K for all the cells until it crossed the firing threshold (taken arbitrarily to be 0 mV). The cell was then assumed to have produced an action potential, and its potential was reset to a subthreshold value. A value for K of 82.5 µV/ms was used, chosen so that a cell which had an afterhyperpolarization (AHP) of 10 mV below threshold would fire at a frequency of 8.25 Hz. The precise value of the AHP was determined from
AHP = <IT>KI</IT> (1)
where I is the next interspike interval, chosen at random after each threshold crossing from a gamma distribution with a shape parameter of 4 and a mean equal to the reciprocal of the desired average firing rate. Inputs to the cells produced excitatory postsynaptic potentials (EPSPs) with a sigmoidal rising phase of variable duration and an exponentially decaying falling phase with a time constant of 4.8 ms. These summed linearly with the cell membrane potential and were therefore able on some occasions to cause a threshold crossing to occur slightly earlier than it would have done without the input. The advantage of this model neuron is that it is computationally simple so that simulation time is reduced. It is an obvious simplification of the complex behavior of a cortical neuron but suffices to produce a realistic spike train having the desired stochastic properties.

The motoneurons (Fig. 1, bullet ) were simulated with the use of a published model (model 9-2; see Powers 1993). Briefly, this is a single-compartment biophysical model, having a fast and slow potassium and a low- and high-threshold calcium voltage-gated conductance, whose combined action produces a realistic AHP after cell firing. The membrane time constant was 6 ms. Action potentials were not simulated explicitly; instead, whenever the membrane potential crossed threshold, it was clamped to a depolarization of +20 mV for 1 ms and then reset to the threshold level. After this, there was a 1-ms refractory period. The threshold for spike initiation was not fixed but varied with a fast and slow time course. This is a far more realistic model of the motoneuron than that according to Ashby and Zilm (1982), but considerably more computationally intensive to simulate.

The model of Powers (1993) was extended for the current study as follows. Inputs to the cell produced activation of a synaptic conductance having a reversal potential of 0 mV with an alpha  function-shaped time course given by
<IT>g<SUB>s</SUB></IT>(<IT>t</IT>) = (<IT>g</IT><SUB>max</SUB><IT>te</IT><SUP>1−<IT>t</IT>/τ</SUP>)/τ (2)
where gmax = 0.015 µS and tau  = 0.2 ms. These parameters produced an EPSP with a height of 70 µV and rise time of 1 ms. Synaptic noise was simulated in the motoneurons by adding Gaussian noise to the simulated membrane potential with an SD of 2 mV and time constant of 4 ms, in accordance with experimental data from Calvin and Stevens (1968). The cells were made to fire tonically by activating the synaptic conductance by a constant amount independent of the activity of simulated inputs. The level of this constant shunt conductance was chosen to produce the desired firing rate in preliminary, calibration runs, with no cortical inputs to the motoneurons. One advantage of this model over more simplistic treatments is that the membrane potential rise after a spike is not quite linear but has a decreasing slope as the threshold is approached. This is in accord with findings in human motoneurons by Olivier et al. (1995) and Matthews (1996) and may bias the cells to respond somewhat preferentially to synchronous inputs (Matthews 1996).

Motoneuron pool

The mean firing rate for each motoneuron and the number of motoneurons simulated were determined from the model of a tonically firing motoneuron pool proposed by Wani and Guha (1975), based on the results of Milner-Brown et al. (1973a,b) regarding motoneuron recruitment order and firing rate. Their model was solved numerically for the first dorsal interosseus (1DI) muscle contracting at 2% of its maximum voluntary contraction (a strength typical of the activation of this muscle during many skilled tasks). This indicated that 103 of the total of 377 motoneurons would be tonically firing. The first recruited motor unit would be firing at 8.5 Hz, and the last would be firing at 8.0 Hz.

Surface EMG

After the discharge of each motoneuron in the model, a motor unit action potential (MUAP) was generated. The form of this potential was intended to approximate closely the potential which would be recorded by a surface EMG electrode. The potentials were based on data gathered in two human subjects specifically for this investigation.

Single motor unit discharges were recorded from the 1DI muscle using intramuscular needle electrodes in two human subjects. The subject was given auditory and visual feedback of the rate of discharge of the unit and was required to maintain firing at a low rate (~10 Hz) by abducting the index finger against a force transducer. EMG was recorded differentially from two adhesive surface electrodes (H59P, ARBO Medical), one over the knuckle of the index finger and one over the belly of the muscle, and 30-Hz, high-pass filtered. The output of the force transducer together with the EMG and unit discharge were recorded on magnetic tape for subsequent analysis.

Off-line, the single unit action potentials were converted to standard transistor-transistor logic (TTL) pulses by a double time-amplitude window discriminator. A computer digitized the time of occurrence of these pulses and sampled the surface EMG (5-kHz sampling rate) and the force (100 Hz). Accepted potentials were monitored during data capture to ensure that only potentials from one unit triggered the discriminator; additionally, the interval histogram was monitored for the presence of unreasonably short intervals.

Averages of the surface EMG were then compiled with the TTL accept pulses as triggers. These motor unit averages (Lemon et al. 1990; Yemm 1977) provide a representation of the MUAP as it is recorded by the surface electrodes. To determine the twitch tension of the unit, averages of the force record were also produced with only those spikes where the succeeding interspike interval was >100 ms. This spike selection reduced contamination from the effects of succeeding spikes, which could artificially increase the measurement of twitch tension (see Calancie and Bawa 1986; Nordstrom et al. 1989). The model of Wani and Guha (1975) gives expressions for the twitch tension of 1DI motor units as a function of unit recruitment order. It was therefore possible to determine the approximate recruitment order of each unit from its twitch tension.

A total of 12 units from subject NA and 20 from subject MM were clearly discriminable and yielded sufficient triggers to produce a high signal-to-noise ratio in the averages; Fig. 2 shows examples of MUAPs for some of these units. A priori, it might be expected that the amplitude of the surface MUAP would be positively correlated with the unit's twitch tension (Lemon et al. 1990) and that the waveform shape should change smoothly as the motoneuron recruitment number increased. An interpolation algorithm was therefore used to produce a surface potential for each of the 103 motor units in the model.


View larger version (24K):
[in this window]
[in a new window]
 
FIG. 2. Examples of the MUAPs used in this study. On the left are examples of actual potentials recorded from 2 human subjects and extracted from the surface EMG by motor unit triggered averaging. Motor unit recruitment number is shown beside the traces, estimated from the unit twitch tension using the model of Wani and Guha (1975). up-arrow , time of the intramuscular action potential used to trigger the averages. On the right are examples of interpolated MUAPs, generated by the algorithm described in the text. These traces are plotted starting at the potential onset. A: subject NA; B: subject MM. Note the smooth change in size and shape of the interpolated potentials in A not evident in B.

For each subject, the surface MUAPs were aligned by eye at the first point where there was a clear deviation in the signal from the background. A section of the average following this point was used, with duration chosen to be sure that even late components of the surface potential were included. This was 35 ms for subject NA and 30 ms for MM. At the sampling rate of 5 kHz used, this yielded a set of N values of potential for each unit (N = 150 or 175). These were treated as vectors sj and subjected to principle component analysis with the aim of describing the variation in shape of the MUAPs by a small number of parameters (Manly 1986). Each unit's MUAP was expressed in terms of the mean µ and principle components pi to yield the coefficients cji.
<B>s</B><SUP><IT>j</IT></SUP>= <B>μ</B>+ <LIM><OP>∑</OP><LL><SUB><IT>i=</IT>1</SUB></LL><UL><SUP><IT>i</IT>=<IT>N</IT></SUP></UL></LIM><IT>c</IT><SUP><IT>j</IT></SUP><SUB><IT>i</IT></SUB><B>p</B><SUB><IT>i</IT></SUB> (3)
Here j is an index identifying a particular experimentally recorded MUAP, and i indexes the principle components. The equality of Eq. 3 is exact if the sum is taken over all n principle components pi; however, it is a feature of principle component analysis that use of only the first few terms of the summation yields a good approximation to sj. The coefficients cji then provide a representation of sj in a small number of parameters.

Plots of ci for i = 1-4 were then made versus the motoneuron's recruitment number n. For subject NA, c1, c2, and c4 showed a deterministic relationship with n; this was expected, as the amplitude of the individual MUAPs showed a high correlation with n in this subject (r = 0.73, P < 0.01). These parameters were therefore fitted by straight lines (c1 and c2) or a quartic (c4) as judged appropriate by eye. The line c3 did not show any systematic relationship with n. A set of 103 surface unit potentials were then generated from the first four principle components, allowing coefficients c1, c2, and c4 to vary in accordance with the fits and fixing c3 to make the following relationship true
[<B>μ</B>] + <LIM><OP>∑</OP><LL><SUB><IT>i</IT>=1</SUB></LL><UL><SUP><IT>i</IT>=4</SUP></UL></LIM><IT>c<SUB>i</SUB></IT>[<B>p</B><SUB><IT>i</IT></SUB>] = 0 (4)
where the bracket notation [x] has been used to indicate the mean of the points in vector x. This constrains the average voltage level of each interpolated waveform to be zero, so ensuring that there are no small differences between the negative and positive phases, which on average would lead to a DC offset in the final simulated EMG.

For subject MM, there was no significant correlation between MUAP amplitude and recruitment number, and no deterministic relationships between the ci and n were seen. In this case, the mean and SD of each ci were determined, and c1-c3 selected from a normal distribution having these parameters. The parameter c4 was adjusted according to Eq. 4. The failure to find a deterministic relationship for the ci in this subject implies that the surface MUAP amplitude is uncorrelated with recruitment order, contrary to the findings of Lemon et al. (1990). Thomas et al. (1987) however showed a similar lack of correlation, and Lemon et al. (1990) noted that the extent of correlation varied between subjects. Examples of the interpolated MUAPs from each subject's data are shown in Fig. 2.

It should be noted that, although the interpolation method used attempted to reproduce the pattern of MUAP variation with recruitment order seen in the recordings, the recruitment order of the units had little influence on the outcome of the simulations. As stated above, the tonic firing rates of the simulated units only varied from 8 to 8.5 Hz between the lowest and highest threshold unit. The motor units were therefore almost identical, apart from differences in their MUAPs. Details of the correlation between MUAP amplitude and unit recruitment order will not therefore critically influence the output of the model.

In the simulation, whenever a motoneuron fired it caused its MUAP to occur after a delay tau , given by
τ = <IT>D</IT>/(55 + 30<IT>N</IT>/377) (5)
where D is the peripheral conduction distance, estimated as 0.5 m for a hand muscle in the adult macaque monkey. Equation 5 allows for the small amount of peripheral dispersion that occurs and makes the simplifying assumption that conduction velocity for the 377 motoneuron axons is perfectly correlated with recruitment number n (Burke 1981) over the range of 55-85 m/s (Lemon et al. 1986). The range of conduction times tau  obtained for the pool of 103 motoneurons used here was 8.0-9.2 ms.

Network

The three main classes of cortical cell whose firing was used to produce STAs of the simulated rectified EMG are shown in Fig. 1. One nonsynchronized CM cell made monosynaptic connections to all the motoneurons of the pool. It had a spike train that was independent of the activity of other cortical cells. The synchronized CM (SCM) cells also made monosynaptic connections to the motoneurons, but in contrast their discharges were synchronized with each other because they received synaptic input from a "common input" cell firing tonically at 40 Hz. This input produced an EPSP in the SCM cells. In initial simulations, the size and rise time of this EPSP was adjusted so that the cross-correlation peak between two of the SCM cells would have a peak duration of 15 ms and a peak strength A = 0.06, where
<IT>A</IT>= 2<IT>P</IT>/(<IT>n + n</IT>) (6)
and where P is the area of the cross-correlation peak above baseline, nT is the number of cross-correlation trigger spikes, and nr is the number of cross-correlation response spikes.

These figures for cross-correlation peak width and height are in accordance with the mean values found when pairs of identified CM cells in the monkey motor cortex are cross-correlated (Baker 1995). Subsequent simulations investigated the effect of varying the width of the cortical synchrony peak (see RESULTS). It is obviously unlikely that this observed cortical synchrony is produced by a single common input cell; rather it is the result of common input from a population of cells. This is reflected in the slow rise time of the EPSP from the model common input cell. However, for the purposes of the modeling described here, it is only important that the SCM cells in the model show synchronous firing which is quantitatively similar to that seen in experimental data; the exact means by which this is achieved is unimportant. The number of SCM cells present in the network could be varied, thereby altering the amount of synchrony in the total input to the motoneuron pool.

A third class of cell was one synchrony only (S) cell shown in Fig. 1. This made no monosynaptic connections to the motoneurons, but received the common input that also projected to the SCM cells. It was thus synchronized with the SCM cells but could exert no causal effects on the motoneuron pool itself.

The S, CM, and SCM cells all fired with a tonic firing rate of 10 Hz. This is lower than the maximum burst rate of cells recorded in the monkey primary motor cortex. However, the results will not be changed greatly by altering this rate since the measure used to quantify synchrony between cortical cells (A; Eq. 6) is independent of firing rate.

Corticospinal conduction delay

The pyramidal tract contains fibers of widely varying diameters, implying a range of conduction velocities from 5 to 85 m/s (Humphrey and Corrie 1978; Mediratta and Nicoll 1983; Tan et al. 1979). Because monosynaptic connections from the cortex are only present to any great extent in primates, it has been assumed that it must be only the very large corticospinal fibers (another primate feature) which make such connections. However, more recent research has shown that, although CM connections are less frequently seen from slow PTNs, they do exist (Fetz and Cheney 1980; Lemon et al. 1986; Porter and Lemon 1993).

Inclusion of PTNs with different conduction velocities would not be expected to alter the form of the CM cell PSFs in the model described here, other than to shift their latency to remain appropriate for the conduction velocity of the particular cell under consideration. However, if the SCM cell population has a variable conduction velocity, then the PSFs from both SCM and S cells will be affected, as the synchrony of spike activity in the corticospinal terminals just before making contact with the motoneurons will be more dispersed than at the cortical level. It was therefore considered important to include a realistic corticospinal conduction velocity distribution in the model. An analysis concerning the effect of conduction velocity dispersion on cross-correlation peaks in the bulbospinal system has been carried out by Davies et al. (1985).

To produce a realistic distribution of conduction velocities for the SCM cell population, the antidromic latencies of 132 identified CM cells which responded antidromically to a PT stimulus in three different Macaca nemestrina monkeys were used (R. N. Lemon, K.M.B. Bennett, and D. Flament, unpublished observations). Assuming a conduction distance of 50 mm from PT to cortex (Humphrey and Corrie 1978) and a utilization time of 0.1 ms (Lemon 1984), this allowed a conduction velocity to be calculated for each cell. A cumulative distribution graph for these data is plotted in Fig. 3 ().


View larger version (16K):
[in this window]
[in a new window]
 
FIG. 3. Cumulative probability curve showing the distribution of conduction velocity estimated from antidromic responses of 132 identified corticomotoneuronal (CM) cells to pyramidal tract stimulation in 3 monkeys. , observed distribution; - - -, a correction after Humphrey and Corrie (1978) to remove the effects of electrode sampling bias. Based on unpublished data from R. N. Lemon, K.M.B. Bennett, and D. Flament.

The theoretical expressions developed by Humphrey and Corrie (1978) to correct PTN conduction velocity distributions for the sampling bias toward large cells in the recordings were applied to the data of Fig. 3. The result of such a correction is shown in Fig. 3 (- - -).

In the simulations, the conduction velocity of each SCM cell in the population was chosen by reading off conduction velocity values from the corrected plot of Fig. 3 at equally spaced intervals along the cumulative probability axis. Thus if the cumulative distribution is represented as V(P), then the conduction velocity of the nth cell &cjs1726;n out of a total of n was determined by
&cjs1726;<SUB><IT>n</IT></SUB>= <IT>V</IT>([(2<IT>n</IT>+ 1]/[2<IT>n</IT>+ 2]) (7)
The conduction time from the cortex to the spinal segment of termination for each SCM cell was derived assuming a conduction distance of 100 mm, producing a variation in conduction time from 1.6 ms to 10 ms. Discharges from the fastest and slowest conducting SCM cells were used in the STA analysis.

The unsynchronized CM cell was simulated in all cases with a conduction time of 1.4 ms, making it the fastest conducting cell in the population. Changes in this parameter would simply shift effects in time, without changing their form.

Before making synaptic contacts with motoneurons, corticospinal axons branch and become narrower (Lawrence et al. 1985; Shinoda et al. 1976, 1979). Shinoda et al. (1986) estimated the terminal conduction velocity to be as low as 1 m/s. In addition, motoneurons in a single motor nucleus can extend over two to three spinal segments (Jenny and Inukai 1983), implying a variable conduction distance along the stem axon from cortex to individual motoneurons. The delay between cortical cell firing and the production of an EPSP in a motoneuron could therefore vary considerably from one motoneuron to the next. To simulate the effect of such terminal slowing of impulse conduction, an additional delay was added to the conduction time from each CM or SCM cell to each motoneuron. The extra delay was randomly determined from a uniform distribution between 0 and 1 ms. Shinoda et al. (1986) found <= 1-ms conduction time in the terminal collaterals. The use of a uniform distribution <= 1 ms may somewhat overestimate the variability in conduction time; however, the cross-correlation peak between a CM cell and motor unit which these assumptions produced had a similar width to those published by Lemon et al. (1990).

Implementation

Simulations were carried out with a 0.2-ms time step. In one pass, each cell was first updated; for the motoneurons this required solution of the seven differential equations given by Powers (1993), by using an exponential integration scheme (MacGregor 1987). If any cortical cell produced an action potential, a EPSP was set up in its target motoneurons after a preset delay. Whenever a motoneuron fired, it was followed by the production of its surface MUAP. The instantaneous value of the total surface EMG potential was determined by summing the appropriate points of all currently occurring motor unit potentials. This voltage value was then written to disk together with the times of any cell action potentials that had occurred.

Simulations were run on 486 and Pentium PC-compatible computers using a program custom written in Pascal (Borland). Approximately 8 h of computer time were required to simulate 1,000 cortical cell spikes on a 50-MHz 486 computer; this was reduced to <2 h on a 133-MHz Pentium-based machine. To produce results in a reasonable time, simulations were run simultaneously on <= 25 computers, with a different initial seed for the random number generator in each case. The separate data files were then combined end to end to produce the final data for analysis. Spikes from the initial and final one second of each file so combined were deleted so that only periods when the network had reached a stable activity pattern were included in the analysis. The files so produced were analyzed using a commercially available neurophysiological data analysis package (Spike2, Cambridge Electronic Design). STAs of the simulated rectified EMG were compiled using spikes from each of the three classes of cortical cells (CM, SCM, and S).

    RESULTS
Abstract
Introduction
Methods
Results
Discussion
References

Figure 4 presents analysis of three simulation runs. In Fig. 4A, a cross-correlation histogram is presented between two SCM cells; this will be identical for all SCM-SCM pairs and for each SCM cell cross-correlated with the S cell. This shows the central peak caused by common synaptic input, with a peak width of ~15 ms. In Fig. 4, parts B and C, the interspike interval histograms are shown for a cortical and motoneuronal cell, respectively; in each case they resemble those seen in experimental data. Note that the cortical cell model was designed to produce an interval histogram with a gamma distribution. By contrast, the variability in the interspike intervals in the motoneuron model of Fig. 4C resulted from the simulated noise, in agreement with the proposal of Calvin and Stevens (1968).


View larger version (33K):
[in this window]
[in a new window]
 
FIG. 4. Analysis of 3 simulations, illustrating features of the model and the effect of changing the number of synchronized CM (SCM) cells in the network. A: cross-correlation between 2 SCM cells, showing a central peak 15-ms wide. Size of this peak agrees with experimental observations. B and C: interval histograms of cortical cells and motoneurons, respectively, compiled with 1-ms binwidth. D-F: spike-triggered averaging (STAs) of rectified EMG triggered from 4 different simulated cortical cells. Constructed using MUAPs from subject NA. SCMs and SCMf are the slowest and fastest SCM cells present in each network. Different panels result from simulations with different numbers of SCM cells in the network (10, 20, and 30 cells for D-F, respectively). Arrows, time of the triggering spike; dotted lines, earliest possible onset latency for a causal, monosynaptic effect from a fast CM cell. G: cross-correlation histograms, constructed with a cortical cell as the trigger and all 103 motoneurons in the pool as the response. Same simulation as F, with 30 SCM cells in the network. Motoneuron firing was registered at the spinal cord; hence, latency of the effects does not include the peripheral conduction delay, unlike the STAs. Dotted lines, earliest possible onset latency for a causal effect (1.4 ms). All measures compiled from 10,000 trigger spikes.

Figure 4, D-F, presents STAs of the rectified EMG compiled from four different cortical cells: the unsynchronized CM cell, the slowest and fastest conducting SCM cells in the simulation, and the synchrony only cell (S). The moment of spike discharge is indicated by the arrow at time 0. The results of three simulations are shown, differing only in the number of SCM cells present in the cortical network. This was 10 cells in D, 20 in E, and 30 in F. Averages were compiled from 10,000 spike triggers, a number often quoted in the STA literature as sufficient to show most PSFs (Buys et al. 1986; Fetz and Cheney 1980; Lemon et al. 1986).

It is apparent from Fig. 4 that all cell types produced a PSF which exceeded the noise level in the averages. This was true even when only 10 SCM cells were synchronized together (Fig. 4D). However, the effects produced purely by synchrony (S), with no causal connections to the motoneuron pool, were qualitatively different. The pure CM PSF was brief, and discharge of the motor units was sufficiently synchronized to the triggering spike that the PSF appeared bifid, representing the biphasic nature of the MUAPs coming through in the average of rectified EMG. By contrast, the S cell PSFs were broader and smoothed out. They were the same amplitude as the CM cell PSF only at the highest level of synchrony investigated (30 SCM cells in the network; Fig. 4F). The SCM cell effects, as would be expected, resembled the sum of the CM and S effects, with the CM effect shifted to be appropriate for the conduction delay of the particular SCM cell used as the trigger. Figure 4G shows cross-correlation histograms between each cortical cell class and the entire motoneuron pool (discharge recorded at the spinal cord so that the peripheral conduction delay is not included) for the same simulation as in Fig. 4F, with 30 SCM cells present in the network. Here the summation of the brief peak of the CM cell with the broad peak of the S cell in the SCM cell histograms is readily apparent.

There has been much emphasis in the literature on measurement of the latency of a PSF as a means of excluding effects that are purely due to synchrony. For each average in Fig. 4, the dotted line indicates the earliest latency at which a cortical cell can have a causal effect on the EMG. This is formed by summing the conduction time from cortex to motoneuron with the shortest peripheral conduction time (8 ms). As expected, the CM cell effect occurred shortly after the earliest possible causal latency. However, so too did the effects from S and the fastest SCM (SCMf) cells, at least as judged from these relatively noisy averages of 10,000 spikes.

These measurements of the latency of PSF from STAs of 10,000 spikes should be compared with those from the cross-correlation histograms in Fig. 4G. These histograms were calculated using the cortical cell spikes as a trigger and the multiple unit activity of all motoneuron spikes simulated as the response. The signal-to-noise ratio should therefore be comparable with that of the STA because the same number of spikes contributed to both measures. In the cross-correlation, the peaks for the three synchronized cortical cell classes all began before the earliest causal onset latency time for a fast CM connection (dotted lines), unambiguously indicating a synchrony component to the generation of these peaks. The reason for the difference probably lies in the different sensitivities of the two methods. Because of cancellation effects, averages of rectified EMG can markedly underestimate small potentials (Baker and Lemon 1995). Thus the slowly rising initial phase of the synchrony peaks could be easily distinguished in the cross-correlation histograms but was lost in the noise of the STAs. The earliest discernable component therefore appeared to have a latency consistent with a wholly causal effect.

To investigate this effect more thoroughly, the simulation presented in Fig. 4F, with 30 SCM cells in the network, was rerun for a much longer period of time, until 470,000 spikes of each of the cortical cells were generated. In Fig. 5A are shown STAs compiled from all 470,000 spikes for the three different cell classes (only the fastest-conducting SCM cell, SCMf, is illustrated). With the now much higher signal-to-noise ratio in these averages, the PSFs from the SCM and S cells can be seen to have a slowly rising component at their onset, which clearly began before the earliest causal onset latency line (dotted line). Figure 5B presents a number of averages compiled from epochs of 10,000 spikes. With the eye guided by the averages above in Fig. 5A, it is possible to discern earlier components in some of these epochs for the SCM and S cells. However, these features consistently failed to rise above the noise.


View larger version (22K):
[in this window]
[in a new window]
 
FIG. 5. A: STAs from the same simulations shown in Fig. 4F (30 SCM cells present), triggered from many more spikes (n = 470,000). ···, earliest possible onset latency for a monosynaptic effect; - - -, peak latency, measured from these relatively noise-free averages; up-arrow , time of the triggering spike. B: averages compiled from epochs of 10,000 spikes with the same data as in A. A considerable variability is present in the apparent onset latency of the post-spike facilitations (PSFs). Peak latency is somewhat more constant for the CM and SCMf cells than for the S cells.

In addition to measurements of onset latency, Fig. 5 also presents measurements of the latency of the peak of the PSFs (dashed lines). This showed a marked difference between cell classes, being 18.2 ms for CM and 18.4 ms for the fastest SCM, but 21.8 ms for the S cell. The peak of the PSF is, by definition, the point at which the signal-to-noise ratio is highest, and in the 10,000-spike epochs of Fig. 5B the peak latency shows little variability around the value determined from the all spike averages in Fig. 5A. This suggests that the peak latency may prove to be of more use than the onset latency in excluding PSFs generated purely by cortical synchrony. However, a disadvantage of peak latency is that it changes with the conduction velocity of the CM cell, such that slower conducting cells than the ones shown in Fig. 5 will have a latency which overlaps with that of the S cell.

Quantitative determination of the separation of causal and synchrony PSFs

To provide an objective and quantitative assessment of whether PSFs generated by pure synchrony can be distinguished from those caused by monosynaptic connections, the 47 epochs of 10,000 spikes resulting from the simulation presented in Fig. 5 were subjected to a further analysis which automatically measured various parameters from each average. It was important that this analysis should be entirely automated, so that no experimenter bias could enter the measurement process, a concern when working with noisy averages of this kind.

The latency of the PSFs was difficult to measure automatically, as the PSF was often preceded by a dip in the average (see Fig. 5A), reflecting the autocorrelation of the triggering unit (Lemon and Mantel 1989; Perkel et al. 1967). To account for this, a regression line can be fitted to a background region of the average and then subtracted (Bennett and Lemon 1994) (cf. Fig. 6B). The SD of the average around the regression fit is then determined, and the point at which the PSF exceeds the fit +2 SD is taken as the onset latency (Fetz and Cheney 1980).


View larger version (17K):
[in this window]
[in a new window]
 
FIG. 6. A: data used to define a significance level for a PSF; 60-ms-long averages of rectified EMG with respect to 10,000 spikes were compiled, with the triggers shifted in time by 1-2 s, thereby destroying any causal link between spike train and EMG. A regression line was fitted to the first 30 ms of the average and subtracted; SD of the 1st 30 ms was then found, together with the maximum voltage in the period from 30-60 ms. This maximum was divided by the SD, forming the measure plotted on the abscissa. Ordinate lists the proportion of averages where it was smaller than the value plotted. In 95% of trials, the 30-ms test period did not exceed 5.70 times SD, providing an appropriate test level to determine the 1st point of a PSF, which is significantly different from the background. Results are illustrated from simulations by using motor unit data from subject NA; similar results were obtained with data from subject MM, with a similar significance level (5.77). B: example of the automatic latency detection method. The trace is an STA of EMG produced with MUAPs from subject NA, triggered by 10,000 spikes of the CM cell (up-arrow ). black-square, 30-ms control region used to fit a regression line (lower - - -); above this is shown regression line +2, and +5.7 times SD of the background region. square , region which was searched for a crossing of the detection threshold.

The validity of this procedure was investigated as follows. Control averages were constructed from the simulated data of Fig. 5, which were triggered by the spikes of a cortical cell, shifted in time by 1-2 s such that any temporal relationship with the EMG was lost. Fluctuations in such averages will be caused only by noise. A regression fit was then made to the first 30 ms of these averages (the "background") and subtracted. The maximum value of the average in the 30-ms period after the background region was determined and expressed as a fraction of the SD of the background region. This procedure was repeated for a large number of such control averages, and a distribution of this corrected maximum found. This is plotted as a cumulative distribution in Fig. 6A, which therefore shows the proportion of averages that did not exceed a certain level. Use of the standard criterion level of mean + 2 SD would result in 60.2% of these control averages exceeding the detection threshold. If a 95% confidence level is required, a criterion of 5.7 SD above the regression fit must instead be used.

Two measurements of onset latency were therefore made from the PSF epochs of 10,000 spikes. The section of the average from 40 to 10 ms before the trigger point was used as a baseline, and its best fit line was calculated. The point at which the average in the period from 10 ms before to 20 ms after the trigger first exceeded this line plus two, or plus 5.7 times the SD of the regression corrected baseline region, was used as the automatic latency measurement. This permitted a comparison of the results by using the commonly used criterion level and a more rigorous, statistically defined level. Figure 6B shows an example of this method applied to a CM cell STA.

In addition to these onset latency measurements, the time of the peak, and the peak width at half-maximum (PWHM) were also determined from each 10,000-spike average. For the PWHM, the algorithm first found the mean level in the baseline region from 40 to 10 ms before the trigger. Working out from the peak of the PSF, it then determined the first points either side of the peak where the average fell below one-half the height of the peak above the baseline. The time separation of these points defined the PWHM.

Figure 7 plots histograms of the distribution of these various measures derived from the PSF for each of the different categories of cortical cell. Figure 7, A and B, shows the two latency measures, derived from a 5.7 times SD and 2 times SD criterion, respectively. In the latter case, a considerable number of averages had an onset latency which is unreasonably early. These were presumably caused by noise crossings of the detection threshold and confirm that such a low criterion is inappropriate. Use of the higher criterion produced less variability, but for the CM cell led to a bimodal latency distribution caused by the slightly bifid nature of the CM cell PSF. However, for both sets of latency measurements, the synchrony-only cell S had a latency distribution which overlapped that of the other cells which had monosynaptic connections to motoneurons. Thus on the basis of comparatively noisy averages of 10,000 spikes, measurements of latency could not distinguish genuine from synchrony only PSFs.


View larger version (34K):
[in this window]
[in a new window]
 
FIG. 7. Distribution of various measures made from STAs of 10,000 spikes for 4 different cortical cells. A and B: PSF onset latency, determined by a 5.7 or 2 times SD threshold, respectively. CM and SCM cell distributions overlap with the S cell, preventing reliable separation of the effects by onset latency. C: peak width at half-maximum (PWHM). D: PSF peak latency. Bimodal distributions for the CM cell in A and slowest SCM (SCMs) in C result from the bifid shape of the monosynaptic component of the PSF in these simulations.

Figure 7C shows the distribution of the PWHM for these cells. This measure provided good separation of the cell classes. The CM cell had the narrowest peak; this was appreciably broadened for the SCM cells where a monosynaptic and synchrony effect combined. The slow SCM cell yielded a bimodal width distribution because the cusp of the bifid monosynaptic part of the PSF fell close to the half-maximum detection level. The synchrony-only S cell had the widest peak of all, and there was little overlap with the peaks of the other cell classes, such that this measure could form the basis of a reliable classification scheme.

Figure 7D shows distributions of peak latency for each cell type. These values were very tightly clustered for the CM cells, both with and without synchrony, reflecting the rapid rise time of the monosynaptic part of the PSF. The peak time measured from the PSF produced by the S cell was by contrast highly variable from epoch to epoch. It was later than that for the fast CM and SCM cells and did not overlap substantially with their distributions. However, the S cell peak time was earlier than that for the slowest SCM cell. Hence, in these simulations, peak time only permitted accurate discrimination of synchrony versus genuine PSFs for fast conducting PTNs.

Effect of different surface motor unit action potentials

Measurements made of timing of effects in STAs differ from those taken from cross-correlation histograms in one crucial respect; they include a component due to the shape and duration of the MUAPs. A concern with the classification presented by Fig. 7 must therefore be that it is highly dependent on the MUAPs used in the simulations. Different sets of potentials could result in altered criteria being needed to separate synchrony from genuine PSFs such that the generality of the results would be lost. To investigate this, the EMG was recalculated from the stored motoneuron firing times of the simulation illustrated in Fig. 7 using diverse sets of MUAPs.

These potentials were formed by taking the standard sets of potentials computed as described in METHODS from the two human subjects and stretching or compressing them in time. The temporal scale factors used were 25-200%, in 25% steps. Thus with a 25% scale factor, the potentials were four times shorter and resembled the fast MUAPs which might be recorded with highly selective intramuscular electrodes. With a 200% scale factor, the potentials lasted twice as long, mimicking a less-selective surface recording from a large forearm muscle. Using eight different scale factors on each of the 2 sets of MUAPs from the different subjects gave 16 sets of potentials in total.

Figure 8 illustrates STAs computed from simulations where these different potentials were used, with the CM or S cells as trigger. Figure 8A presents data using MUAPs from subject NA, and Fig. 8B presents data from subject MM. As the scale factor increased, the PSFs for each cell became wider. Using briefer potentials, the bifid nature of the CM cell PSF was lost, presumably because the interval between each phase of the MUAPs was brief compared with the jitter in motor unit firing relative to the triggering cortical cell. By contrast, averages constructed with the broadest potentials in Fig. 8B showed three phases of the CM cell PSF. However, in spite of these trends, the PSFs from CM and S cells remained clearly distinctive.


View larger version (18K):
[in this window]
[in a new window]
 
FIG. 8. Effect of different MUAPs on PSFs from CM and S cells. Averages in A were constructed from potentials from subject NA, and those in B were constructed from potentials from subject MM. Numbers on the left give the extent of the temporal scaling of the interpolated MUAPs of Fig. 2, with 100% representing no scaling. The briefer potentials lead to monophasic CM cell PSFs, whereas with longer potentials 2 or 3 phases are evident, reflecting the complexity of the waveforms of the underlying motor unit potentials. In all cases, the CM and S effects are clearly different. up-arrow , time of triggering spike. Compiled from 300,000 spikes.

Figure 9 presents the results of automated measurements from epochs of 10,000 spikes, combined across STAs computed with all 16 possible MUAPs. Unsurprisingly, the latency, shown in Fig. 9, A and B, measured with either of the two threshold criteria, failed to discriminate between the genuine and synchrony PSFs, just as it did when measurements were made using the original set of MUAPs in Fig. 7, A and B. Figure 9C shows the PWHM distributions combined across all MUAPs. Although these were broader for each cell class compared with Fig. 7, the CM and S cells were still essentially nonoverlapping. The reason for this is that the width of the synchrony PSFs was predominately determined by the width of the cortical synchrony, together with the dispersion of the synchronized corticospinal axons caused by the range of conduction velocities present. The duration of the MUAP therefore had relatively little effect on the PWHM.


View larger version (41K):
[in this window]
[in a new window]
 
FIG. 9. Histograms of the distribution of measures made from PSFs for 4 different cortical cells. Thirty epochs of 10,000 spikes were measured for simulations by using each of 16 different sets of MUAPs, and all the measures were combined for this figure. Hence, differences seen between cells here will be robust across MUAPs. A and B: PSF onset latency, measured with criteria of 5.7 and 2 times SD. C: PSF peak width. D: PSF peak time.

Figure 9D shows the peak latency measures combined across motor unit potentials. In this case, the combination was disastrous; whereas CM and S cells did not overlap on this measure when only one set of motor unit potentials was used, the distributions were now heavily overlapping. However, each histogram showed multiple peaks. These were caused by the summation of several very narrow peaks from each motor unit potential type. This suggests that, for all classes of MUAPs, PSFs from S and CM cells had reliably different peak times. However, there was no one criterion level which could be used to separate the two categories; the optimum criterion changed with different sets of MUAPs.

Sensitivity to the width of cortical synchrony

One of the most important parameters of the model presented above is the width of the cross-correlation peak between two synchronized cortical cells. The chosen value of 15 ms is the best estimate available from a small data set gathered in our laboratory (Baker 1995). Smith and Fetz (1989) found a somewhat larger mean width of 23 ms. These values are means; the width of correlation peaks between any two cells can vary widely, from as short as a few milliseconds to >30 ms. Additionally, the magnitude of synchronous activity in motor cortex can vary in a task-dependent manner (Baker et al. 1997). Because experimental data on this point are currently so sparse and because the findings may be expected to depend importantly on it, it was decided to investigate a range of different cortical synchrony widths and to determine whether the PWHM remains capable of distinguishing genuine from synchrony PSFs.

Figure 10 shows the results of these simulations. Figure 10A illustrates the cross-correlation peak between two synchronized cortical cells. In these different simulations, the common input to the cells was chosen to make the cross-correlation peak duration ~5, 15, 25, and 35 ms (15 ms was the standard simulation width, used until this point). The size of the EPSP from the common input cell and hence the area of the cross-correlation peak were held constant. Finally, because our recent report has highlighted the importance of oscillations in the motor cortex at ~25 Hz (Baker et al. 1997), we simulated oscillatory synchrony between the cortical cells by adding a 25-Hz, 300-µV sine wave to the membrane potentials of the S and SCM cell. A cross-correlation produced after this manipulation is shown in the final column of Fig. 10.


View larger version (29K):
[in this window]
[in a new window]
 
FIG. 10. Effect of changing width and type of synchrony between cortical cells. A: cross-correlation histograms between 2 SCM cells, triggered by 300,000 spikes. Each column shows results from a different simulation, where the width and type of synchrony was changed. Note that the ordinate scale does not begin at 0. B: STA of rectified EMG triggered by CM (thin lines) and S cells (thick lines). Constructed with unscaled MUAPs from subject NA. Arrows, time of the trigger. n = 300,000, except for the simulation with 15-ms width, where n = 470,000. Calibration bars, 1 µV, 20 ms. C: distribution of the PWHM measure. Distributions are summed over all sets of MUAPs used in this study. For the S cell, the distributions are shown as histograms. For the other 3 cell classes, the distributions are given by bars, representing the region containing 95% of the measurements. CM, solid bars; fastest SCM, shaded bars; slowest SCM, open bars. Two bars are shown for the PWHM of the CM cell in the 5-ms simulation, reflecting the bimodal distribution for this cell.

Figure 10B illustrates STAs of the rectified EMG in each of these simulations. Unscaled MUAPs were used from subject NA to generate these traces. As expected, the STA triggered by the CM cell (thin line) did not vary with changes in synchrony present in the remainder of the cortical network. The STA triggered by the S cell (thick line) shows that increasing the width of cortical synchrony caused the PSF from this cell to widen slightly. Because the cortical cell cross-correlation peak area was conserved, the wider synchrony time course had a smaller peak amplitude in the cross-correlation histograms (Fig. 10A), and this was reflected in the smaller amplitude of the synchrony PSF in Fig. 10B as the cortical synchrony width increased.

For the simulation with oscillatory synchrony in the final column of Fig. 10B, the STA showed oscillations of similar amplitude to the peaks produced by the nonoscillatory synchrony. Although experimentally such oscillatory STAs are rarely encountered, they have been observed (Baker et al. 1997). In experimental recordings, the subpeaks become smaller the further away they are from the main, central peak. This presumably reflects the fact that in vivo the oscillations cover a broadband of frequencies, such that averages of distant cycles accumulate phase errors and lead to cancellation. This was not seen in the simulated STA of Fig. 10B because the oscillations were at a fixed frequency of 25 Hz.

Figure 10C presents the distribution of the PWHM measure for each of the different cortical synchrony widths. As in Fig. 9, these are summed distributions from simulations using the full range of MUAPs available. For simplicity, the distributions of these measures from the CM, fastest SCM, and slowest SCM cells are summarized by bars which encompass 95% of the distribution of the measures for these cells. For the S cell, the distribution histograms are presented in full. The PWHM permitted good separation of CM from S cell effects up to a cortical synchrony width of 25 ms. However, at a width of 35 ms, the histograms overlapped for the two cell classes. At this width the peak became small relative to the noise in averages of 10,000 spikes, such that its PWHM was often underestimated. Extensive overlap was also seen in the oscillatory synchrony case. This also resulted partly from the small size of the PSF but partly from the difficulty in defining a "baseline" for the oscillatory S cell PSF. However, inspection of Fig. 10B shows that the peak width of the S cell oscillatory PSF was not greatly different from that of the CM cell PSF.

The results presented above allow the effect of different criteria of PWHM to be tested. Figure 11A plots the errors made in classifying genuine versus synchrony PSFs using different PWHM criteria for the four classes of CM, fastest and slowest SCM, and S cells. For the CM and SCM cells, the errors are false-negatives, resulting from rejection of a PSF as synchrony generated when in fact it had a monosynaptic component. For the S cells, they are false-positives, resulting from incorrect acceptance of synchrony effects. The data in this figure were derived from Fig. 9 and therefore assumed a cortical synchrony width of 15 ms. It is evident from this plot that a criterion of PWHM <9 ms will provide a significance test for the presence of a CM connection, with an approximate 95% confidence level, since only 4.5% of the S cell PSFs has a PWHM less than this. This criterion leads to only 0.4% of the unsynchronized CM cell PSFs being falsely rejected. However, between 37 and 40% of the SCM cells will be falsely rejected, depending on the conduction velocity of their corticospinal axons.


View larger version (29K):
[in this window]
[in a new window]
 
FIG. 11. A: performance of attempts to classify cell types on the basis of PWHM with different criteria. Data are derived from simulations with a cortical synchrony width of 15 ms, combining across all MUAPs as in Fig. 9. Error rates for CM, fastest SCM (SCMf), and slowest SCM (SCMs) are false-negatives, the proportion of PSFs that was incorrectly rejected as resulting from synchrony effects. Error rates for S cells are false-positives, the proportion of PSFs incorrectly accepted as caused by CM connections. B and C: dependence of classification errors on the parameters of cortical synchrony assumed in the model, with the use of a fixed criterion of PWHM <9 ms (B) or <7 ms (C). D: distribution of PWHM in 76 PSFs measured experimentally from STAs of hand and forearm muscle EMG, triggered by 43 motor cortical cells. Bins shaded black are those excepted to have a monosynaptic component because their PWHM is <7 ms. Those shaded gray would be accepted on the more relaxed criterion of PWHM <9 ms.

A more relaxed criterion can of course be used, which will include more of the SCM effects, but at the expense of an increased false acceptance rate of the synchrony PSFs from the S cells. Figure 11B plots the classification error rates with the use of a 9-ms PWHM criterion for the different simulations described above, with differing widths of cortical synchrony. For synchrony more narrow, or more broad, than the 15 ms cross-correlation width used for Fig. 11A, substantially more S cell PSFs are incorrectly accepted. Given the uncertainty surrounding the correct parameters for the time course of cortical synchrony, the 9-ms criterion therefore seems unduly relaxed as a significance test for whether a PSF contains a monosynaptic component.

Figure 11C illustrates the classification errors using a stricter criterion of PWHM <7 ms. The false acceptance rates of the S cell PSFs are now <5% for cortical synchrony widths from 5 to 25 ms. PSFs which pass this test can therefore be assigned as monosynaptic, nonsynchrony effects with high confidence. However, this is gained at the expense of a high false rejection rate of the SCM cell PSFs.

Comparison with experimental data

The main result of the simulations presented above is that PWHM can be used to separate PSFs due purely to synchrony from those with a monosynaptic component. It is of immediate interest to compare these simulated PWHM distributions with those measured from experimental data, to determine what fraction of PSFs would be rejected by this method. The PWHM was therefore measured from 76 PSFs using data previously obtained in this laboratory (Baker 1995; Bennett and Lemon 1994). These PSFs were seen in STAs of hand and forearm EMG triggered by 43 neurons recorded from the primary motor cortex in 4 M. nemestrina monkeys performing the precision grip task of Lemon et al. (1986). Figure 11D shows a distribution histogram for these PWHM values. Bins shaded in black (56/76, 74%) are those which would be accepted on the criterion developed above of PWHM <7 ms; those shaded in gray (8/76, 10%) are those which would also be accepted if a more relaxed criterion (PWHM <9 ms) were adopted. The unshaded bins (12/76, 16%) would be rejected as possibly caused by synchrony effects.

It can therefore be seen that the majority of PSFs measured experimentally meet the conditions for containing a monosynaptic contribution, even using the stricter of the two criteria illustrated in Fig. 11.

The distribution of antidromic latency for the cells shown in Fig. 11D was calculated for the two categories separated on the basis of the 7-ms PWHM criterion. These two distributions showed no significant differences (t-test, P > 0.8), indicating that there was no bias for acceptance or rejection of a PSF depending on the conduction velocity of the fibers.

    DISCUSSION
Abstract
Introduction
Methods
Results
Discussion
References

Model assumptions

In constructing a model such as the one presented here, it is inevitably necessary to make assumptions where experimental data are incomplete. Before conclusions can be drawn from this model, it is necessary to examine these assumptions.

First, a number of assumptions were made concerning the parameters of the cortical synchrony. For most of the simulations presented, the strength and duration of the synchrony was set in accordance with the mean values found from our own experimental data (Baker 1995). However, there is considerable variability in strength and duration of cross-correlation peaks between pairs of cortical cells. There are often multiple peaks, with a narrow component superimposed on the broad central peak. Nevertheless, the general conclusions of the current work appear unaltered by different widths of cortical synchrony (Figs. 10 and 11), and a good separation of CM from S cell PSFs can still be achieved. The result is therefore likely to be robust in the face of synchrony consisting of a mixture of different cross-correlation widths among different cell pairs in the population.

The cause of cortical synchrony was here assumed to be a single source of common input (represented by one cell in Fig. 1). However, it is likely that there are a number of sources of common input, which connect to different subpopulations of the cortical cells. A large number of different patterns of connectivity could result in the same cross-correlation peak width and size as seen experimentally between pairs of neurons. These effects will act to reduce the mean amount of synchrony between a given cell and the rest of the synchronized cortical population below that which occurs when all cells receive the same common input. Hence the present simulations investigated the pattern of connectivity which would give the largest synchrony PSFs while remaining true to the experimental data on cortical cross-correlations. Different connectivity patterns would give smaller synchrony PSFs, leading to less likelihood of falsely classifying them as caused by monosynaptic connections. However, synchrony can vary significantly during different phases of task performance (Baker et al. 1997), and it is possible that cortical synchronisation could be considerably greater in certain highly specific circumstances.

It was also assumed that synchrony is equally prevalent in CM cells with different conduction velocities. There is currently no experimental evidence concerning how common input is distributed among slow and fast CM cells. It is conceivable that the slow CM cells could receive common input from a different source than the fast CM cells. This is given some plausibility, since slow and fast CM cells projecting to the same muscle appear to fire differently during a movement, consistent with a different role in movement control for the two subpopulations (Maier et al. 1993). In addition, there is evidence that the intracortical connections of fast and slow PTNs differ (Kang et al. 1991; Tsukahara et al. 1968).

At the spinal level, it was assumed that each CM cell evokes the same EPSP in each motoneuron in the pool. Precise quantitative details are not available on this, but the assumption is supported by the finding of Mantel and Lemon (1987) that a CM cell which produces a PSF in a muscle has a cross-correlation peak with most single motor units which could be recorded in that muscle. It is possible that there are consistent changes in EPSP size with motoneuron recruitment order, as discussed by Heckman and Binder (1993). However, these would be expected to have less effect on the situation modeled here, where only the lower 27% of the entire motoneuron pool is active (103/377 cells).

There is some evidence that slowly conducting CM cells produce smaller PSFs, suggesting that they may make weaker connections to motoneurons (Fetz and Cheney 1980; Lemon et al. 1993). However, without data from single unit cross-correlations, no precise values for this effect can be included in the model.

We did not include in out model facilitation of motoneurons by corticospinal cells mediated by nonmonosynaptic linkages. Although interneurons do produce such facilitation (Fetz et al. 1996), it is not known if they are intercalated in the corticospinal pathway. However, in the anaesthetized monkey, there is a conspicuous lack of nonmonosynaptic EPSPs in motoneurons following pyramidal tract stimulation (Maier et al. 1996), suggesting that such connections may be weak.

Possibility of PSF caused by cortical synchrony

One of the most important predictions of the simulations presented here is that a cortical cell which makes no causal connections to motoneurons but is synchronized with other cells that do could produce a significant PSF of rectified EMG. This disagrees with the reports in the literature which have addressed this question to date (Lemon et al. 1985; Porter and Lemon 1993; Smith and Fetz 1989), which concluded that cortical synchrony played little role in the production of PSF.

The current simulation results show that, if the synchronized CM colony contains as few as 10 cells, then another neuron which is also synchronized with these but is not a CM cell can produce a significant PSF. No data exist on the number of CM cells which may be synchronized together in this way because recordings of no more than two CM cells have been reported. However, three observations are relevant here.

First, in recordings from monkey motor cortex, pairs of cells have been observed which are synchronized together as judged from a cross-correlation histogram; however, one cell produced a strong PSF, and the other cell produced no PSF at all (Porter and Lemon 1993). Second, the yield of PTNs showing PSF in an experiment is relatively small (55% from a highly selected population of PTNs) (Fetz and Cheney 1980; Lemon et al. 1993). Finally, the peak amplitude of most experimental PSFs relative to background (percentage modulation) is <15% (Lemon et al. 1986). In the simulation presented in Fig. 5, the fast SCM cell PSF had a percentage modulation of 15.2%. If >30 CM cells were synchronized together, the size of the PSF would be larger than that seen experimentally. This would be true independent of how many sources of common input projected to the SCM cell pool. The current model predicts that no CM cell can have a cross-correlation peak of the same amplitude as assumed here with >30 other CM cells projecting to the same muscle without producing a PSF larger than that encountered experimentally. We chose to vary the strength of cortical synchrony by keeping the strength of correlation between individual cells fixed and altering the size of the synchronized colony. This is equivalent to changing the magnitude of the correlation peak with a fixed colony; both effect the number of action potentials in the corticospinal input synchronized with the trigger.

Taken together, the implication of these findings is that the level of synchrony used here is at the upper limit of that present in the cortex.

Differences between genuine causal and synchrony PSFs

In all simulations illustrated, there were clear qualitative differences in the PSFs produced by the CM, SCM, and S cells. The S cell peaks were broader than the CM and SCM peaks and had a later time to peak. Both S and SCM cells generated PSFs with earlier onset latencies than reasonable for a causal effect (Fig. 5). These differences are clearest in the relatively noise-free averages of Figs. 5 and 8.

One important difference between CM and S cell PSFs was that the CM effects were often bifid in appearance (e.g., Figs. 4 and 5). However, Fig. 8 shows that when MUAPs were used which were shorter than those interpolated from the human recordings, the bifid nature of the CM cell PSFs was lost. It is interesting that, in experimental studies in the monkey, PSFs usually only exhibit a single peak. The small size of muscles in the monkey compared with the human, with correspondingly smaller interelectrode distances, may mean that MUAPs are briefer (see data of Mantel and Lemon 1987); the data of Fig. 8 with compressed MUAPs then most accurately reflects the experimental situation, with corresponding monophasic PSFs. A second factor contributing to the production of monophasic PSFs in experimental data may be the presence of a superimposed synchrony PSF; for example, the SCMf effects for the two higher levels of synchrony investigated in Fig. 4 show no inflections on their PSFs larger than the noise level. Finally, it is possible that MUAPs show more variation in their shape and duration than accounted for by the present method of interpolating from a limited sample of recorded potentials; this increased variability could then act to smooth out the bifid cusp of the CM cell PSFs.

As described above, the use of PSF latency to exclude S cells from being classified as CM has considerable practical difficulties due to the relatively high noise level present in averages of 10,000 spikes, a number commonly considered sufficient to detect PSFs (Lemon et al. 1986). With this number of triggering spikes, onset latency cannot be used reliably to discriminate synchrony from genuine PSFs (Figs. 7, A and B, and 9, A and B). Errors can be made in both directions. Hence an onset latency which is too early to be produced by monosynaptic effects must indicate a synchrony component to the peak; however, it does not preclude the presence of a superimposed genuine PSF. Conversely, a PSF onset at appropriate latency could be produced by purely synchronous effects. Given the very shallow initial part of the S cell PSF and the tendency of averages of rectified EMG to underestimate small effects (Baker and Lemon 1995), an enormous number of spikes need to be averaged to measure these latency differences accurately. Recording in the order of 300,000 spikes is impractical, given the limitations on the length of experimental sessions when working with chronic animals and the difficulties in holding cells in the record for long periods. Any other measure that depends on accurate knowledge of the onset or offset of the PSF, such as rise time or total duration, will suffer from the same problems.

A measure that achieved good separation of genuine from pure synchrony PSFs was the PWHM. The same criterion of PWHM <7.0 ms correctly classified CM and S cells for cortical synchrony widths from 5 to 25 ms, making it an impressively robust measure. The reason for the relative invariance of PWHM with cortical synchrony width is that the synchrony between inputs to motoneurons is wider than that between cortical cells because of the modeled distribution of corticospinal conduction delays. Hence, even with relatively tight cortical synchrony (e.g., 5-ms cross-correlation width), the dispersion in the synchronous input to the motoneurons is still sufficient to broaden the S cell PSF significantly compared with the CM cell PSF (cf. Davies et al. 1985). As noted above, the dispersion in corticospinal conduction time is one of the main assumptions of this model. However, it may safely be assumed that the model would be robust to changes in this dispersion so long as the cortical synchrony width was not also allowed to become too narrow.

The PWHM criterion produced a high proportion of false-negatives for the SCM cells, as illustrated by Fig. 11, A-C. This is because a considerable proportion of the PSF amplitude for these cells is generated by synchrony. However, as argued above, the levels of cortical synchrony used here are likely to be at the upper end of those consistent with experimental data. Hence, in practice, many SCM cells would be expected to have PSFs with a far smaller synchrony component and a PWHM which falls well within the criterion for a genuine monosynaptic effect.

The PWHM of CM cell PSFs shows some slight sensitivity to the MUAP duration, as can be seen by comparing the distribution histograms of Fig. 7C (for a single set of MUAPs) and Fig. 9C (measures from many different MUAP sets combined). However, for the cortical cells which are synchronized, the PWHM measurements do not appear greatly affected by the size or shape of the MUAPs. This is a useful property of the PWHM measure compared with the peak time, since in an experimental situation it is difficult to quantify MUAP duration.

A PSF with PWHM <7 ms is therefore unlikely to be produced purely by synchrony effects. This is the case for a range of assumed synchrony time courses at the cortex (Fig. 11C). As shown by Fig. 11D, the majority of experimental PSFs meet this criterion, adding credence to the STA technique. The status of experimental PSFs which are wider than this remains uncertain. It is not possible to be sure that they have a monosynaptic component, as cell type S in the simulation was clearly capable of producing such PSFs. However, as noted above, the level of synchrony used here is probably at the upper limit of that which is compatible with experimental data. At weaker levels of cortical synchrony, S cell PSFs would be smaller and possibly not rise above the noise. Figure 11C makes it clear that a high percentage of SCM cells with both synchronized and monosynaptic components to their PSFs will be rejected by the PWHM criterion. Future experimental studies would therefore be advised to treat PSFs with PWHM <7 ms as monosynaptic in origin, while allocating the remainder to a different category whose origin which remains uncertain.

Quantification of PSF

An important finding of this study is that the size of a PSF from a cell which does make connections to motoneurons can be heavily influenced by the extent to which it is synchronized with other CM cells. Thus in Fig. 4, whereas the fast SCM cell PSF is only slightly larger than that of the CM cell when 10 SCM cells are synchronized together (Fig. 4D), it is more than twice the height when 30 SCM cells are present (Fig. 4F). Measures made from the area under the peak (e.g., mean percentage increase) (Cope et al. 1987) will be even more susceptible, given the broad time course of the synchrony contribution. Little can be done to correct for this synchrony contribution in an experimental PSF where the size of the synchronized CM cell pool is unknown, and the only solution would seem to be the exercise of great caution in the use of STA to quantify the strength of the connection from a CM cell to the motoneuron pool innervating a particular muscle (see Bennett and Lemon 1994).

    ACKNOWLEDGEMENTS

  The authors thank Dr. Peter Kirkwood for critical discussion of this work.

  This work was funded by the Wellcome Trust and the Medical Research Council.

    FOOTNOTES

  Address reprint requests to S. N. Baker.

  Received 22 December 1997; accepted in final form 4 June 1998.

    REFERENCES
Abstract
Introduction
Methods
Results
Discussion
References

0022-3077/98 $5.00 Copyright ©1998 The American Physiological Society