Previous Article | Next Article 
The Journal of Neuroscience, September 15, 1999, 19(18):8083-8093
Neuronal Interactions Improve Cortical Population Coding of
Movement Direction
E. M.
Maynard1,
N. G.
Hatsopoulos2,
C. L.
Ojakangas2,
B. D.
Acuna2,
J. N.
Sanes2,
R. A.
Normann1, and
J. P.
Donoghue2
1 Department of Bioengineering, University of Utah,
Salt Lake City, Utah 84112, and 2 Department of
Neuroscience, Brown University, Providence, Rhode Island 02912
 |
ABSTRACT |
Interactions among groups of neurons in primary motor cortex (MI)
may convey information about motor behavior. We investigated the
information carried by interactions in MI of macaque monkeys using a
novel multielectrode array to record simultaneously from 12-16 neurons
during an arm-reaching task. Pairs of simultaneously recorded cells
revealed significant correlations in their trial-to-trial firing rate
variation when estimated over broad (600 msec) time intervals. This
covariation was only weakly related to the preferred directions of the
individual MI neurons estimated from the firing rate and did not vary
significantly with interelectrode distance. Most significantly, in a
portion of cell pairs, correlation strength varied with the direction
of the arm movement. We evaluated to what extent correlated activity
provided additional information about movement direction beyond that
available in single neuron firing rate. A multivariate statistical
model successfully classified direction from single trials of neural
data. However, classification was consistently better when correlations
were incorporated into the model as compared to one in which neurons
were treated as independent encoders. Information-theoretic analysis
demonstrated that interactions caused by correlated activity carry
additional information about movement direction beyond that based on
the firing rates of independently acting neurons. These results also show that cortical representations incorporating higher order features
of population activity would be richer than codes based solely on
firing rate, if such information can exploited by the nervous system.
Key words:
motor cortex; monkey; population coding; movement
direction; correlation; cell assemblies
 |
INTRODUCTION |
A major problem in understanding
higher brain function is determining how complex neuronal
representations of sensory and motor information are constructed from
the activities of individual neurons. Higher order sensory and motor
representations appear to emerge from the firing of neuronal
assemblies, but it has yet to be determined whether spatial and
temporal interactions contribute to these representations. Interactions
clearly occur within the widespread network of synaptically modifiable
horizontal connections in which cortical neurons are embedded (Huntley
and Jones, 1991
; Hess and Donoghue, 1994
; Rioult-Pedotti et al., 1998
).
Although controversial, various forms of interactions among cortical
neurons, including firing synchrony and local oscillations, have been
considered as key to forming higher level representations (Singer,
1995
; Vaadia et al., 1995
; DeCharms and Merzenich, 1996
; Riehle et al., 1997
; Hatsopoulos et al., 1998b
).
The representation of information by populations of neurons has been
extensively studied in both sensory (Nicolelis et al., 1995
; Meister,
1996
) and motor areas of cortex (Georgopoulos et al., 1988
;
Schwartz, 1993
; Tanaka, 1994
). A common observation is that the
firing rates of single cortical neurons during repeated presentations
of the same stimulus vary considerably (Henry et al., 1973
; Tomko and
Crapper, 1974
; Heggelund and Albus, 1978
; Rose, 1979
; Dean, 1981
; Lee
et al., 1998
; Shadlen and Newsome, 1998
). Whether this variability is
meaningless noise, reflecting the stochastic nature of neurons, or
carries information is widely debated (Richmond and Optican, 1990
; Bair
and Koch, 1996
; Buracas et al., 1998
; Lee et al., 1998
; Oram et al.,
1998
; Shadlen and Newsome, 1998
; Zhang et al., 1998
).
The firing rate of primary motor cortex (MI) neurons broadly reflects
the direction of arm-reaching movements (Georgopoulos et al., 1986
).
The broad tuning of individual cells to movement direction coupled with
their firing rate variability (Lee et al., 1998
) have inspired theories
of directional coding in which an accurate representation emerges as a
result of averaging over a large number of cells (Georgopoulos et al.,
1992
; Salinas and Abbott, 1994
; Sanger, 1994
; Tanaka, 1994
). Although
these methods are quite successful in deriving movement direction,
response variability was treated as a simple random noise process. This assumption is called into question both because individual MI neurons
nearly always represent additional motor variables (Kalaska, 1988
; Fu
et al., 1993
) and because weak correlations have been encountered in
the firing rate variability of pairs of cortical neurons (Gawne and
Richmond, 1993
; Zohary et al., 1994
; Lee et al., 1998
). The
contribution of correlated discharge has largely been dismissed as
reflecting coding redundancy in a noisy system. However, correlations
in response variability can provide additional information that is
available only from the assembly (Deadwyler and Hampson, 1996
;
Hatsopoulos et al., 1998b
; Oram et al., 1998
). In the present study we
demonstrate that statistical dependencies in the trial-to-trial firing
rates of pairs of neurons can be an additional source of information
about motor behavior.
 |
