Abstract
Precise spatiotemporal patterns in neural discharge are a possible mechanism for information encoding in the brain. Previous studies have found that such patterns repeat and appear to relate to key behavioral events. Whether these patterns occur above chance levels remains controversial. To address this question, we have made simultaneous recordings from between two and nine neurons in the primary motor cortex and supplementary motor area of three monkeys while they performed a precision grip task. Out of a total of 67 neurons, 46 were antidromically identified as pyramidal tract neurons. Sections of recordings 60 s long were searched for patterns involving three or more spikes that repeated at least twice. The allowed jitter for pattern repetition was 3 ms, and the pattern length was limited to 192 ms. In all 11 recordings analyzed, large numbers of repeating patterns were found. To assess the expected chance level of patterns, “surrogate” datasets were generated. These had the same momentbymoment modulation in firing rate as the experimental spike trains, and matched their interspike interval distribution, but did not preserve the precise timing of individual spikes. The number of repeating patterns in 10 randomly generated surrogates was used to form 99% confidence limits on the repeating pattern count expected by chance. There was close agreement between these confidence limits and the number of patterns seen in the experimental data. Analysis of high complexity patterns was carried out in four long recordings (mean duration 23.2 min, mean number of neurons simultaneously recorded 7.5). This analysis logged only patterns composed of a larger number (7–11) of spikes. The number of patterns seen in the surrogate datasets showed a small but significant excess over those seen in the original experimental data; this is discussed in the context of surrogate generation. The occurrence of repeating patterns in the experimental data were strongly associated with particular phases of the precision grip task; however, a similar task dependence was seen for the surrogate data. When a repeating pattern was used as a template to find inexact matches, in which up to half of the component spikes could be missing, similar numbers of matches were found in experimental and surrogate data, and the time of occurrence of such matches showed the same task dependence. We conclude that the existence of precise repeating patterns in our data are not due to cortical mechanisms that favor this form of coding, since as many, if not more, patterns are produced by spike trains constructed only to modulate their firing rate in the same way as the experimental data, and to match the interspike interval histograms. The task dependence of pattern occurrence is explicable as an artifact of the modulation of neural firing rate. The consequences for theories of temporal coding in the cortex are discussed.
INTRODUCTION
The nature of the code used by the brain to store and process information is an important issue in neuroscience. Conventionally, it has been assumed that the rate at which single neurons discharge action potentials is the relevant code, but in the last 10 years, there has been considerable interest in “temporal” codes, in which the precise timing of individual action potentials is the important variable.
Several different variants of temporal coding have been proposed. Most assume that synchronous firing of multiple neurons carries information (Gerstein et al. 1989). However, the “jitter” in the firing allowed before discharges are no longer considered to be synchronous varies between different models. In the cerebral cortex, many workers have provided evidence for peaks in the crosscorrelation histogram between cells that are around 20 ms wide (Fetz et al. 1990; Vaadia et al. 1995; S. N. Baker, R. Spinks, A. Jackson, and R. N. Lemon, unpublished observations); these may be accompanied by side lobes that indicate an oscillatory nature to the synchrony, at frequencies close to 40 Hz in the visual system (Eckhorn 1994; Singer and Gray 1995) or 25 Hz in the motor system (Donoghue et al. 1998;Murthy and Fetz 1996; Baker et al., unpublished observations; E. M. Pinches, S. N. Baker, and S. N. Lemon, unpublished observations).
Other work on temporal coding has emphasized synchrony at considerably higher precision, with around 1 ms allowed jitter. The existence of such precise synchronization is predicted by the “synfire chain” theory of Abeles (1991). The technique of unitary event analysis (Grün 1996) is designed to investigate synchrony between pairs of cells on this time scale, and some evidence for its existence in experimental recordings has been obtained (Grammont and Riehle 1999; Riehle et al. 1997; but see Pauluis and Baker 2000). However, an even stronger prediction of the synfire chain theory is that spike trains recorded from multiple neurons simultaneously will show patterns in their discharges. Such patterns are postulated to span relatively large time intervals (hundreds of milliseconds), but the location of spikes within them is precise to within a few milliseconds. An example of such a pattern might be a spike from neuron 1, followed 57 ms later by one from neuron 2, and 124 ms later by a spike from neuron 3. Since such patterns are defined both across time, and across multiple neurons (space), they are referred to as spatiotemporal patterns. The synfire chain theory proposes that a given pattern will be generated each time the brain enters a particular state. If these states are repeated, for example during repetitive performance of the same behavioral task, the pattern should also recur. The existence of repeating spatiotemporal patterns at above chance levels would therefore be an important piece of evidence that this type of temporal coding is used by the brain.
A number of previous studies have claimed to find such patterns in experimental data (Abeles et al. 1993; Dayhoff and Gerstein 1983; Frostig et al. 1990;Prut et al. 1998; Strehler and Lestienne 1986; Villa and Abeles 1990; Villa et al. 1999). However, it is a challenging statistical problem to determine whether the number of repeating patterns seen is in excess of those expected by chance. As pointed out by Abeles and Gerstein (1988), the number of patterns expected can depend on a high power of the firing rates of the neurons analyzed, such that it will be very sensitive to changes in firing rate over the analysis period. Most studies to date have assumed a constant firing rate in assessing the chance pattern expectancy.
Recently, Oram et al. (1999) presented an analysis of patterns in spikes recorded from the visual cortex, and the lateral geniculate nucleus (LGN), which showed that if both the stimulusinduced modulation in firing rate and the nonPoisson variability in the number of discharges per stimulus was accounted for, the number of patterns seen was close to that expected by chance. However, the work of Oram et al. (1999) had a number of limitations compared with previous studies that had found the converse result. Only patterns generated by a single neuron were investigated, and the total duration of the patterns was limited to 50 or 75 ms, compared with the few hundred milliseconds often used previously. TheOram et al. (1999) study investigated only triplet or quadruplet patterns, and the length of data used for analysis was short (400 ms, compared with around 1 min in most previous studies). Finally, the section of data used to search for repeating patterns was that following a single visual stimulus. If repeating patterns are encoding stimulusspecific information, it might be expected that they would repeat at a similar poststimulus latency with successive stimuli; in this case, excess repeating patterns would not be expected when analyzing only single trials, but only when the responses to multiple stimuli are compared.
This study presents an analysis of the number of repeating patterns observed in multiple singleunit recordings made from the primary motor cortex (M1) and the supplementary motor area (SMA) in awake behaving macaque monkeys performing a precision grip task. We show that both high and low complexity patterns can be seen in these data, involving single or multiple neurons, and that patterns repeat at preferred points of the behavioral task. However, the number of such patterns does not exceed those seen in “surrogate” data, constructed to match the firing rate modulation and interspike interval distribution of the experimental data. We conclude that there is no evidence in our data that precise repeating spatiotemporal patterns are used by the cortex to encode information.
METHODS
Recordings
Multiple singleunit recordings were gathered from three adult purposebred female macaque monkeys (monkey 29, m. nemestrina; monkeys 30 and 33, m. mulatta), according to methods described in full in previous publications (Baker et al. 1999). Briefly, the animals were trained to perform the precision grip task of Lemon et al. (1986) for a food reward. This required the monkey to squeeze two springloaded levers between finger and thumb into electronically defined position windows, hold for a period of around 1 s, and then release. Performance on the task was cued to the animal by three auditory tones indicating that the levers were in target, that the animal had held for long enough, and finally that the levers had been released and a reward would follow. The onset of the middle tone (End Hold) was used for alignment of some of the data analysis. Once the monkeys were fully trained, they were implanted with a stainless steel headpiece for head fixation (see Lemon 1984). Two electrodes were also chronically implanted in the pyramidal tract at the level of the medulla to permit antidromic identification of pyramidal tract neurons (PTNs) (see Baker et al. 1999). After a period of retraining with head fixation, a stainless steel recording chamber was mounted over a craniotomy made over M1 or SMA. All surgeries were performed under general anesthesia (induced by ketamine hydrochloride 10 mg/kg im and maintained by inhalation of Isoflurane in 50:50 O_{2}:N_{2}O) and aseptic conditions, and were followed by a full course of antibiotic (Terramycin/LA 20 mg/kg, Pfizer) and analgesic (Vetergesic 10 μg/kg, Reckitt and Colman Products) treatment. Following recovery from the chamber implant surgery, recordings were made from the cortex using an Eckhorn microdrive (Thomas Recording, Marburg, Germany) (Baker et al. 1999; Eckhorn and Thomas 1993). This allowed up to 16 independently moveable microelectrodes to be placed in the cortex on a 4 × 4 grid (interelectrode spacing 305 μm). Great care was taken to ensure that recordings were made only from single neurons, completely uncontaminated by spikes from other units. This included the use of collision tests to establish the identity of PTNs online (seeBaker et al. 1999). Continuous waveform data from each electrode were recorded on tape (24kHz sampling rate) for later offline discrimination into the times of singleunit discharge using custom written “clustercutting” software (Eggermont 1990). Discriminated units that contained interspike intervals shorter than 1 ms were discarded, as this is usually a sign of poor unit discrimination.
After sufficient recordings had been made from one cortical area (e.g., M1) underlying the initial craniotomy, the recording chamber was reimplanted under general anesthesia over a different area (e.g., SMA), or to the same area of the other hemisphere, and further recordings were obtained. The animal was retrained to perform the task with the other hand in cases where the chamber was moved to the opposite side, so that task performance was always with the contralateral hand. Following the completion of recordings, the monkeys were killed by an overdose of anesthetic and perfused through the heart; the brain was then removed for histological confirmation of the pyramidal tract electrode location. All procedures were carried out in accordance with the UK Animals (Scientific Procedures) Act 1986.
In monkey 29, the Eckhorn system was not used; neurons were recorded using a simpler twoelectrode drive (interelectrode spacing 1.5 mm) (see Baker et al. 1997).
Analysis
Datasets were searched for all spatiotemporal patterns that were a maximum of 192 ms long, and repeated at least twice to a precision of 3 ms, using the algorithm of Abeles and Gerstein (1988). The precision of 3 ms was chosen according to the data of Abeles et al. (1993). The algorithm was run in two stages. The first stage detected all patterns that occurred twice; this ran on a PC computer. The second stage detected patterns that had been so found more than once (i.e., those that occurred 3, 4, etc., times). The large numbers of patterns involved here meant that this second stage of theAbeles and Gerstein (1988) algorithm was by far the most computationally intensive; it was accordingly run on 64 nodes of a parallel supercomputer (Hitachi SR2201) at the University of Cambridge High Performance Supercomputing Facility. Patterns were characterized by their complexity (i.e., the number of spikes forming the patterns), and by their occurrence (the exact number of times that the pattern was detected).
Application of the repeating pattern search to the full datasets was not feasible, since the very large numbers of patterns found would have exceeded our capacity to store and process them by many times. Two separate analyses were therefore carried out. One sought to investigate repeating patterns of low complexity (triplets and higher) and used sections of data 60 s long extracted from the available recordings, which were normally much longer. The section for analysis was chosen to include a period when the animal was performing the behavioral task. Eleven datasets were included in this analysis. The second used the full available dataset, but did not search for repeating patterns that had a complexity lower than a preset limit. The level of this complexity threshold was set individually for each dataset to keep the files listing patterns to manageable size (<40 MB), and ranged from 7 to 11. Four datasets were included in the high complexity analysis.
Significance assessment
The method used to determine the number of repeating patterns that would be expected by chance is illustrated in Fig.1. Beginning from an experimentally recorded spike train (Fig. 1 A), the interspike interval (ISI) histogram was constructed (Fig. 1 F, solid line), and the modal interspike interval determined. The instantaneous firing rate of the neuron was then estimated by convolving the spike train with a Gaussian kernel, whose standard deviation was set equal to the modal interspike interval (Fig. 1 B) (Silverman 1986). This provided a representation of changes in firing rate on a momentbymoment basis during task performance. The neuron illustrated in Fig. 1 showed a higher firing rate just before movement onset [see finger (F) and thumb (Th) lever displacement signals in Fig. 1 C], a lower rate during the hold period (the arrowhead indicates the End Hold signal), and then again at a higher rate during the removal of the digits from the manipulandum.
A spike train was then simulated as follows. First, an inhomogeneous Poisson process was generated, with rate equal to n times the determined instantaneous firing rate (Fig. 1 D). This spike train was nfold decimated, such that only everynth spike was used, with the remainder being deleted. The starting point for the decimation was chosen at random from the firstn spikes. This produced a surrogate spike train that had an instantaneous gamma distribution of ISIs (with order parameter n), but a firing rate modulation equal to that of the experimentally recorded spikes (Baker and Gerstein 2000). A gamma process is a better representation of a neuron spike train than a Poisson process, since it allows for the presence of a relatively refractory period (Kuffler et al. 1957; Stein 1965). The procedure was repeated for values ofn = 1–30, and the timeaveraged ISI histogram for each of these surrogates was determined.
Because the underlying firing rate used to simulate each of these surrogates was continuously modulated, none of the surrogate ISI histograms should resemble that expected for a gamma process; instead, they will be a complex mixture of the histograms expected at multiple different rates. However, if it is assumed that the experimental spike train can be modeled as a nonhomogeneous gamma process of constant order, then a similar mixing in the timeaveraged ISI histogram should occur for both experimental and surrogate spike trains, since they share the same rate modulation time course. To estimate the best gamma order to use, the surrogate histogram that fitted the experimental ISI distribution best in the leastsquares sense was determined (Fig.1 F, gray histogram), and this value of n was used to generate the surrogate spike trains for this cell that were used for all subsequent analysis (Fig. 1 E; n = 16 for this neuron).
The above procedure was repeated for all neurons within a given dataset. Ten surrogates for each cell were then simulated, using the same firing rate profile and order parameter n, but different values of the seed for the random number generator used to produce the inhomogenous Poisson processes. Each of the 10 surrogate datasets was then subjected to the same search for repeating patterns as the original data. The mean μ and standard deviation ς of the number of patterns at a given complexity and number of occurrences was determined from these 10 surrogates.
The problem now arises of how to place appropriate confidence limits on the number of patterns expected in the experimental data, given the results from the 10 surrogates. Within a particular complexity and number of occurrences, there are many different possible patterns, involving different neurons and different intervals between the component spikes. The number of each of these individual, identified patterns in a dataset should follow a Poisson distribution (seeAbeles and Gerstein 1988, Appendix), in which the variance equals the mean. When all of these different possible patterns are combined and counted together to obtain a total count for this complexity and number of occurrences, the distribution of the total will not simply be the convolution of the individual distributions, because the component pattern counts are not independent (for example, if fewer than average spikes from neuron number 1 occur, this will reduce the probability of all patterns involving that cell). We therefore should expect that the variance of the total pattern count will not in general be equal to the mean. However, since the total pattern count is the sum of a large number of random variables, then by the central limit theorem it should follow an approximately normal distribution. We therefore placed 99% confidence limits on the number of patterns expected as μ ± 2.58ς, which assumes a normal distribution and uses the mean and standard deviation determined empirically from the 10 surrogates simulated. This gave the range of the number of patterns that would be expected by chance from spike trains that simply showed the same modulation in firing rate and ISI distribution as the original data.
The surrogate spike trains not only matched the ISI distribution of the experimental spike trains, but also the serial correlation of successive ISIs. This is illustrated in Fig. 1 G, which shows the correlation coefficient between ISI_{i} and ISI_{i+j} for different interval lagsj for the same cell as used in the earlier parts of this figure. The black circles show the correlation for the experimental spike train, and the vertical extent of the boxes shows the correlation range over the 10 simulated surrogates. There is good agreement here, indicating that the serial ISI correlation arises from the modulation in firing rate of the neuron, which is mimicked by the surrogate data. It should be noted that we investigated different methods of estimating instantaneous firing rate. Methods based around the reciprocal of the ISI (e.g., Pauluis and Baker 2000) generated surrogates with ISI histograms close to the experimental distributions but often failed to produce a similar profile of serial ISI correlation as in the original spike data. The kernel estimate of instantaneous firing rate described above usually produced surrogates that matched the serial ISI correlation more accurately, and we therefore chose it to generate the surrogates for all of the analysis described here.
The experimental spike trains often showed peaks in their crosscorrelation histograms (see Baker et al. 2000 for further analysis of these). Since the surrogates were generated for one cell at a time, the discharges of pairs of surrogate spike trains had no shortterm synchronization (Sears and Stagg 1976). This is illustrated in Fig. 1 H, which shows a crosscorrelation peak between two experimental spike trains (thick line), which is not present in correlation between the two corresponding surrogate trains (thin line).
Figure 1 I shows the modulation in firing rate of the neuron featured in panels A–G aligned to the End Hold behavioral marker of the precision grip task, which signaled the moment when the animal had held within the desired target zone for long enough, and could release to obtain a reward. The thick line shows the profile of the experimentally recorded neuron and the thin line that for a surrogate. As expected, the taskrelated changes in firing rate are well matched; there is some slight smoothing in the surrogate compared with the experimental rate profile, resulting from the kernel estimator used to determine the instantaneous firing rate.
Thus the surrogate spike trains matched the instantaneous firing rate modulation and ISI statistics of the experimental spike trains, but did not reproduce the precise timing of the spikes. They therefore represent a “null hypothesis”; any patterns seen in the surrogates can be assumed to occur by chance, and not by any neural mechanism (e.g., synfire chains) that favors pattern generation.
It has been previously noted that certain pseudorandom number generators can produce an excessive number of patterns in surrogate data (Abeles and Gerstein 1988). All analyses reported here used the generator supplied with MATLAB version 5.1 (Mathworks), which has 35 seeds and a theoretical repeat cycle of 2^{1492} numbers. An analysis was also carried out using an alternative generator (“ran1”) (Press et al. 1989), and a CD ROM of numbers generated by combining the output of two “good” pseudorandom generators with physically produced noise (Marsaglia 1996); similar results were obtained.
RESULTS
The analysis described in this paper was based on 12 recording sessions in 3 animals; 3 sessions were recorded from SMA, the remainder from the hand representation of M1. Nine sessions came frommonkey 33, two from monkey 30, and one frommonkey 29. A total of 67 cells were recorded, of which 46 were antidromically identified as PTNs; between 2 and 9 neurons were simultaneously recorded (mean 5.8). An average of 25.7 min of data (range 6.6–62 min) were recorded for each session, with 1,173–91,333 spikes per neuron (mean 24,000).
Patterns of low complexity
To investigate repeating patterns of low complexity, 60slong sections of data were analyzed. Even in these short files, large numbers of repeating patterns were detected. Across the 11 datasets (3 from SMA) analyzed in this way, between 999 and 219,414 triplets were seen that occurred twice. The maximum occurrence of a triplet was 83, and the highest complexity seen was a pattern involving 23 spikes, which occurred twice.
Figure 2 A shows the number of patterns at different levels of complexity and number of occurrences seen in one recording, which included nine neurons (6 PTNs) recorded simultaneously. The black bars show the number of patterns seen in the experimental data. For clarity of display, the histograms only show a subset of the pattern categories that were detected. At the extremes of complexity and occurrence in this dataset, 1 triplet that occurred 32 times was found, and 1 pattern of complexity 23, which occurred twice.
The other bars in Fig. 2 A show the number of patterns seen in surrogate datasets constructed to reproduce different facets of the experimental data. The open bars show the mean number of patterns seen in the surrogates that reproduced the firing rate modulation of the data, and matched the ISI distribution by choosing the appropriate order of the gamma distribution (see methods); the error bars show 99% confidence limits on these. The remaining bars show the effect of forcing the gamma distribution to have ordern = 1 (thereby making it a Poisson process), or of using surrogates that had a constant firing rate equal to the mean experimental firing rate, without the momentbymoment rate modulation.
A very large number of repeating patterns was found (see ordinate scales in Fig. 2 A). In this 60s length of recording, there were a total of 10,688 spikes from 9 neurons. Across all complexities and occurrences, there were 1.11 × 10^{6}repeating patterns, which contained 1.08 × 10^{7} spikes. On average, therefore each spike participated in around 1,000 patterns.
Figure 2 A shows that the number of patterns expected by chance on the basis of the most realistic surrogates encompasses the numbers seen in the experimental dataset, so that there is no evidence here for an excess of repeating patterns above the chance level. At low complexity and number of occurrences, the differences between the different surrogates are small, but become more pronounced for the higher complexity and occurrence patterns. Of especial note is the failure of the constant firing rate surrogates to predict the number of patterns accurately: they tend to underestimate the number expected at high complexity, whereas for the triplets they produce an overestimate. Surrogates with Poisson distributed ISIs but variable firing rate give numbers of patterns close to the experimental values, but are not as accurate in predicting the number seen as the surrogates that match both ISI distribution and firing rate modulation.
In two cases in Fig. 2 A the number of patterns seen experimentally lay outside the limits generated from the surrogates. However, for each of these (complexity 4, 8 occurrences, and complexity 7, 3 occurrences), the number of patterns expected was very low (2 or 3). The assumption of a normal distribution for the pattern count is unlikely to be valid for such small expected numbers, since the count is constrained to be nonnegative. In subsequent quantitative analysis, therefore, combinations with 10 or fewer patterns expected were excluded. Over the 11 datasets analyzed, there were 480 combinations of complexity and number of occurrences with more than 10 patterns detected in the mean of the 10 surrogates. In seven instances the number of patterns seen experimentally was above the 99% confidence limits, and in nine it was below. The probability of 16/480 instances lying outside the 99% confidence limits if the null hypothesis is correct was computed as 1.0 × 10^{−} ^{5} (Binomial distribution). However, since on a similar number of occasions the number of patterns was below the confidence limits as above, it is unlikely that there is an excess of patterns in the experimental data compared with the surrogates. Instead, we suggest that the discrepancy arises from a slight departure of the pattern count from the normal distribution, which was assumed (see discussion).
Following this analysis, it remained possible that there was a weak trend in the experimental data toward an excess of patterns, which failed to reach significance when tested on a pointbypoint basis. To check for this, Fig. 2 B shows a scatter plot of the number of patterns seen in experimental versus surrogate data (the surrogate value is the mean of 10 different surrogates), for all 480 different combinations of complexity and number of occurrences with more than 10 expected patterns. The points lie very close to the identity line, and the regression line fitted to them has a slope not significantly different from one and intercept not significantly different from zero. This analysis therefore also fails to show an excess number of patterns in the experimental data.
Some of the neurons that we recorded had been antidromically identified as PTNs, so that it was possible to investigate the number of patterns produced by PTNs versus unidentified (UID) cells. Patterns were therefore categorized according to whether they involved only PTNs, only UIDs, or a mixture of both, and the numbers found at each complexity and occurrence compared between experimental and surrogate data. An alternative subdivision that we tried categorized the patterns according to whether they involved one or more than one neuron. The results of these five analyses were similar to those for the grouped data described above and provided no evidence for an excess of patterns in the experimental data in any category.
Patterns of high complexity
Patterns of high complexity were investigated by searching for patterns in the entire length of available data for four recording sessions (1 from SMA), all in monkey 33. These sessions lasted between 15.7 and 31.6 min. Figure3 A shows the number of patterns seen in different categories in one dataset, using a similar display to Fig. 2 A. Analysis here was confined to patterns of complexity 9 or higher to limit computation time and storage space. The histograms only show patterns with two occurrences; 19 examples of patterns with complexity 9 that occurred 3 times were also seen.
The results are similar to those for the low complexity analysis of Fig. 2 A. At these high complexities, the surrogates constructed with constant firing rate considerably underestimate the number of patterns expected. However, the surrogates with ISI distributions and firing rate modulation matched to the experimental spike trains produced similar numbers of patterns to those observed experimentally; in no case did the experimental number of patterns lie outside the 99% confidence limits.
Over the four datasets so analyzed, 38 combinations of complexity and number of occurrences had more than 10 patterns in the mean of 10 surrogates. In 10/38 instances, the number of patterns in the experimental data were lower than the 99% confidence limits determined from the surrogates; in no case was it higher. This was a significant departure from the null hypothesis (P < 10^{−} ^{13}, binomial test), but in the direction of fewer patterns seen experimentally than expected by chance.
Figure 3 B shows a scatter plot of the number of patterns expected from the surrogates and seen experimentally for the high complexity analysis. The line marks the identity line y= x; although the points lie close to this line, 34/38 are below it, confirming a trend for there to be fewer patterns in the experimental compared with the surrogate data.
Effect of different kernel widths in firing rate estimation
To generate the surrogate datasets, it was necessary to estimate the instantaneous firing rate (Fig. 1 B). When a kernel estimator is employed for this, a decision must be made on the width of the kernel. For the main analyses here described, the width was set equal to the mode of the interspike interval histogram for the neuron under consideration. To investigate the sensitivity of the results to changes in the kernel width, surrogates from one dataset were generated using a kernel width that was a constant k times the modal interspike interval; k was varied from 0.5 to 5.
If a narrow kernel function is used, a high proportion of its density will lie close to the time of the original spike. This would lead to the preservation of part of the precise timing information in the surrogate spike train. Figure4 A shows the fraction of the kernel area that lay within the central 3 ms (illustrated schematically by the inset in Fig. 4 A), as a function of the parameter k; this fraction has been averaged over the seven neurons recorded in the dataset used for this analysis. The curve is rapidly decreasing, showing that for values of k above approximately 0.5, little timing information is retained to within the 3ms precision that we used to search for repeating patterns.
Figure 4 B shows how the number of patterns of complexity 11 that occurred twice varied in the surrogate datasets generated using different values of k. For all of these surrogates, all but one of the patterns of complexity 11 occurred twice; the remaining one occurred three times. The pattern count shows little change up tok = 1 and is close to the value seen in the experimental data (dashed line); thereafter it declines slowly. This contrasts with the rapid decline in preservation of temporal information as k increases (Fig. 4 A). We therefore conclude that repeating patterns are not generated by neural mechanisms related to precise spike timing. The gradual decline in the number of patterns with high k is probably due to oversmoothing of the firing rate estimate, so that it is unable to follow the rapid changes in rate that occur in some cells. This reduces the number of repeating patterns present in the surrogates. The very low number of patterns seen in surrogates with constant rate (Fig.3 A) represents this process taken to its extreme (infinite kernel width).
Composition of patterns
There has been some controversy in the literature on whether patterns involving one (Abeles et al. 1993), or more than one (Prut et al. 1998; Villa et al. 1999) neuron are most common. Figure5 A presents the results of an analysis that addresses this issue. For a dataset used in the “high complexity” analysis above, each pattern of complexity 9 that had been found was examined and classified according to the number of neurons participating in it. The fraction of patterns so classified is shown by the black bars for the experimental data, and for the most realistic surrogates (gamma ISI distribution, variable firing rate) by the open bars. The error bars mark 99% confidence limits on these proportions. The hatched bars show the results using a surrogate with Poissondistributed ISIs and variable firing rate. Only 0.21% of the 1,010,065 patterns analyzed involved 1 neuron; most included the discharge of 4 or 5 of the recorded cells. The distribution of the proportion was well matched by both types of surrogate, although the most realistic surrogate (variable rate gamma) gave the closest fit to the distribution seen experimentally. Similar results were seen among the other three datasets used in the high complexity analysis: the great majority of patterns contained spikes from more than one neuron.
It is of interest to see whether spikes involved in repeating patterns have a preferred distribution of intervals; this could indicate the involvement of a regular process in the generation of the patterns (Lestienne and Strehler 1987) and would provide evidence for a special mechanism for generating repeating patterns even though they are not more numerous than expected by chance. Figure5 B shows the distribution of intervals between spikes in the patterns of complexity 9 analyzed in Fig. 5 A. This distribution has been computed by considering all possible intervals, and not just first intervals. Thus if a pattern comprised the cells (1, 2, 3) firing at times (0, 10, 25) ms, intervals of 10, 25, and 15 ms would be counted in the histogram of Fig. 5 B. Note that since intervals between the spikes of different neurons can enter this analysis, some intervals shorter than the refractory period of a single cell can occur; these simply indicate the near synchronous discharge of two cells.
The thin line of Fig. 5 B shows the interval distribution seen in patterns found in the experimental data, and the gray shaded region shows the 99% confidence limits of the distribution seen in the surrogate datasets that had a modulated firing rate and gamma distribution of ISIs. The fit between these is very good, and there is no evidence for a nonchance distribution of intervals within the patterns. The thick line in Fig. 5 B shows the interval distribution seen when the surrogate with modulated firing rate but Poisson ISIs is used. This has an excess of short intervals above those seen in the experimental data. There is a dip in the interval distribution at 0 ms for the Poisson curve. This results from the 3ms binwidth at which all of this analysis was carried out. For this surrogate, a given neuron could generate a 3ms interval with the spikes of any other neuron in the dataset, including itself. However, a 0ms interval could be formed only in combination with another cell; this reduces the probability of a 0ms interval compared with a 3ms one, producing the dip at 0 ms seen in the thick line in Fig.5 B.
Figure 5 C shows the ISI distribution histogram summed over all nine neurons in this dataset. The mode of the ISI distribution lies close to the peak in the pattern spike interval distribution, indicating that the ISI statistics of the individual cells can affect the timing of the repeating patterns. It is clearly important, therefore when examining patterns for intervals preferred above chance to compare the experimental data with surrogates that replicate the observed ISI distribution.
Task relationship of patterns
An important argument for precise repeating patterns encoding information in the cortex is the fact that they appear preferentially in association with particular behavioral events (Abeles et al. 1993). The analysis described above simply looked at the total number of patterns in the entire dataset analyzed. Although the overall number of patterns was not in excess, it remained possible that they could show a preferential distribution over the task performance. The following analysis was designed to test this idea.
For each of the “low complexity” analyses described above, we first determined the onset time of each occurrence of a repeating triplet. The rate of these events was then determined in 100ms bins aligned to the End Hold behavioral marker; this marker is shown as time 0 on the abscissa in Fig. 6 and signifies the successful completion of a 1s hold of the precision grip within the displacement target limits. The thick line in Fig.6 A shows the rate of the triplets seen in one experimental dataset. The high rate (note the kHz scale) reflects the fact, mentioned already, that on average each spike participated in a large number of patterns. It is clear that there is a profound modulation in the number of repeating patterns during the task, with a peak just before and after the hold period (the 1 s prior to the alignment marker). The gray shaded region shows the 99% confidence limits of the pattern rate assessed from the surrogate data with realistic ISI distribution and firing rate modulation. This shows the same pattern of taskdependent modulation. Figure 6 B presents the total discharge rate of the three neurons that were recorded in this dataset. The overall firing rate is greatest for these cells before and after the hold period (−2 to −1 s and 0 to 1 s, respectively); this modulation in discharge rate is clearly responsible for the changes in the pattern rate seen. Since the surrogate data have a similar firing rate modulation as the experimental data, they too show a change in the number of repeating patterns with different phases of the task.
Figure 6 C shows a similar analysis for patterns of complexity 9 in one of the datasets subjected to the high complexity analysis. Once again, a profound modulation in the number of patterns is seen with task phase; however, this is closely mirrored by the behavior of the surrogates (99% confidence limits marked in gray). Although the summed firing rate of the nine neurons in this dataset (Fig. 6 D) was greatest at the time when most patterns occurred, the depth of modulation in the pattern rate was considerably greater. This probably reflects the highly nonlinear dependence of high complexity patterns on firing rate (Abeles and Gerstein 1988).
The datasets illustrated in Fig. 6 showed considerable modulation of the neuronal firing rate with task performance. The extent of such modulation varied over the datasets analyzed; similar results were, however, seen to those in Fig. 6 for all other datasets, independent of how deep the modulation of neuronal firing rate was.
Patterns with a fraction of the component spikes absent
The synfire chain theory (Abeles 1991) is expressed in probabilistic terms, and it suggests that patterns may not repeat exactly on every repetition. In this case, it is important to investigate the occurrence of patterns where not all of the constituent spikes are required to be present for a repetition to be counted. Figure 7 shows the results of such an analysis. In one of the datasets subjected to high complexity analysis, a single pattern of complexity 20 that occurred in its exact form twice was chosen. The two exact repetitions are shown in Fig. 7 Aas rasters of the four neurons that contributed to this pattern; the spikes that formed the pattern are shown with thick lines.
Figure 7 B shows that this pattern had a highly specific relationship to the task. The displacement signals from the finger and thumb levers of the precision grip manipulandum are shown aligned to the onset of the two pattern occurrences; the gray shaded region shows the pattern duration. On each repetition, shortly after the last spike in the pattern, the animal began to squeeze the levers to move them into the target zone.
A search was carried out to find instances of this pattern in which up to half (10/20) of the component spikes were allowed to be missing. Sixtytwo such inexact matches were found and are illustrated in raster form in Fig. 7 C. Remarkably, these showed a similar highly specific relationship to the task as seen for the exact form. This is shown in Fig. 7 D as overlain lever displacement signals; in almost all cases, the inexact patterns also occurred just before the levers were squeezed to initiate a trial.
However, analysis of the surrogate data indicated that this result was likely to be simply due to the taskdependent modulation in firing rate of the neurons recorded. Over the 10 surrogates, between 57 and 71 (mean 62.8) inexact matches to the template of Fig. 7 A were found. Figure 7 E shows the manipulandum lever position signals aligned to the onset of the matches found in one of these surrogate datasets. A very similar pattern of task dependency is present as in Fig. 7 D for the matches made with the original experimental data.
Finally, Fig. 7 F plots the distribution of the occurrence times of the inexact pattern matches relative to the End Hold behavioral marker. The thick line shows the distribution in the experimental data, and the gray shaded region the range in the 10 surrogates. There is good agreement between these. The inexact pattern repetitions and their task dependency would also therefore seem to be explicable simply by the modulation in firing rate of the neurons whose spikes are analyzed.
DISCUSSION
Significance of repeating patterns
The results presented in this study extend the recent findings ofOram et al. (1999) in visual cortex and LGN in a number of directions. The number of repeating patterns found in our data from monkey motor cortex was no more than that expected in spike trains which simply mirrored the modulation of firing rate and ISI statistics of the experimental data (Figs. 2 and 3). This was the case even though we searched for patterns in longer sections of data, with greater allowed maximum pattern length, and at a wider range of complexities than Oram et al. (1999). Additionally, while Oram et al. (1999) analyzed only the discharge of a single neuron, we have shown that patterns involving more than one cell are also not evident in numbers greater than expected by chance. As well as a good match between the overall numbers of patterns in experimental and surrogate data, we have also shown similar task dependence in the preferred times of pattern occurrence, both for exact (Fig. 6) and inexact (Fig. 7) pattern repetitions. Finally, detailed descriptive statistics of the patterns themselves, such as the number of neurons participating in a pattern (Fig. 5 A) and the intervals between spikes in a pattern (Fig. 5 B) were very similar between experimental data and the surrogates. We therefore conclude that repeating patterns in data from the motor cortex are nothing remarkable, but arise solely as an epiphenomenon of the modulation in neural firing rate.
Two features of the surrogate datasets that we used allowed them to match the number of repeating patterns found in the experimental data. The first feature, modulation in firing frequency, was the most important. Surrogates with a constant rate did not in general provide an accurate estimate of the number of patterns seen in the experimental data. The number appeared to be overestimated at low complexity (triplets, some quadruplet combinations), but underestimated at higher complexities (Figs. 2 and 3). This difference is probably related to the precise details of the rate modulation of the units analyzed, and to the way in which pattern counts vary with rate. Triplet counts will depend on a product of three firing rates, while, for example, complexity 10 counts will depend on a product of 10 rates (Abeles and Gerstein 1988). If mean rates are substituted for the actual, timevarying rate profiles, this could therefore have very different effects on the estimated number of high versus low complexity patterns. Although several previous reports on pattern significance have not attempted to quantify the effects of firing rate modulation, the data presented have often shown extensive changes in rate (e.g., Fig. 1 of Villa and Abeles 1990; Figs. 3 and 4 of Abeles et al. 1993).
The surrogates used here were generated with a momentbymoment estimate of the instantaneous firing rate, rather than relying on a trial averaged firing rate profile. In work involving awake behaving animals, there is often considerable variation in neuron response from one trial to the next, which can make the latter unreliable (Brody 1999a,b; Pauluis and Baker 2000; Baker et al., unpublished observations). This was noted by Oram et al. (1999), who found that it was necessary to generate surrogates with the number of spikes in a given trial fixed to be equal to that seen experimentally to match the observed pattern count. However, not only the magnitude, but also the latency or shape of the neuronal response can vary from trial to trial (see Brody 1999a,b; S. N. Baker and G. L. Gerstein, unpublished observations); this may have been less important in the Oram et al. (1999) study, due to the short duration of their trials and the lower rate modulation compared with the present work. In certain situations therefore, it can be inadequate to use simply the trialaveraged firing rate profile, scaled according to the number of spikes observed on that trial. An instantaneous firing rate measure as used here is capable of capturing such subtle changes in cell responses, which could have an important effect on the number of patterns generated.
The second feature mimicked by the surrogate spike train was the interspike interval statistics of the experimental spike trains (Fig.1 F). This permitted a better fit to the number of repeating patterns seen at a given complexity and occurrence between experimental and surrogate data. In addition, inclusion of a realistic ISI distribution, incorporating a relative refractory period, allowed the patterns in the surrogates to match the preferred interval between spikes in patterns found in experimental data (Fig. 5 B).Oram et al. (1999) have also noted an improved match between surrogate and experimental pattern counts when a refractory period was imposed on the surrogate spikes.
While there were not large differences between repeating patterns in surrogate and experimental data such as would be expected if patterns formed a temporal code in the cortex, two small discrepancies are worthy of note. First, for the low complexity analysis, in 16/480 combinations of complexity and occurrence the experimental pattern count lay outside the 99% confidence limits calculated from the surrogates; this was significantly different from the 1% (4.8 combinations) expected by chance. The finding that almost equal numbers of these departures were above as below the confidence limits argues against them reflecting a weak excess of patterns in the data. Rather, it seems likely that the distribution of number of patterns does not entirely follow the Gaussian density that was assumed for the calculation of the confidence limits. A distribution that is heavy tailed compared with the normal curve would produce the observed result. As noted in methods, the pattern count at a given level of complexity and number of occurrences is the sum of a large number of Poisson distributed counts of individually identified patterns. The Poisson distribution is heavy tailed compared with its Gaussian approximation (Press et al. 1989, p. 550), so that this explanation for the discrepancy does not seem unreasonable. Full characterization of this distribution would have required analysis of a considerably higher number of surrogates per dataset (e.g., 10,000 rather than the 10 used here); this was not carried out due to the greatly increased computation time that it would entail.
A different explanation must be sought for the finding in the high complexity analysis that 18% of the combinations of complexity and occurrence analyzed showed significantly fewer patterns in the experimental data compared with the surrogates; this trend was confirmed by inspection of the remaining combinations, which almost all had slightly lower experimental pattern counts than the average of the surrogates (Fig. 3 B). The most likely explanation for this discrepancy lies in the rather simple nature of the surrogates themselves. First, the firing rate modulation was assumed to occur in a gradual manner, such that it could be represented realistically by the kernel firing rate estimator. However, neurons may on occasion increase their firing rate abruptly (see Pauluis and Baker 2000). Smoothing across such rate change points would lead to an overestimate of the firing rate just before the change, and an underestimate after it. It is conceivable that these errors could lead to the slight excess of repeating patterns seen in the surrogates. Second, it is possible that the underlying distribution of the ISIs did not exactly follow the gamma distribution assumed for the surrogates. Finally, surrogates were simulated with the same order of the gamma ISI distribution throughout the dataset. It is possible that the experimentally recorded neurons could have changed the regularity of their firing at different times in the task performance, and that this would lead to the minor discrepancies observed. It is interesting to note that only in the high complexity analysis is the discrepancy between experimental and surrogate counts seen. The rate of such repeating patterns is known to depend on high powers of the neurons' firing rates (Abeles and Gerstein 1988); their expectancy is likely to be exquisitely sensitive to minor inaccuracies in the surrogates.
We made no attempt to produce more realistic surrogates that matched the experimental data better, since the current analysis had already answered the issue we set out to address. If such simple surrogates as used here can produce as many if not more repeating patterns as in the experimental data, the hypothesis that repeating patterns are generated by cortical circuits specialized for that purpose and used to encode information is unlikely.
Consequences for theories of temporal coding
We conclude that there is no evidence in our experimental data to support the idea that neurons specifically generate precise (millisecond range) temporal codes. It is possible, as suggested byLestienne and Strehler (1987) and Lestienne (1999), that neurons are nevertheless sensitive to particular spatiotemporal patterns of inputs. Since the rate of patterns depends on a power of the neuron firing rate, this could lead to a form of amplification, by which cells effectively respond to a power of the firing rate of their inputs instead of linearly. The utility of raising a noisy signal to a power to improve the signaltonoise ratio has been explored by Baker and Gerstein (2000). Larkum et al. (1999) have recently provided evidence for a biophysical mechanism that would make single cortical cells sensitive to spatiotemporal patterns in their inputs. As yet, it is not clear to what extent cell firing in the conscious animal is governed by such mechanisms.
Two possibilities for cortical information coding remain following the present study, and the related work of Oram et al. (1999). It is possible that neurons simply transmit information in their firing rate, which must be averaged over a suitable time interval and population of similarly responding cells to obtain a low noise signal (Lee et al. 1998; Shadlen and Newsome 1998). An additional possibility, which does not necessarily exclude the first, is that information is contained in the broader time course synchronous firing of neurons that is commonly observed as ∼20ms peaks in crosscorrelation histograms of cortical cells (Fetz et al. 1990) (see Fig. 1 H), and may be accompanied by oscillatory activity (Murthy and Fetz 1996; Singer and Gray 1995). Such synchrony is different both from the millisecond level timing investigated here and proposed by the synfire chain theory (Abeles 1991), and from a broader timecourse comodulation in firing rate (Miller et al. 1993; Roelfsema et al. 1998). Synchrony with medium jitter has been shown to modulate with behavioral state (Baker et al. 1997; Singer and Gray 1995;Vaadia et al. 1995; Baker et al., unpublished observations), making it a plausible substrate for information processing in the cortex.
Acknowledgments
The authors thank N. Philbin, R. Spinks, and E. Pinches for assistance with the recordings, and G. Gerstein, P. Kirkwood, and D. Wolpert for critical comments on this work. The use of the supercomputer at the Cambridge High Performance Computing Facility is gratefully acknowledged.
This work was funded by the Medical Research Council, the Wellcome Trust, and the Brain Research Trust.
Footnotes

Present address and address for reprint requests: S. N. Baker, Dept. of Anatomy, University of Cambridge, Downing Street, Cambridge CB2 3DY, UK (Email: snb11{at}cam.ac.uk).

The costs of publication of this article were defrayed in part by the payment of page charges. The article must therefore be hereby marked "advertisement" in accordance with 18 U.S.C. Section 1734 solely to indicate this fact.
 Copyright © 2000 The American Physiological Society