MATERIALS AND METHODS |
Task. Two male monkeys (macaca
fasicularis) were operantly conditioned to perform an
instructed-delay task consisting of visually guided planar reaching
movements from a central holding position to a radially located target
position (Georgopoulos et al., 1986
; Hatsopoulos et al., 1998a
).
Animals moved a two-joint manipulandum in the horizontal plane to
direct a cursor from a central hold position to one of eight possible
radially positioned targets (1 cm2, 6 cm
away from start position) viewed on a computer monitor. A trial was
composed of three epochs: a hold period during which time the monkey
had to maintain the cursor at the hold position for 0.5 sec, a random
1-1.5 sec "instructed delay" period during which the target for
the forthcoming action appeared but movement was withheld, and a
"go" period initiated by target blinking (mean reaction time,
~365 msec). Manipulandum position was monitored using a digitizing
tablet (Numonics) sampled at 72 Hz. Monkeys generally completed >30
movements in each of the eight directions during a recording session.
Array recordings. When >85% correct trials were attained,
an array of 100 microelectrodes was implanted in the MI arm
representation, medial to the spur of the arcuate sulcus in the
mediolateral dimension, abutting the central sulcus in the
anteropoterior dimension. The microelectrode arrays (Bionic
Technologies, Salt Lake City, UT) consisted of 100, 1.5-mm-long
platinized tip, silicon probes arranged in a square grid (400 µm on
center) with probe impedances between 200 and 500 k
(1 nA, 1 kHz
sine wave; Nordhausen et al., 1996
). Twenty-two in the central matrix
of the 100 available electrodes were selected in advance for electrical
connection to one of two percutaneous connectors using 25.4 µm (1 mil) Teflon-coated platinum-iridium (10%) wire (10IR1T;
MedWire, New York, NY). An extra wire (50.8 µm diameter, Pt-Ir
20%) with ~5.0 mm of the terminal insulation removed was inserted
subdurally and used as a reference. For added ease of handling, the 12 wires were held together in four or five spots with silicone elastomer
(MDX4-4210; Dow Corning, Midland, MI). Each percutaneous connector
consisted of a 12-pin connector (FR-12S-6; Microtech,
Pittsburgh, PA) set in a custom-designed titanium pedestal. The
back of the array and percutaneous connectors were coated with silicone
elastomer to mechanically protect the wires and maintain electrical
insulation at the bond pad sites.
Animal care and surgical procedures used in the identification of the
arm area in monkey MI are described elsewhere (Donoghue et al., 1998
),
as are the procedures for implanting the array and closing the surgical
defect (Rousche and Normann, 1992
). All surgical procedures conducted
for this study were in accordance with Brown University Institutional
Animal Care and Use Committee-approved protocols and the Guide for the
Care and Use of Laboratory Animals (National Institute of Health
publication no. 85-23, revised 1985).
Recordings began after a 1 week period of postsurgical recovery. Custom
shielded cabling carried signals from either 11 or 22 electrodes to a
custom 16-channel amplifier (gain, 25,000;
3 bB points 250 Hz and 7.5 kHz; 5th order Butterworth anti-aliasing filter; input impedance,
109
). Signals were simultaneously
digitized at 20 kHz (Experimenter's Workbench; DataWave Technologies,
Longmont, CO) and processed to extract action potentials. All
candidate action potential waveforms (
0.5 to +1.5 msec around
threshold crossing) exceeding ~2× background noise were stored for
off-line processing and analysis. A given data set included the
simultaneously recorded neural responses from 16 electrodes and
kinematic data in the form of the digitized x-y position of the manipulandum.
Array electrodes often recorded action potentials originating from more
than one neuron (Nordhausen et al., 1996
; Maynard et al., 1997
) and
ranged from multiple nonclassifiable units to well isolated potentials.
Classification of the action potential waveforms was accomplished using
a semiautomatic spike-sorting algorithm based on eight waveform
parameters (Autocut; DataWave Technologies). The classification
templates resulting from this process were subsequently verified and
adjusted manually. In a small number of the identified neurons, it was
likely that a few appropriate waveforms were excluded in the process of
removing incorrectly classified waveforms. Units with stable,
identifiable waveform shapes and a refractory period >1 msec
demonstrated in autocorrelograms were considered to be well-isolated
single units. Action potentials were reduced to an event time series
for further analysis. Based on the lack of a clear refractory period,
13 of the 85 (15.3%) recorded neurons were excluded from further
analysis. Some of the remaining neurons showing a few events in the
refractory period were nevertheless included in the analyses to
maximize sample size and labeled as multiunit neurons (Table
1). As will be shown below, results from
data sets that included these multiunits were similar to those sets
with only single units.
For each identified unit, the mean firing rate during each movement
trial was calculated from the number of action potentials present in
the 600 msec interval beginning 200 msec before movement onset and
ending 400 msec afterwards. This interval most commonly showed the
greatest firing rate change from background for the recorded units.
Movement onset was defined as the time that hand tangential velocity
exceeded 1.0 cm/sec.
Statistical analyses. The statistical significance of an
observed correlation coefficient for a pair of neurons was assessed by
comparing it with a correlation coefficient computed after the sequence
of trials associated with one neuron was randomly shuffled relative to
that of the other neuron. Shuffling removes trial specific correlations
and creates an "independent" condition similar to that which would
be obtained in serially recorded (i.e., one neuron at a time) data. For
each cell pair in each movement direction, 10 unshuffled subsets of
data (>29 movement trials per subset) were selected from the complete
data set and used to compute 10 correlation coefficients and then
averaged. The temporal sequence of these trials was then randomly
shuffled, and a set of 10 control correlation coefficients was
computed. The observed average correlation coefficient for each pair of cells in each of the eight movement directions was assessed for significance by comparing it to the mean of the 10 control correlation coefficients using a two-sample t test
(p < 0.01).
In practice, the covariance matrices that characterized the
second-order multivariate Gaussian model used for classification (see
Fig. 8) for more than two neurons were not always invertible because
singularities result from small numbers of trials. To avoid this
problem, we resorted to generating realizations of a multivariate
Gaussian distribution with a covariance structure (and mean) given by
the actual neural data using a pseudorandom number generator without
having to invert the covariance matrix (Johnson, 1994
). One hundred
realizations were generated for each of the eight movement directions.
Each trial of the actual neural data were then compared to the
k-nearest neighbor realizations (using a Euclidean distance
measure) and assigned to the movement direction that generated the
majority of the k-nearest neighbor realizations. The nearest
neighbor classifier error rate never exceeds twice that of the optimal
Bayesian classifier and will approach the Bayes' error rate as the
number of realizations and k go to infinity (Duda and Hart,
1973
). A variety of values for k were used and were found to
generate qualitatively similar results. The results reported here were
based on a value of 10. To ensure that the classifier was not
overfitting the neural data used to construct it, the classifier
was cross-validated by testing it on neural data that was not used to
build it. Cross-validation was implemented by leaving one trial
out and generating the mean and covariance parameters based on the
remaining trials. This was repeated for all trials in the data set.
 |
RESULTS |
The data used for this study were obtained from five recording
sessions in two monkeys made over a 1 year period. Most neurons demonstrated a marked modulation in their average firing rate around
arm movement but did not modulate their firing rate strongly if at all
during the instructed-delay period (Fig.
1A). To characterize the directional tuning of each neuron, a linear regression was used to
fit a cosine function to the observed firing rates in the eight
movement directions (Schwartz et al., 1988
). The directional tuning
curves of 64 of the 72 (89%) units were fit by a cosine function
(p < 0.05). The peak of this cosine function
identified the preferred direction of that cell (Schwartz et
al., 1988
). The aggregate of preferred directions did not uniformly
represent the movement space (Fig. 1B) in either
monkey, particularly between 0 and 90°, over all five data sets (Fig.
1B, bottom plot). This nonuniformity is
most likely attributable to limited sample size because others have
found directional tuning spanning the entire two-dimensional workspace
(Georgopoulos et al., 1986
).

View larger version (42K):
[in this window]
[in a new window]
|
Figure 1.
A, Raster plots and perievent
histograms of four simultaneously recorded neurons in monkey MI for
movements in all eight directions. Each vertical mark in
the raster trace is a single action potential attributed to a single
neuron. Each row of the raster is a single trial in
which the monkey moved from a central hold position to a radially
displaced target position. Rasters and histograms are aligned to the
onset of movement (dotted line) and show the interval
from 1 sec before to 500 msec after the onset of movement. The
instruction signal (solid square), go cue (solid
triangle), and reward (solid circle) are also
shown in the rasters. Each perievent histogram was constructed with a
bin width of 20 msec and then smoothed. B, Distributions
of the preferred directions for neurons recorded from one data set
(top) and from all data sets combined
(bottom). The preferred direction of each unit was
determined by fitting a cosine curve to the directional tuning function
of the unit.
|
|
We recorded data from the same location in one animal in sessions
separated by 23, 26, and 329 days. A preliminary stability analysis based on the inspection of waveforms, autocorrelograms, and
perievent histograms in one monkey suggested that >50% of the sites
recorded action potentials from different neurons after a 1 month
period and 100% after 14 months. Therefore, some of the recording
sessions may have contained partially overlapping populations.
Variability in neural discharge
As is evident in Figure 1A, the modulation in MI
firing during repeated movements in the same direction vary somewhat
from trial to trial. Response variability for the set of MI neurons was
quantified by comparing the SD of the observed firing rates for
all of the movements performed in a given direction to the mean firing
rate for that direction. The relationship between the trial-to-trial
variability of the response and the mean firing rate is best fit by a
power function (Fig. 2). The dependence of the variability (SD) on the mean firing rate (M) is described by a
power function of the form: SD = 1.44 M0.50
(r2 = 0.85). This relationship
is consistent with the findings of Lee et al. (1998)
who derived a
relationship of SD = 1.21 M0.59
(r2 = 0.77). The variance
exceeds that predicted by a Poisson count in which the mean equals the
variance. Although the origin of this excess variance is unclear, we
will show that some of this variability as revealed from pairwise
correlations can be accounted for by the covariation between cells.

View larger version (22K):
[in this window]
[in a new window]
|
Figure 2.
Trial-to-trial variability in the firing rates of
individual neurons for each of the eight movement directions. The SD in
the firing rate of each cell for each of the eight movement directions
(open circles) is plotted as a function of the average
firing rate of the cell in each direction. The dashed
line is the result of fitting a power function to this
distribution. The solid line is the relationship between
the SD and mean predicted from a Poisson process.
|
|
Covarying neural discharge
To study correlations among MI neurons, the following equation
(Eq. 1) was used to normalize responses to their respective z-scores on a trial-by-trial basis (Crow et al., 1960
; Lee
et al., 1998
):
|
(1)
|
In Equation 1, i indexes a particular neuron,
k specifies the trial number, d specifies the
movement direction, md is the mean, and
d is the SD of the firing rates for movements in
a given direction, d.
Although responses of single neurons were variable from trial to trial,
this variability was correlated for many pairs of MI neurons. Figure
6A (left) depicts an exceptional example
of correlated activity on a trial-to-trial basis for a pair of neurons simultaneously recorded from MI. The correlation strength across MI
neuron pairs was evaluated by computing the correlation coefficient,
ijd,
between the responses of neuron i and j for each
movement direction d, separately:
|
(2)
|
where E[... ] represents the expected value of the
product of the firing rate z scores for neurons i
and j for movement direction d. This expected
value was estimated from the sample mean over all trials for one
movement direction. A Fisher's z transform was used to
convert the correlation coefficients so that they became normally
distributed and so that populations of correlation coefficients could
be more easily compared (Lee et al., 1998
):
|
(3)
|
Correlation coefficients were computed for each pair of cells,
separately for each direction. The average correlation for all pairs of
cells in the eight movement directions was 0.21 ± 0.34. The
slight bias toward positive correlation in the trial-to-trial variability is consistent with other reported correlations in cortical
neurons [e.g., 0.12 by Zohary et al. (1994)
and Lee et al. (1998)
].
In the entire population of simultaneously recorded MI neurons for each
of the eight movement directions, significant correlations
(t test; p < 0.01) were found in 2538 of
the 3976 (63.8%) pairs. Figure 3 shows
the broad distribution of the 2502 significant correlation coefficients
between pairs of simultaneously recorded MI neurons. Note that a few of
the interactions present in the pairs of MI neurons for movements in a
particular direction were quite strong (i.e., correlation coefficients
of ~0.8 or higher).

View larger version (45K):
[in this window]
[in a new window]
|
Figure 3.
Distribution of response correlations between
pairs of MI neurons that were significantly different
(p < 0.01) from the shuffled responses. The
distribution is constructed from 2502 correlation coefficients
calculated for each cell pair in each of the eight movement directions.
The bimodal distribution is a result of the removal of the
nonsignificant correlation coefficients.
|
|
Correlation and directional preference
The dependence of the pairwise response correlation on direction
preference was investigated. A strong relationship between the
similarity of directional tuning and response correlation would be
expected if the correlation resulted from redundant information about
movement for cell pairs. Directional tuning similarity for each pair of
simultaneously recorded MI neurons was evaluated by computing the dot
product between the preferred directions of the two cells. This
procedure gives a quantitative measure of similarity between two vector
quantities such as direction because the dot product is proportional to
the cosine of the angle subtended by the two vectors. Thus, pairs of
cells with a dot product near 1.0 have very similar preferred
directions. Correlation coefficients between the responses of neuron
pairs as a function of the similarity in their directional tuning
curves are shown in Figure 4. Although a
linear (slope = 0.09; r2 = 0.74) relationship was found between the mean correlation coefficient and the similarity of the preferred directions, the variability in
correlation strength was very large over different cell pairs. That is,
different cell pairs could exhibit very different correlation coefficients even when the similarity in their preferred directions was
common. Thus, there appears to be a weak relationship between the
directional preference and correlation strength when correlations are
measured on a broad time scale (i.e., 600 msec).

View larger version (28K):
[in this window]
[in a new window]
|
Figure 4.
Mean response correlation coefficients plotted as
a function of the similarity between the directional tuning curves.
Each mean correlation value is the average of all correlation
coefficients within a 0.1 bin (over the x-axis). A weak
relationship exists between the strength of correlation and the
similarity of the directional tuning of the cells. Error bars indicate
SDs.
|
|
Correlation and interelectrode distance
Based on the distribution of MI afferents from other regions of
cortex, cells that are close to each other have a higher likelihood of
receiving common inputs, thus modulation of the shared input could
result in the observed correlation in firing of neuronal pairs. Shared
input has been a common explanation for coding redundancy when cell
pairs are recorded on single electrodes (Gawne and Richmond, 1993
). The
possibility that correlations between MI neurons were the result of
similar covarying inputs was investigated by plotting the magnitude of
the correlation coefficients as a function of the cortical distance
between the electrodes on which the neurons were recorded (Fig.
5). For the significant correlations in
our population of simultaneously recorded neurons, there was no
significant relationship between the strength of the correlations and
the cortical distance between the neurons (linear regression,
p < 0.36). As is evident in the plot, pairs of neurons
recorded from the same electrode (distance = 0) were no more
likely to be correlated than neurons separated by distances as much as
2 mm. This result suggests that these correlations emerge from
spatially broad interactions and not from local common input.

View larger version (25K):
[in this window]
[in a new window]
|
Figure 5.
Mean significant correlation coefficient values
(averaged over each discrete distance) between the responses of two MI
neurons as a function of the distance between the electrodes on which
the neurons were recorded. No significant relationship between the
strength of the correlation and the interelectrode spacing was
observed. Error bars indicate SDs. Note that, unlike Figure 4, mean
correlation coefficient values are based on only those that were
significantly different from zero. When all correlation coefficient
values are considered, the relationship between the mean correlation
and electrode distance remains flat and nonsignificant.
|
|
Directional variations in correlation
Although correlated activity between cortical neurons has been
previously observed, we discovered that correlation coefficients between pairs of MI neurons could vary with movement direction [see
Hatsopoulos et al. (1998b)
, for directional changes in synchrony measured on a fine time scale]. Notice how the same cell pair exhibits
strong correlated activity on a trial-by-trial basis for one movement
direction (Fig. 6A,
top) but much weaker correlated activity for a different movement
direction (Fig. 6A, bottom). Figure
6B provides an example of directional variations in
correlation between two MI neurons over all eight directions. As is
evident from the firing rate tuning curves of the two cells, they have similar preferred directions of ~135°; yet they engage in
significant correlated firing for movements in the 45 and 270°
directions. A one-way ANOVA for direction revealed that 78%
(387/497) of the cell pairs exhibited significant differences in
correlation with movement direction (p < 0.05).
This finding indicates that correlation strength is not necessarily a
fixed property of neuronal pairs but rather can depend on movement
direction.

View larger version (20K):
[in this window]
[in a new window]
|
Figure 6.
Variations in correlated activity with movement
direction. A, A cell pair showing strong trial-by-trial
covariation in firing for movements in the 135° direction
(left) but much weaker covariation for movements in the
315° direction. The normalized firing rates of one cell (open
squares and dashed line) are plotted together
with those of another cell (solid circles and
solid line) over all trials in the recording session.
Whereas the normalized firing rates are significantly correlated for
135° movements (rij = 0.85;
p < 0.001), they are not significantly correlated
for 315° movements (rij = 0.05;
p < 0.79). B, Variations in the
firing rates of two cells (top) and in their correlation
coefficient (bottom) over eight movement directions. The
two cells fire maximally in the 135° direction, whereas their peak
correlation occurs in the 45 and 270° directions.
|
|
Covarying neural discharge and directional coding
We evaluated whether correlations could contribute to directional
coding by comparing the performance of statistical models that
incorporate second-order relations with those that assume independence
of neurons in the population. These statistical models were used to
classify single trials of neuronal data into one of eight classes
corresponding to each of the eight movement directions. Bayesian
estimation techniques (Sanger, 1996
; Oram et al., 1998
; Zhang et al.,
1998
) as well as approximations to Bayesian estimation were used to
determine the direction of the arm movement. Bayesian techniques use
probabilistic descriptions of the relationship between the firing rate
of a neuron and the direction of arm movement to return the probability
of a movement direction given the firing characteristics of a
population of cells within a trial. This probabilistic description was
constructed by assuming that the number of action potentials across all
the neurons in the population,
for movements in a particular
direction, d, followed a multivariate Gaussian distribution:
|
(4)
|
where N is the number of neurons in the population,
is the mean
population response, and Cd is
the covariance matrix for movements in direction, d. The
vertical lines surrounding the covariance matrix denote its
determinant. The superscript T denotes the transpose
operator. The population response on a given trial,
, and the mean population
response,
d,
were obtained directly from the experimental data.
To determine the direction of the arm movement from the response of the
ensemble of MI neurons, Bayes' rule was used to transform the
probability of observing a response in a particular movement direction,
P(
|d), into a probability that a movement was performed in a particular direction given a particular pattern of neural firing (i.e., the neural
"response"):
|
(5)
|
where P(d) is the prior probability of a particular
movement direction and
P(
) is the prior
probability of a population firing pattern.
For each trial, these eight probabilities were collapsed into a single
estimate for the direction of movement by choosing the direction with
the highest probability of occurrence (i.e., the maximum a
posteriori estimate or MAP estimator). If the direction chosen was
the same as the actual direction in which the monkey moved, then the
estimation was successful (chance, 12.5%). Because of the limited data
available for these analyses, cross-validation of the model was
performed by removing each trial from the data set, forming the various
probability distributions on the remaining data, and estimating the
direction of movement for the firing rates in the removed trial.
The role that the response correlations might play in directional
coding was tested by comparing two versions of the Gaussian model
defined in Equation 4. The first model assumed that the responses of
each of the neurons were conditionally independent. This was
accomplished by setting the off-diagonal terms of the covariance
matrix, Cd, to zero and only
estimating on-diagonal terms, i.e., the variances of each of the
neurons. In other words, the responses of all cells were assumed to be
uncorrelated with each other and, therefore, independent of each other
given that they are normally distributed. It is important to stress
that this is a conditional independence assumption. That is,
conditional on a particular movement direction, the probability that
one neuron gives a particular response is independent of the
probability that any other neuron gives a particular response. Clearly,
cells with similar preferred directions are unconditionally dependent.
The second model, on the other hand, incorporated the observed
second-order correlations by estimating the off-diagonal elements of
the covariance matrix. To illustrate the difference between the two
models, data based on the activity of two neurons were classified into
one of two possible movement directions using either the independent
model (Fig. 7A) or the second-order model (Fig. 7B). Trials from one data set are
correctly classified (solid symbols) or incorrectly classified (hollow
symbols) as 180° (circles) or 270° (triangles) movements. In
addition, contour plots show the isoprobability lines representing the
joint probability of activity of two cells conditional on either the 180° direction (dashed lines) or 270° direction (sold lines). Notice how the second-order model (Fig. 7B) fits the data
better resulting in fewer incorrectly classified trials.

View larger version (27K):
[in this window]
[in a new window]
|
Figure 7.
Classification of trials from one data set into
one of two movement direction categories using a pair of simultaneously
recorded neurons. Trials are classified correctly
(solid) or incorrectly (hollow) as either
180° movements (circles) or 270° movements
(triangles). Isoprobability lines represent the joint
probability of activity in the two cells conditional on either 180°
movements (dashed lines) or 270° movements
(solid lines). A, Trials are classified
using the independent model in which the two cells are assumed to be
uncorrelated. B, Trials are classified using the
second-order model in which the covariances between the cells are
incorporated into the model.
|
|
The classification results for the independent and second-order models
based on the complete ensembles of recorded cells (instead of just
pairs of neurons) are shown in Figure
8A. The shaded bars in
the figure show the success of the model that treated the neurons as
independent encoders of movement direction in determining the direction
of movement from a population of simultaneously recorded MI neurons.
The dark bars show the success of the second-order model that
incorporated the pairwise covariances in predicting the direction of
movement from a simultaneously recorded population of MI neurons. In
all cases, incorporating the second-order interactions between neurons
improved classification of movement direction (Fig.
8A). On average, incorporating the covariance
structure resulted in correct determination of movement direction in
72% ± 17% (±1 SD) of the trials; the model without the correlations correctly estimated the direction of movement in 61% ± 12% of the
trials. The difference between these two populations (11%) was
significant (paired t test, t(4) = 3.97, one-tailed, p < 0.05); this test demonstrates that
these correlations provide additional information about direction not
present in the firing rates of individual neurons treated as
independent coding elements. In both cases, the models performed
significantly better than the chance rate of 12.5%.

View larger version (31K):
[in this window]
[in a new window]
|
Figure 8.
A, Results of applying two models
to decode the direction of movement from the population of neural
responses. Dark bars (second order) are the percentage
of trials in which a model incorporating correlations in the
trial-to-trial responses of pairs of neurons correctly estimated the
movement direction. Hatched bars (independent) are the
numbers of trials for which a model that did not incorporate the
correlations between cells correctly predicted the direction of
movement. These models were applied to five experimental data sets
denoted by a seven -numeric code. B, A comparison in
performance of two versions of the second-order model. The improved
classification performance of the original second-order model as
compared to the independent model is labeled as the "variable
covariance" model (black bars). In the other version,
the covariance matrix used was fixed and was set to the average over
the eight movement directions (gray bars). Notice
how the original second-order model does consistently better than the
fixed covariance model over all data sets.
|
|
We next investigated to what extent the observed changes in correlation
with direction could account for the 11% improvement in the
performance of the second-order model. If we treat the correlation
strength as composed of a fixed DC component and a variable AC
component, the bulk of the improvement in directional estimation may be
accounted for by incorporating the fixed correlational structure into
the second-order model. On the other hand, the observed changes in
correlation with movement direction may be contributing substantially
to the improvement of the second-order model. To tease apart these two
contributions, we compared the performance of the original second-order
model having variable correlations to a second-order model whose
covariance matrices were fixed to the average of their values over all
eight movement directions. Figure 8B shows that the
original second-order model does consistently better than the model
with a fixed covariance structure. The model with fixed covariance
performed only 5% better than the independent model compared to an
11% improvement with the original second-order model. Thus, 55% of
the improvement [(11
5)/11] in performance of the
second-order model is because of variations in correlation with
movement direction.
The nature of information contained in covarying
neural discharge
Correlations between the responses of neurons can either encode
information redundantly, independently, or synergistically (Gawne and
Richmond, 1993
). When information is redundantly encoded in a pair of
neurons, there is less information in the simultaneous firing rates of
the two cells than if the two neurons independently represented
information (i.e., if one neuron indicates a leftward direction with
100% reliability, having a second neuron that indicates left is
redundant). Likewise, when two neurons synergistically encode
information, the simultaneously recorded firing rates of the two
neurons contain more information than if the neurons were taken
independently (i.e., the coactivation of a left direction and a right
direction neuron may indicate a new movement direction). Because
accurate estimates of information content in rate codes requires more
data than we were able to collect (Golomb et al., 1997
), we used
statistical models that were constrained by the experimental data to
determine the nature of information available from pairs of neurons.
The model was constructed using (1) the cosine directional tuning
functions for neurons whose firing rate modulation was well-described by a cosine, (2) a model for the dependence of the response variability of a single cell on its mean firing rate (Fig. 2), and (3) the observed
correlations in the responses of simultaneously recorded cell pairs.
These three parameters were derived from the data and then used to
create discrete bivariate Gaussian models,
P(r|d), for pairs of neurons as well univariate
Gaussian models for individual neurons. These models were similar to
those used in the Bayesian classification analysis described above
except that they used discrete instead of continuous normal
distributions. Using Equations 5 and 6, it was possible to compute the
average mutual information I(D;R)
between the joint response of the pair of neurons and the direction of
the arm movement. It was compared to the sum of the information from
the two constituent neurons (Gawne and Richmond, 1993
). If the
correlations represent directional information redundantly, the ratio
of the joint mutual information to the sum of the mutual information
from the two constituent neurons treated separately is <1.0. If the
ratio is >1.0, the correlations between the simultaneously recorded
neurons contribute more information than is available from the two
cells treated independently.
|
(6)
|
where the angled brackets represent an average over neuronal
responses. Figure 9 shows the
distribution of the information fraction for 342 pairs of neurons. This
distribution indicates that covarying discharges in a pair of MI
neurons can synergistically encode information about movement
direction, although for many cases the coding is redundant. The
significant shift in the mean of this distribution (1.064 ± 0.006; p < 0.001; Wilcoxon signed-rank test) suggests
that that correlations in the trial-to-trial response variability of MI
pairs as a group are more likely to increase than decrease information
content. In some cases, the presence of correlations can increase the
representational capacity of the pair of neurons in excess of 30% over
the independent case.

View larger version (33K):
[in this window]
[in a new window]
|
Figure 9.
Distribution of the information fraction for pairs
of MI neurons. The information fraction measures whether the
simultaneously recorded responses of a pair of neurons contain more
information than the sum of the information extracted from the two
neurons treated separately. Fractions >1.0 indicate that the joint
response of the neurons contains more information than the two neurons
taken independently. Fractions <1.0 indicate redundancy in the
representation of information in the pair of neurons.
|
|
 |
DISCUSSION |
The present study demonstrates the capability of neuronal
populations to represent information in their joint activity. The second-order covariations in firing between MI neurons provide directional information above that found in single neuron firing rates.
Directional classification of data on a trial-to-trial basis was
consistently better when these interactions were considered. Information-theoretic analysis demonstrated that the joint activity of
certain neuron pairs could provide >30% additional information about
direction (Fig. 9). These findings suggest that higher-order ensemble
codes that incorporate the interactions of a population of neurons can
make important contributions to cortical representations [see
Hatsopoulos et al. (1998a)
for a comparison between ensemble and
standard population codes] if these codes can be exploited. We have
not demonstrated that these codes are actually used but that additional
information is available when these interactions are considered.
Noise in single neuron firing rates
Although it has been shown that single MI neurons carry
information about movement direction in their average firing rate (Georgopoulos et al., 1982
), the firing rate variability is higher than
would be expected from a random Poisson process [Fig. 2; see also Fig.
2 in Lee et al. (1998)
]. However, we know that precisely timed and
controlled movements can be generated despite this variability. To
overcome the problem of random noise, precise movements could be
represented by averaging the responses of a population of MI neurons
(Georgopoulos et al., 1986
). When noise sources are statistically independent across different neurons, the effect of noise on the population average decreases by a factor of one over the square root of
the number of neurons, thereby increasing the signal-to-noise ratio.
This explanation becomes less tenable as noise becomes more correlated
across cells (Zohary et al., 1994
). We have shown that this noise is
partially correlated, in agreement with others (van Kan et al., 1985
;
Gawne and Richmond, 1993
; Zohary et al., 1994
; Lee et al., 1998
), but
unlike these studies we have emphasized that this correlated noise is
not fixed but rather varies with behavioral condition, i.e., movement direction.
Our approach considered correlated noise as part of the signal and
demonstrated that it carried additional information about movement
direction above that obtained from the firing rates of individual
neurons acting as independent encoders. Thus, rather than limiting the
information capabilities of neural populations, correlations may be
able to extend them as an additional type of signal. A key feature of
this approach is the finding that the firing covariation between
neurons varies across different movement directions. As Figure
8B demonstrates, a major contribution to the
improvement in classification performance of the second-order model was
because of the fact that the directional variations in firing
covariance among motor cortical neurons were incorporated into the
model. That is, the firing covariances were preserved for each movement
direction unlike other studies, which pooled the covariances over
different stimulus or movement directions (Zohary et al., 1994
; Lee et
al., 1998
).
These findings complement other results from MI that demonstrate that
information can be obtained from correlated activity between neurons on
a finer time scale. Riehle et al. (1997)
have shown that synchronous
pairwise firing on a millisecond time scale can signal the expectation
of visual cues instructing the monkey to initiate a movement in the
absence of rate modulation. In addition, results from our laboratory
have demonstrated that directional information is available in fine
synchronous interactions of neurons that cannot be accounted for by
independent rate modulation (Hatsopoulos et al., 1998b
). We have not
fully examined the relationship, if any, between broad correlations in
spike count measured on the time scale of hundreds of milliseconds and
fine synchrony measured on the millisecond time scale. However,
preliminary results suggest that they are relatively independent
phenomena. The flat relationship between broad correlation and
interelectrode distance (Fig. 5) is inconsistent with the relatively
strong dependence of fine synchrony with distance that has been
observed previously (Kwan et al., 1987
; Hatsopoulos et al., 1998b
,
their Fig. 1D). This also suggests that the two
processes may be independent.
Our results suggest that the observed firing variability may be
considered a potential signal when viewed in the context of an ensemble
of neurons. This may be true not only in MI but throughout cortex.
Studies in the hippocampus and somatosensory cortex using different
ensemble-based approaches have also found that classification of
behavior or sensory localization improves when higher dimensional features of a neural ensemble are incorporated (Deadwyler et al., 1996
;
Nicolelis et al., 1998
).
Source of correlation
Although correlated activity between cortical neurons has been
widely detected, the source of these interactions has not been established. One possibility is that correlations are a consequence of
common input from other neurons. Nearby neurons, which have been the
focus of most recording studies, might be expected to be correlated
because they would share thalamic projections (in MI, Shinoda and
Kakei, 1989
). However, the bulk of our neurons would not likely be
driven by common thalamic sources because our electrodes were spaced by
at least 400 µm, which is a considerable distance for two cells to
receive strong, common input. Moreover, the observation that the
correlation strength does not vary with distance argues further that
common input is not a major source of correlation.
An alternative hypothesis is that correlations arise from large-scale
reciprocal connections among cortical neurons (Juergens and
Eckhorn,1997
). Neurons in motor cortex form a network of connections with each other with a high degree of lateral interconnectivity (Ichikawa et al., 1985
; Huntley and Jones, 1991
; Hess and Donoghue, 1994
). Moreover, these lateral connections can become stronger with the
acquisition of new motor behaviors (Rioult-Pedotti et al., 1998
).
Functional studies have also shown that a large, distributed area of
the MI arm area activates when simple reaching movements are made,
suggesting that motor cortex is not organized as a set of independent
columns (Sanes et al., 1995
).
Functional role of correlations
It is unclear what functional role these correlations may be
playing, if any. One possibility is that they are the consequence of
anatomical connections between cortical neurons that may be important
in promoting population codes. Recent theoretical work has suggested
that lateral connections could be useful in refining a coarse code
among a population of cortical neurons (Pouget et al., 1998
). A
particular arrangement of lateral connections between directionally
selective neurons would reduce the effects of neural noise by creating
a smooth "hill" of activity over the set of cortical neurons.
According to the neural model they develop, the correlation between two
directionally selective neurons is predicted to be proportional to the
product of the derivative of their tuning curves. Although this
prediction has not been systematically tested with our data, the
variations in correlation demonstrated in Figure 6B
are somewhat consistent with it.
Another hypothesis specific to motor cortex is that these correlations
are used to create temporary synergies among muscle groups. Layer 5 pyramidal neurons in primary motor cortex form part of the
corticospinal tract and, thus, directly and indirectly synapse onto
motor neuron pools that innervate skeletal muscles. Hepp-Raymond et al.
(1996)
have recently made electromyographic recordings from hand
muscles that are active during precision grip movements. It was
observed that particular muscles became correlated and, more
importantly, the correlation strength varied with the overall level of
force. These variable muscle synergies may be created by linking
different neurons in motor cortex.
Methodological limitations
This work has demonstrated that the presence of second-order
interactions in a population of MI neurons can significantly add to the
representation of directional information. This is an initial step
toward examining higher-order relationships among ensembles of neurons.
We restricted ourselves to a relatively simple multivariate model
(i.e., the Gaussian model) that considered only pairwise correlations
between neurons. This was, in part, because of the limited amount of
data that could be acquired from a behaving animal and is a major
challenge in the endeavor of discovering neural ensemble
representations. It is unknown, however, whether even higher-order
moments may account for more of the response variability and,
therefore, provide more information about direction or any other
parameter. As multielectrode recording becomes more prevalent, novel
statistical methods will be needed that are mathematically sound and
tractable, and which incorporate multidimensional interactions.
Although we have shown that higher-order representations provide more
information about movement direction, these experiments do not
demonstrate that these higher-order representations are actually being
used by the animal to control motor behavior. In fact, most if not all
studies that postulate a coding scheme based on the firing rates of
single neurons have the same problem. A recent study of olfactory
coding in an invertebrate preparation has attempted to show that
correlated neural activity as revealed by field potential oscillations
is necessary for fine olfactory discrimination (Stopfer et al., 1997
).
This was accomplished by pharmacologically desynchronizing a population
of neurons without affecting their individual firing properties and
showing that olfactory discrimination was impaired. In the future,
similar perturbation experiments could be done in the cortex of mammals to determine whether higher-order interactions among neurons are necessary to guide behavior.
 |
FOOTNOTES |
Received Jan. 12, 1999; revised June 28, 1999; accepted June 28, 1999.
This work was supported by National Institutes of Health Grants NS25074
and NS10269-02, National Science Foundation Grant SBR-N720368, the Utah
Center of Excellence in Neural Interfaces, and the W. M. Keck
Foundation. E.M.M., N.G.H., and J.P.D. are responsible for the concept,
analysis, and writing of the this manuscript. All other authors were
involved in the surgical and technical aspects of this work, as well as
in the data acquisition. We thank Elie Bienenstock and Stu Geman for
helpful comments on this manuscript. We also thank Beth Travers and
Gerhard Friehs for their technical assistance.
The first two authors are the primary authors of this paper.
Correspondence should be addressed to Nicholas Hatsopoulos, Box 1953, Department of Neuroscience, Brown University, Providence, RI 02912.
 |
APPENDIX |
The probability of observing a number of action potentials across
an ensemble of N motor cortical neurons,
, given the monkey's movement
direction, d, was modeled as a multivariate Gaussian
distribution:
|
(4)
|
Eight of these distributions were specified, one for each
direction. The parameters that fully specify each distribution are: (1)
the mean spike count over the ensemble of neurons for a particular
movement direction,
d. This
vector of mean values is N-dimensional; (2) the covariance matrix, Cd, computed for a
particular movement direction. This matrix is N × N.
An element of the covariance matrix is computed as follows:
where Md is the number of
trials for direction d. The inverse of the covariance matrix
is defined as follows:
CdCd
1
where I is the identity matrix. The vertical lines in the
denominator of the above equation represents the determinant of the matrix.
The mutual information, I(D;R), between the
movement direction, D, and the response of the neural
ensemble, R, is defined as follows:
|
(6)
|
The a priori probability of a movement direction,
P(d), is determined by the experimenter and was specified to
be uniform. The a posteriori probability of a movement
direction given the response of the ensemble,
P(d|
), was estimated by using Bayes' rule (Eq. 5). This equation can be
decomposed into three parts: (1) the logarithm of the ratio of the
a posteriori to a priori probability
distributions (the term in the brackets). If the firing rates of the
ensemble do not alter one's estimation of the probability of a
particular direction, then the ensemble response doesn't provide
information about direction. In this case, the ratio of distributions
in Equation 6 will be 1, and the log of 1 is zero, which confirms our
intuition. On the other hand, if the ensemble response does alter
one's estimation of the probability of a particular direction, then
the ratio is not 1. If the a posteriori probability is
large, the ratio will be >1, which will give a positive log. (2) The
expectation or average of the information over all movement directions
conditional on the ensemble response (the sum sign and its terms). The
information is weighted more for directions that are conditionally more
likely. (3) The average over all ensemble responses (is represented
with the angular brackets).
 |
REFERENCES |
-
Abbott LF,
Dayan P
(1998)
The effect of correlated variability on the accuracy of a population code.
Neural Comput
11:91-101[Abstract/Free Full Text].
-
Bair W,
Koch C
(1996)
Temporal precision of spike trains in extrastriate cortex of the behaving macaque monkey.
Neural Comput
8:1185-1202[Abstract].
-
Buracas GT,
Zador AM,
DeWeese MR,
Albright TD
(1998)
Efficient discrimination of temporal patterns by motion-sensitive neurons in primate visual cortex.
Neuron
20:959-969[ISI][Medline].
-
Crow EL,
Davis FA,
Maxfield MW
(1960)
In: Statistics Manual. New York: Dover.
-
Deadwyler SA,
Bunn T,
Hampson RE
(1996)
Hippocampal ensemble activity during spatial delayed-nonmatch-to-sample performance in rats.
J Neurosci
16:354-372[Abstract/Free Full Text].
-
Dean AF
(1981)
The variability of discharge of simple cells in the cat striate cortex.
Exp Brain Res
44:437-440[ISI][Medline].
-
DeCharms RC,
Merzenich MM
(1996)
Primary cortical representation of sounds by the coordination of action-potential timing.
Nature
381:610-613[Medline].
-
Donoghue JP,
Sanes JN,
Hatsopoulos NG,
Gaal G
(1998)
Neural discharge and local field potential oscillations in primate motor cortex during voluntary movements.
J Neurophysiol
79:159-173[Abstract/Free Full Text].
-
Duda RO,
Hart PE
(1973)
In: Pattern classification and scene analysis. New York: Wiley & Sons.
-
Fu QG,
Suarez JI,
Ebner TJ
(1993)
Neuronal specification of direction and distance during reaching movements in the superior precentral premotor areas and primary motor cortex of monkeys.
J Neurophysiol
70:2097-2116[Abstract/Free Full Text].
-
Gawne TJ,
Richmond BJ
(1993)
How independent are the messages carried by adjacent inferior temporal cortical neurons?
J Neurosci
13:2758-2771[Abstract].
-
Georgopoulos AP,
Kalaska JF,
Caminiti R,
Massy JT
(1982)
On the relations between the direction of two-dimensional arm movements and cell discharge in primate motor cortex.
J Neurosci
2:1527-1537[Abstract].
-
Georgopoulos AP,
Schwartz AB,
Kettner RE
(1986)
Neuronal population coding of movement direction.
Science
233:1416-1419[Abstract/Free Full Text].
-
Georgopoulos AP,
Kettner RE,
Schwartz AB
(1988)
Primate motor cortex and free arm movements to visual targets in three-dimensional space. II. Coding of the direction of movement by a neuronal population.
J Neurosci
8:2928-2937[Abstract].
-
Georgopoulos AP,
Schwartz AB,
Kettner RE
(1992)
Neuronal population coding of movement direction.
Science
233:1416-1419.
-
Golomb D,
Hertz J,
Panzeri S,
Treves A,
Richmond B
(1997)
How well can we estimate the information carried in neuronal responses from limited samples?
Neural Comput
9:649-665[Abstract].
-
Hatsopoulos NG,
Ojakangas CL,
Maynard EM,
Donoghue JP
(1998a)
Detection and identification of ensemble codes in motor cortex.
In: Neuronal ensembles: strategies for recording and decoding (Eichenbaum H,
Davis J,
eds), pp 161-175. New York: Wiley.
-
Hatsopoulos NG,
Ojakangas CL,
Paninski L,
Donoghue JP
(1998b)
Information about movement direction obtained from synchronous activity of motor cortical neurons.
Proc Natl Acad Sci USA
95:15706-15711[Abstract/Free Full Text].
-
Heggelund P,
Albus K
(1978)
Response variability and orientation discrimination of single cells in striate cortex of cat.
Exp Brain Res
32:197-211[ISI][Medline].
-
Henry GH,
Bishop PO,
Tupper RM,
Dreher B
(1973)
Orientation specificity and response variability of cells in the striate cortex.
Vision Res
1