Abstract
Visual stimuli can evoke waves of neural activity that propagate across the surface of visual cortical areas. The relevance of these waves for visual processing is unknown. Here, we measured the phase and amplitude of local field potentials (LFPs) in electrode array recordings from the motion-processing medial temporal (MT) area of anesthetized male marmosets. Animals viewed grating or dot-field stimuli drifting in different directions. We found that, on individual trials, the direction of LFP wave propagation is sensitive to the direction of stimulus motion. Propagating LFP patterns are also detectable in trial-averaged activity, but the trial-averaged patterns exhibit different dynamics and behaviors from those in single trials and are similar across motion directions. We show that this difference arises because stimulus-sensitive propagating patterns are present in the phase of single-trial oscillations, whereas the trial-averaged signal is dominated by additive amplitude effects. Our results demonstrate that propagating LFP patterns can represent sensory inputs at timescales relevant to visually guided behaviors and raise the possibility that propagating activity patterns serve neural information processing in area MT and other cortical areas.
SIGNIFICANCE STATEMENT Propagating wave patterns are widely observed in the cortex, but their functional relevance remains unknown. We show here that visual stimuli generate propagating wave patterns in local field potentials (LFPs) in a movement-sensitive area of the primate cortex and that the propagation direction of these patterns is sensitive to stimulus motion direction. We also show that averaging LFP signals across multiple stimulus presentations (trial averaging) yields propagating patterns that capture different dynamic properties of the LFP response and show negligible direction sensitivity. Our results demonstrate that sensory stimuli can modulate propagating wave patterns reliably in the cortex. The relevant dynamics are normally masked by trial averaging, which is a conventional step in LFP signal processing.
- cerebral cortex
- cortical oscillations
- cortical waves
- local field potentials
- spatiotemporal dynamics
- visual processing
Introduction
Neural activity in cerebral cortex can self-organize into spatiotemporal wave patterns that propagate across the cortical surface and link the activity of distinct neuron populations. These patterns have been observed in recordings from multielectrode arrays (Freeman and Barrie, 2000; Rubino et al., 2006; Townsend et al., 2015; Zanos et al., 2015) and from voltage-sensitive imaging (Prechtl et al., 1997; Wu et al., 2008; Huang et al., 2010; Muller et al., 2014). Studies of spontaneous (intrinsic) activity patterns showed that their temporal dynamics are not random, but rather follow repeated motifs (Mohajerani et al., 2013; Townsend et al., 2015). Sensory stimulation also produces reproducible wave patterns (Sato et al., 2012; Muller et al., 2014), but the central questions of whether and how these waves can represent incoming sensory signals remains unsolved.
Here, we used multielectrode array recordings to measure propagating activity patterns in local field potentials (LFPs) in the motion-processing medial temporal (MT) area of anesthetized marmoset monkeys. Visual stimuli comprised drifting grating and dot-field stimuli. We found that, on a trial-by-trial basis, the propagation direction of LFP waves is modulated by stimulus motion direction. After averaging recordings over multiple trials, the LFP signals instead display a characteristic “source-to-sink” pattern of amplitude, which is less sensitive to stimulus motion direction. We show that this difference arises because the averaged LFP pattern is primarily generated by additive signals (Sauseng et al., 2007; Becker et al., 2008) and is insensitive to the phase of ongoing cortical oscillations (Gruber et al., 2014). Our fine-scale spatiotemporal analysis thus reveals cortical dynamics that may support visual representations that are masked by the standard process of cross-trial signal averaging.
Materials and Methods
Experimental design.
Recordings were made from the MT area of four adult male marmosets (Callithrix jacchus). Details of preparation have been described previously (McDonald et al., 2014; Solomon et al., 2015). Anesthesia and analgesia were maintained by intravenous infusion of sufentanil citrate (6–30 μg kg−1 h−1) and an inspired 70:30 mixture of N2O and carbogen (5% CO2, 95% O2). Dominance of low frequencies (1–5 Hz) in the electroencephalogram (EEG) recording and absence of EEG or electrocardiogram changes under noxious stimulus (tail pinch) were taken as the chief signs of an adequate level of anesthesia. We found that low anesthetic dose rates in the range cited above were always very effective during the first 24 h of recordings. Thereafter, drifts toward higher frequencies (5–10 Hz) in the EEG record were counteracted by increasing the rate of venous infusion or the concentration of anesthetic. The typical duration of a recording session was 48–72 h.
Stimuli were presented on a cathode-ray-tube monitor (Sony G500, refreshed at 100 Hz, viewing distance 45 cm, mean luminance 45–55 cd m−2). Stimuli comprised drifting sine-wave gratings (Michelson contrast 0.5, spatial frequency 0.2–0.5 cycles/degree, temporal frequency 4–5 Hz) or fields of drifting circular white dots (Weber contrast 1.0; dot diameter 0.4°; drift velocity 20 deg/s). For both gratings and dot fields, different motion directions (90° steps) were presented in a large, stationary circular window (30°) for 2 s and then the screen was held at mean luminance for 2 s. The procedure was repeated until 100 repetitions were made for each of the four directions. In three recordings, the contralateral eye was occluded; in the others, the ipsilateral eye was occluded. Data were recorded using multielectrode arrays (10 × 10 electrodes, 1.5 mm length, electrode spacing 400 μm; Blackrock Microsystems). Recording surface insertion depth was targeted to 1 mm.
We denote the full LFP recording from one animal for each stimulus type (dot field or grating) as x(d, r, s, t) for stimulus direction d (of D = 4 perpendicular motion directions), trial repetition r (of R = 100 repetitions), recording site s (of S = 96 electrode sites), and sampling time t (of a 4 s trial, sampled at 1024 Hz). For some analyses, we considered the trial-averaged LFPs, which we denote by xav(d, s, t). For the analyses reported here, responses to gratings and dot fields did not differ significantly, so responses to both stimuli were pooled for further analysis. Analyses of simultaneously recorded single- and multiunit extracellular activity have been described previously (Solomon et al., 2015, 2017).
Oscillation analysis.
Complex wavelet coefficients X(f, d, r, s, t) were calculated from x(d, r, s, t) in each recording for oscillations at center frequencies f from 2 to 20 Hz at 1 Hz intervals using 7 cycle Morlet wavelets (Torrence and Compo, 1998). Amplitude A = | X | and phase ϕ = arg(X/| X |) were calculated at each frequency from these components. Trial averaged LFPs xav were similarly processed into averaged wavelet coefficients Xav(f, d, s, t). Wavelet coefficients X and Xav were sorted into delta (1–4 Hz), theta (5–8 Hz), alpha (9–13 Hz), or beta (14–20 Hz) ranges. All further statistics were first calculated in original, narrow-band wavelet coefficients and then averaged within oscillation bands to give overall results. Line graphs were smoothed with a 10 ms moving average window. Synchrony between recording sites was quantified by the zero-lag correlation, rsyn(f, d, r, sA, sB, t) calculated from filtered time-series Re(X) between every pair of recording sites sA and sB with a 250 ms sliding window.
Phase and amplitude velocity fields.
Phase maps in single-trial activity were processed to obtain spatiotemporal phase velocity fields (PVFs), following the procedure described by Townsend et al. (2015). The PVFs characterize spatial and temporal changes of phase maps and are constituted of a set of vectors indicating direction and speed of phase change between consecutive time steps at each recording site. We used a fixed time step given by the sampling rate of ∼1 ms between phase maps. The precise choice of time step was unimportant; we found that downsampling LFPs by a factor of 10 resulted in negligible changes to computed PVFs at all frequencies up to 20 Hz (data not shown). We characterized the motion present in each individual PVF v(x, y) comprising N vectors at locations (x, y) by calculating a number of descriptive statistics: the average velocity magnitude, Mean speed = ∑x,y‖v‖/N, overall motion direction, Mean direction = arg(∑x,y v), and wave coherence (also called the average normalized velocity, or order parameter), Coherence = ‖∑x,yv‖/∑x,y‖v‖. The wave coherence ranges from 0 to 1 and represents the degree to which all velocity vectors are aligned, with a value of 1 indicating all vectors aligned to the same direction. We defined plane wave activity to be present when wave coherence was above a threshold of 0.85 for a continuous period of at least 10 ms and other propagating patterns to be present when coherence was 0.5–0.85 for at least 10 ms. These thresholds were chosen because they correspond to clearly visible differences on manual classification of patterns. We found that changing the threshold values by ±0.05 did not influence the number of detected patterns by >10%. To identify stimulus-generated propagations, we identified plane waves and propagating patterns in all trials from 0 to 200 ms after stimulus onset.
To detect spatiotemporal patterns in trial-averaged LFPs, amplitude velocity fields (AVFs) from trial-averaged amplitude maps were calculated in a parallel fashion to PVFs in single trials. Complex wave patterns (sources and sinks) were classified in trial-averaged AVFs by the Jacobian, winding number, and curl around critical points in velocity fields, as described previously (Townsend et al., 2015). To evaluate whether observed complex wave dynamics in AVFs could be an artifact caused by a spatially static but noisy system, we constructed surrogate LFPs, xsur, with separable spatial and temporal scales in the following form:
where As is the spatial envelope, At is the temporal envelope, ft is the target oscillation frequency, and N(0, σ) is normally distributed white noise with mean 0 and standard deviation σ. We chose As to be a Gaussian with random center location and size (full width at half maximum 1–2 mm), At to follow a linear increase and then decrease to generate a source then sink pattern, and σ such that the signal-to-noise ratio of xsyr was equal to 0.1. At higher noise levels, source/sink structure in the data could rarely be recovered. Surrogates xsur were then processed using the same procedure as for xav to find source and sink pattern locations.
Discrimination between stimuli.
In single-trial PVFs, the ability for the mean speed, mean direction, and coherence statistics described above to discriminate stimulus motion direction was quantified using dissimilarity measures. Dissimilarity measures compare the variability across trials of different stimuli with the variability across trials within the same stimulus (Luo and Poeppel, 2007). We also considered the dissimilarity of zero-lag pairwise correlations rsyn in single-trial filtered LFPs. Variabilities for rsyn, mean speed, and coherence were characterized by the coefficient of variance, CV = σ/|μ|, where μ and σ are the mean and SD across trials of the relevant statistic and variability for mean direction Θ was characterized by the angular variance (Philipp, 2009), SΘ = 1 − e̅i̅Θ̅, where the bar denotes a mean taken across trials. Variability was computed in within-stimulus groups comprising 100 trials for each stimulus direction and across-stimulus groups comprising a random selection of 25 trials from each of the four stimulus directions. For each within-stimulus measure, 10 across-stimulus measures were calculated using different random permutations. Dissimilarity was then calculated as the mean across-stimulus measure minus the mean within-stimulus measure. For trial-averaged activity, dissimilarities cannot be computed because the within-stimulus variability is not defined. Instead, differences between stimulus motion directions in trial-averaged amplitudes were quantified at each electrode using the coefficient of variance across all four stimulus directions and then averaged across all recording sites.
In addition to the dissimilarity measures introduced above, we used linear support vector machines (SVMs) (Burges, 1998; Chen et al., 2015) to quantify the capacity of PVF wave statistics to discriminate stimulus direction. For each pair of stimulus directions in a recording, we trained SVMs to classify stimuli by finding an optimal separating hyperplane through all trials of mean speed, mean direction, or coherence taken from 100 to 300 ms after stimulus onset. We selected stimulus pairs from adjacent (90°) or opposite (180°) stimulus directions and calculated the fraction of trials that could be correctly classified by this hyperplane. As a control, we also tested SVMs in shuffled data by repeating the calculation 10 times with all trials randomly assigned to different stimulus directions. We additionally performed a 10-fold cross-validation procedure (Kohavi, 1995) whereby we partitioned all trials from pairs of stimuli into 10 subsamples and then used each subsample as test data for SVMs trained on the other nine subsamples. We repeated this cross-validation five times with different randomly selected partitions in each of the three wave statistics for all pairs of stimulus directions in theta PVFs in all recordings. For PVF mean directions, we used the sine and cosine of direction to remove the circularity of the data. All SVMs were computed using MATLAB's Statistics and Machine Learning Toolbox 2016b (The MathWorks, RRID:SCR_001622).
Mechanisms of trial-averaged activity.
Phase resetting of ongoing oscillations can be characterized using intertrial coherence (also known as the phase-locking factor or phase-locking index), ITC(f, d, s, t) = |e̅i̅θ̅|. Intertrial coherence measures the alignment of phase across trials at each recording site. However, this measure cannot differentiate true phase resetting from additive evoked activity that is time locked to stimulus presentation. This limitation arises because additive activity can also modify the phase of oscillations (Sauseng et al., 2007; Martínez-Montes et al., 2008). A more selective index of phase resetting is obtained by calculating the uniformity of phase across trials after subtracting the sample mean to remove evoked effects (Martínez-Montes et al., 2008). Specifically, uniformity is quantified by the second trigonometric moment of wavelet coefficients after removing the sample mean and confining phase to the range [0, π], R2(f, d, s, t) = |e̅2̅i̅(̅θ̅*̅)̅|, where θ* = arg(X − X̄) is the phase of wavelet coefficients after subtracting the sample mean across trials. Additive activity is measured by the trial averaged wavelet coefficient amplitude, | X̄ |. This measure does not distinguish between activity that is temporally locked or not temporally locked to the stimulus, but it is unaffected by the phase of trials (i.e., unaffected by induced activity). Both measures were compared with trial-averaged amplitude by calculating Pearson's correlation coefficient across space at each time step between R2 and | Xav | and between | X̄ | and | Xav | and then averaging across recordings.
Results
We analyzed measurements of the LFP obtained by implanting planar arrays of electrodes over area MT of marmoset monkeys. A schematic view of the stimulus and recording arrangement is shown in Figure 1A. Responses to multiple repetitions of one stimulus (a moving dot field) are shown for an example electrode in Figure 1B. Each stimulus presentation generated LFP power increases in the range 2–20 Hz, with most of the power in the LFP in delta and low-theta bands (Fig. 1C,D). To characterize spatiotemporal activity patterns, we decomposed the LFP using Morlet wavelets into 19 frequency bands between 2 and 20 Hz and extracted the phase in each band (Fig. 2A,B; see Materials and Methods). We then created phase maps at 1 ms intervals (Fig. 2C) and generated PVFs for each consecutive pair of maps (Fig. 2D). These PVFs characterize the local magnitude and direction of phase propagation and reveal propagating activity patterns, including complex waves, as was observed previously in spontaneous delta-band activity (Townsend et al., 2015). Propagating activity was present in all frequency bands and was qualitatively similar in prestimulus and poststimulus conditions. Initial analysis of the (weaker) oscillations at higher frequencies (25–50 Hz) showed that these frequency bands displayed similar dynamics to beta activity (data not shown). For expediency, we did not consider these higher frequencies further. We characterized the alignment of activity across the recording array by the wave coherence of PVFs (see Materials and Methods). We found that 84% of trials across all recordings and frequencies (24290/28800) had coherence >0.5 for at least 10 ms in the 200 ms period after stimulus onset, indicating the presence of propagating patterns caused by visual stimuli. Only 26% (6224) of these trials exhibited planar traveling waves, showing that the stimulus-driven response instead usually took the form of complex waves, as we showed previously for activity in the absence of patterned visual stimuli (Townsend et al., 2015). Propagating patterns were typically active across 25–100% of the ∼16 mm2 recording area and extended across multiple direction columns in MT (repeating unit <0.3 mm; Solomon et al., 2017) and a large region of visual space (receptive fields in the recording area spanned ∼5–40° eccentricity; Rosa and Elston, 1998; Chen et al., 2015).
Experimental recordings of LFPs. A, Representation of marmoset eye and brain showing approximate position of MT area and electrode array. B, Raw LFPs at one electrode from contralateral eye for 100 repetitions of the same dot-field stimulus. Stimulus onset is at t = 0 and lasted 2 s. Yellow trace indicates single trial used in Figure 2. C, Power spectrogram of LFP calculated separately for each trial and then averaged across all trials, channels, and recordings. Power was concentrated in delta and low-theta bands. Noncausal signal at low frequencies reflects the filtering. D, Power spectrum as a proportion of prestimulus power in the relevant frequency band.
Identification of propagating activity patterns in LFPs. A, Single-trial waveforms (transparent lines) from two channels filtered with Morlet wavelets to extract 10 Hz oscillations (solid lines), with oscillation amplitude shown by the dotted lines. B, Phase of the 10 Hz oscillations shown in A. C, 2D phase maps at each of the time steps indicated in B. Channels with waveforms shown in A and B are indicated. D, Phase velocity fields calculated between consecutive phase maps (typically 1 ms apart; larger time distance is shown here for display purposes). Properties of the PVF are visualized: mean speed is the average vector length, mean direction is the angle of the mean vector, and coherence is close to one because most vectors are aligned in the same direction.
Stimulus dependence of propagating patterns
Individual neurons in area MT show strong selectivity for the direction of movement within their receptive fields (Dubner and Zeki, 1971; Maunsell and Van Essen, 1983; Rosa and Elston, 1998). We therefore investigated whether stimulus motion direction influenced features of the propagating patterns. We compared the variability across trials of one stimulus motion direction, with the variability across the same number of trials drawn at random from the entire set of motion directions. The difference between these two values is their dissimilarity, with greater dissimilarity indicating greater dependence on stimulus motion direction (see Materials and Methods). We first looked for stimulus-dependent patterns of synchronization in activity across the array by calculating zero-lag correlations between all pairs of channels. Pairwise correlations were increased after stimulus onset (Fig. 3A) in all frequency bands, but the dissimilarity index revealed no dependence of correlations on stimulus direction, as shown for theta frequencies in Figure 3B. Delta, alpha, and beta frequencies likewise displayed near-zero dissimilarities. These results are consistent with a previous study showing that LFP coherence in area MT is independent of stimulus motion direction at the oscillation frequencies studied here (Solomon et al., 2017).
Sensitivity of LFP synchrony structure and wave statistics to stimulus motion direction. A, Average zero-lag LFP correlation across all channel pairs averaged across all recordings and frequency bands. Noncausal effects are due to time smoothing of signal filtering. Error bars indicate SEM. B, Average dissimilarity in theta oscillations averaged across all recordings. Synchrony dissimilarity was calculated using the CV of pairwise correlations and averaged across all channel pairs. PVF mean speed and coherence were also calculated using the CV. PVF mean direction was calculated using the circular SD. Shading shows SEM across all trials in all recordings. C, Average dissimilarity in PVF wave mean direction in theta (5–8 Hz) and delta (2–4 Hz) bands. D, Same as C, but for alpha (9–13 Hz) and beta (14–20 Hz) bands.
We next investigated whether stimulus features influenced the properties of propagating activity patterns. We considered three ways in which the waves could be modified. First, stimulus motion might change the type of patterns present (e.g., coherent plane waves vs other complex waves). Such changes would create direction-sensitive differences in the overall coherence of the wave patterns. Second, stimulus direction might change the propagation speed of complex waves in PVFs. Third, stimulus direction might change the propagation direction of complex waves. Our analysis ruled out the first two possibilities; that is, we found no stimulus-dependent change in coherence or average speed of propagating patterns in any frequency band (Fig. 3B). In contrast, however, the average direction of wave propagation showed sensitivity to stimulus direction in all recordings (Fig. 3B). The strongest sensitivity to stimulus direction was present in theta oscillations, with delta (Fig. 3C) and alpha (Fig. 3D) oscillations also displaying direction sensitivity across all recordings.
The direction dependence of propagating patterns is further illustrated in Figure 4. As expected, single-trial phase maps and PVFs taken before stimulus onset (Fig. 4A,C) showed no preferred direction, as signified by the small mean velocity vectors (Fig. 4A,C, white arrows) and the high trial-to-trial variability in mean direction (Fig. 4B,D). After stimulus onset, PVFs showed coherent motion propagating in a direction that depended on stimulus direction (Fig. 4E–H). Figure 4 depicts only trials with high wave coherence, but propagating patterns with lower coherence showed similar dependence on propagation direction. Typically, the stimulus dependence of propagation direction (Fig. 4I) reached a maximum at 200–300 ms after stimulus onset and was present for <500 ms overall, indicating that motion-sensitive propagating patterns are present only transiently.
Influence of drifting dot-field stimulus direction on wave propagation direction. A, Snapshots of phase maps and phase velocity fields before stimulus onset for stimulus drifting in the 90° direction. White arrows show mean velocity vector across all recording sites at double scale. B, Histogram of PVF mean direction before stimulus onset across all trials in this recording (n = 100) for 90° stimulus direction. C, Same as A, but for 180° stimulus direction. D, Same as B, but for 180° stimulus direction. E–H, Same as A–D, but after stimulus onset. I, Circular average of PVF mean direction across all trials in this recording. Shading shows circular SD (Philipp, 2009).
The dissimilarity measure is highly derived, raising the question of whether the direction sensitivity of complex waves could support visual motion discrimination. We addressed this question by training linear SVMs to discriminate stimulus direction using only the mean coherence, mean speed, or mean propagation direction of theta-band PVFs. Optimal SVMs trained on mean direction successfully decoded motion direction in all stimulus pairs across all recordings (60.0% mean correct classification rate vs 50.0% baseline), but those trained on coherence or mean speed did not (Fig. 5). We observed no difference in discriminability between opposite and orthogonal stimulus directions, but both pair types performed better than data with shuffled stimulus labels. We tested the stability of these results using a 10-fold cross-validation procedure (see Materials and Methods). Cross-validated SVMs trained on mean direction performed better than chance in >85% of pairs in each recording, whereas SVMs trained on coherence or mean speed performed at or below chance level in 60–70% of stimulus pairs. These results confirm that propagation direction of complex waves can support visual motion discrimination.
Unsupervised direction discrimination. Bar-and-whisker plots show SVM classification of stimulus for opposite-direction stimulus pairs (n = 56, 2 stimulus pairs per oscillation frequency for 7 recordings), orthogonal-direction stimulus pairs (n = 112, 4 pairs per frequency), and pairs from shuffled stimulus labels (n = 168) using the following theta-band (5–8 Hz) PVF wave statistics: wave coherence, mean propagation speed, and mean propagation direction. Boxes give 25th, 50th, and 75th percentiles of data; whisker tips give maxima and minima; black horizontal lines show median values. Line at 0.5 indicates chance level.
Trial-averaged propagating patterns
The foregoing analyses show stimulus dependence of propagating patterns in area MT on individual trials. We therefore investigated whether stimulus motion direction can likewise be extracted from trial-averaged activity. Because the amplitude of trial-averaged data collapses single-trial amplitudes and single-trial phase, we studied trial-averaged amplitude maps (we explain the relative contributions of amplitude and phase in the following section). Trial-averaged amplitude maps at all frequencies were organized spatially around a central maximum that appeared in a consistent location across stimulus directions (Fig. 6A), likely reflecting the retinotopic organization of the cortical surface in MT and the position of the visual stimulus (Allman and Kaas, 1971; Dubner and Zeki, 1971; Rosa and Elston, 1998). After stimulus onset, amplitude increased and propagated outward over a few millimeters of cortex for a duration of 100–300 ms, forming a source pattern, and then decreased and propagated inward for another 100–300 ms, forming a sink pattern. Shuffling channels before analysis obliterated the source and sink structures, showing that they are not an artifact of signal processing. Similar stimulus-evoked propagating amplitude patterns are present in voltage-sensitive dye imaging in rat cortex with comparable spatial and temporal scales (cf. Fig. 2 in Mohajerani et al., 2013).
Spatiotemporal activity patterns in trial-averaged signals. A, Trial-averaged amplitude maps and AVFs (white arrows) for different dot-field stimuli in the same recording. Source/sink centers are marked with white dots. B, Percentage of trial-averaged amplitude velocity fields showing widespread source and sink patterns when taken across recordings (n = 7) and center frequencies (n = 19). Prestimulus effects are due to signal filtering. C, Histogram of distance between trial-averaged AVF source and sink pattern centers observed across all recordings and center frequencies compared with model data of a spatially static, noisy structure. D, Relative center location of AVF sink pattern with source pattern center fixed at origin for all frequencies in all recordings. Solid orange lines show the trajectory of source and sink centers for selected points; dashed orange lines show the jump in position from last source pattern to first sink pattern. Symbols and colors correspond to recordings. E, Average channel-by-channel coefficient of variation for trial-averaged amplitude between all stimulus directions by frequency band across all recordings.
To characterize the amplitude patterns, we computed AVFs (see Materials and Methods) in each frequency band (Fig. 6A, white arrows). These AVFs show a source pattern at stimulus onset as LFP power increases, with velocity vectors directed away from a central point (Fig. 6A, 0.02 s), and then a sink pattern as power waned a few hundred milliseconds later, with velocity vectors directed toward a central point (Fig. 6A, 0.22 s). We examined AVFs at all frequencies in all recordings to identify when stable, spatially extended source and sink patterns were active (see Materials and Methods) and found that both patterns were present at most frequencies in all recordings, with source activity peaking 40 ms after stimulus onset and sink activity peaking after 180 ms (Fig. 6B). These patterns were typically separated by a period of 100–200 ms, during which there was a stable amplitude map with no propagation visible in AVFs (Fig. 6A, 0.12 s).
We found that the sink patterns could form at up to 2 mm away from their preceding source pattern (Fig. 6C, mean displacement 0.64 mm; electrode separation 0.4 mm) and that the displacement distances and directions between source and sink patterns varied within and across recordings (Fig. 6D). To test whether these measured displacements reflected separate cortical locations or noisy measurements of source and sink patterns with zero displacement, we constructed surrogate LFPs with activity at a fixed center location (see Materials and Methods). We then applied strong Gaussian white noise (signal-to-noise ratio 0.1) and observed source–sink displacements consistently lower than those observed in real data (Fig. 6C, mean displacement 0.28 mm). These results show that sources and sinks were not simply a static structure that expanded and contracted from a single point in cortex. Instead, they could be observed with different dynamics and locations and may therefore reflect separate processes of sensory response.
We tested for stimulus motion direction sensitivity in averaged LFP signals by measuring at each electrode the coefficient of variation (CV) between trial-averaged amplitudes obtained for different stimulus directions. Increased CV after stimulus onset would indicate that trial-averaged activity is primarily dependent on motion direction. We found increases in CV over prestimulus levels in fewer than 15% of electrodes in any of the recordings, and average CV decreased after stimulus onset in all frequency bands (Fig. 6E). These results indicate that trial-averaged amplitudes are not determined by stimulus motion direction. Instead, they primarily reflect the retinotopic organization of MT and the position of stimuli in the visual field.
Evoked and induced contributions to trial-averaged activity
Our observations present a puzzle: in single trials, we see stimulus-dependent propagating patterns with complex wave dynamics and high variability from trial to trial. However, in trial-averaged activity, we see different propagating patterns that are consistent across all stimuli and recordings. This observation recalls the longstanding debate over whether trial-averaged neural responses reflect stimulus-evoked events or stimulus-induced changes in ongoing activity (Makeig et al., 2002; Penny et al., 2002; Fell et al., 2004; Sauseng et al., 2007; Alexander et al., 2013). Trial-averaged activity may be explained either by a stimulus-evoked waveform, additive with background activity, by a stimulus-induced phase resetting of ongoing oscillations, or by some combination of both effects (Min et al., 2007). Many previous studies have used the intertrial coherence (ITC) (see Materials and Methods) as evidence for phase resetting, but ITC is sensitive to both additive and phase effects, making results difficult to interpret (Sauseng et al., 2007; Martínez-Montes et al., 2008). Further, mechanistic studies of the trial-averaged response have generally been limited to measurements from single recording sites without considering the spatial patterns of evoked or induced activity. Here, we investigated the mechanisms of the trial-averaged response by computing evoked and induced activity in individual trials and comparing the spatial patterns of these measures with those in trial-averaged amplitudes.
To distinguish evoked additive effects from induced phase effects, we adapted the methods developed by Martínez-Montes et al. (2008). In the absence of stimulus, oscillations in different trials have random phase, creating a uniform cloud of wavelet coefficients (Fig. 7A). Each wavelet coefficient represents the amplitude and phase of an oscillation at a single time point. Averaging across trials generates an oscillation with amplitude and phase given by the center of the wavelet coefficient cloud (indicated by the red dot in Figure 7A). This center point (in the absence of stimulus) will have an amplitude close to zero when taken across many trials. Stimulus-driven amplitude increases move the cloud of coefficients away from the origin, increasing the ITC and mean amplitude of coefficients without changing the distribution of phase within the coefficient cloud (Fig. 7B). In contrast, phase resetting does not affect single-trial amplitudes but does cluster phase toward one direction, increasing the ITC and the trial-averaged amplitude (Fig. 7C). Phase resetting can be quantified by a statistic known as the R2 phase coherence (Martínez-Montes et al., 2008), which is increased by phase resetting but is unaffected by oscillation amplitude (see Materials and Methods). The mean trial amplitude therefore captures evoked additive effects and the R2 phase coherence captures induced phase resetting.
Comparison of possible generation mechanisms for trial-averaged activity with additive amplitude and phase resetting. A, No stimulus response. Left, Sample single-trial oscillations before and after stimulus onset at t = 0. Center left, Signal averaged across 100 random trials. Center right, Wavelet coefficients for all trials 100 ms after stimulus onset. Red dot indicates center of mass. Right, Cross-trial measures relative to prestimulus period. B, Same as A, but with stimulus evoking a reliable waveform that is added to background oscillations. The cloud of wavelet coefficients is shifted in space without changing in shape, creating an increase in intertrial coherence and mean amplitude measures. C, Same as A, but with stimulus inducing a phase reset of existing oscillations. Wavelet coefficients are aligned to a narrow range of phases, creating an increase in intertrial coherence and R2 phase coherence.
We calculated R2 coherence and mean amplitude relative to prestimulus levels for each electrode at each point in time and then averaged across all electrodes and all recordings. All recordings and all frequency bands exhibited a significant overall increase in both single-trial amplitude (Fig. 8A) and R2 phase coherence (Fig. 8B). We then constructed separate spatial maps of R2 phase coherence and mean amplitude in each frequency band (Fig. 9A). We calculated the correlation coefficient between these maps and the trial-averaged amplitude map at each time step. This procedure allowed us to establish whether induced or evoked activity patterns better explained the trial-averaged amplitude patterns and thus to infer which mechanism was more important in producing trial-averaged activity.
Relative strength of single-trial mechanisms for trial-averaged activity. A, Percentage increase in single-trial amplitude from average prestimulus level averaged across all channels and recordings. Shading indicates SEM. B, Percentage increase in R2 phase coherence from average prestimulus level averaged across all channels and recordings. Lines follow legend in A. Shading shows SEM.
Single-trial mechanisms of trial-averaged amplitude patterns. A, Example of 10 Hz trial-averaged amplitude spatial maps compared with single-trial amplitude and phase coherence maps in one recording. Each row of snapshots is normalized to its maximum value. B, Average spatial correlation coefficients for theta-band oscillations between trial-averaged amplitude and mean single-trial amplitude (Mean amp) maps and between trial-averaged amplitude and R2 phase coherence (R2) maps across all recordings. Shading shows SEM. C, Average spatial correlation coefficients by frequency band between trial-averaged amplitude maps and mean single-trial amplitude maps at the same frequency across all recordings. D, Same as C, but showing correlations between trial-averaged amplitude maps and R2 phase coherence maps.
We first examined theta-band oscillations and found that correlations between trial-averaged spatial activity patterns and mean single-trial amplitudes were greater than those for R2 phase coherence in all recordings (Fig. 9B). These same effects were also present in the other frequency bands tested (Fig. 9C,D). These analyses show that additive activity is the dominant contributor to trial averaged amplitude patterns. We conclude that the source and sink patterns that we see in trial-averaged activity primarily reflect stimulus-evoked power increases. Because phase waves with the same coherent motion direction can nevertheless vary in space and time across trials, they do not consistently summate at individual electrodes, so their direction sensitivity is masked by trial averaging.
Discussion
We have shown that visual stimuli generate propagating activity patterns in area MT of marmosets and that the propagation direction of these patterns is sensitive to stimulus motion direction. Our results therefore reveal a functional link of propagating patterns to visual processing. The direction-sensitive propagating patterns are not present in the trial-averaged response, which instead shows robust expanding (source) and contracting (sink) patterns that are consistent across stimulus directions. This difference arises because the trial-averaged response is dominated by stimulus-evoked increases in power, which do not capture the propagating phase patterns seen in single trials and render the averaged response less sensitive to stimulus direction.
Trial-averaged spatiotemporal patterns
We found that the trial-averaged LFPs reliably showed a spreading source pattern shortly after stimulus onset, followed by a contracting sink pattern centered up to 2 mm away. An evolution from source to sink patterns has also been found in population responses to visual and somatosensory stimuli in mouse cortex using voltage-sensitive dye imaging (Mohajerani et al., 2013). That study found that source/sink patterns occurred on a shorter time scale (∼50 ms duration) and larger spatial scale (up to 4 mm between source and sink, activity across whole hemisphere) than the patterns that we observed, but the fundamental dynamics are similar despite different species and recording modalities. This similarity suggests that spreading source activity followed by contracting sink activity could be a general property of event-related cortical population responses. We observed that the source–sink patterns were largely unaffected by stimulus motion direction and suggest that their structure and dynamics are mostly determined by the topographic structure of sensory cortices. It is also possible that the observed source–sink patterns do not represent true wave propagations, but instead reflect cortical sites responding to stimulus at different latencies, as was shown previously in voltage-sensitive dye imaging in V1 (Sit et al., 2009).
Evoked and induced mechanisms of trial-averaged patterns
There has been a longstanding debate as to whether event-related responses such as evoked potentials are caused by an additive increase in power or by phase resetting of ongoing oscillations (Makeig et al., 2002; Penny et al., 2002; Min et al., 2007). Much of the debate, however, rests on recordings of event-related potentials in single EEG channels. We approach this debate differently, by analyzing the spatiotemporal dynamics of phase resetting and additive amplitudes. Comparing the strength of phase resetting and additive amplitude effects directly is generally not possible without carefully controlled stimuli (Xu et al., 2016) and, indeed, when we used single-electrode measures, we could identify the presence of both mechanisms but could not compare their strength. We then introduced a correlation-based method that allows comparison of the relative impact of evoked and induced effects in any spatially extended recording by incorporating the spatial structure of activity into previous measures (Martínez-Montes et al., 2008). Our new methods allow us to show that both stimulus-evoked additive amplitude and stimulus-induced phase resetting effects contribute to event-related patterns, but with different impact—for amplitudes, the average maximum spatial correlation was rav = 0.6 and, for phase resetting, it was rav = 0.15. Therefore, both mechanisms contribute to the formation of the event-related population response (Min et al., 2007), but trial-averaged patterns are dominated by additive amplitude effects.
Visual discrimination by propagating patterns
We demonstrate a relationship between stimulus motion direction and the propagation direction of LFP activity in cortex. Further, we show by objective classification with SVMs (Fig. 5) that propagation direction can support visual motion discrimination. We note that the propagating pattern measures used to train the SVMs collapse across all electrodes and therefore represent a significant reduction of dimensionality from the original recordings. The successful decoding of stimulus motion using such low-dimensional information further indicates that mean propagation direction is relevant for stimulus motion discrimination. Finally, we observed no difference in discriminability between opposite and orthogonal stimulus directions (Fig. 5), implying that discrimination would become more difficult only at smaller direction differences.
Previous studies had shown that the appearance of a stimulus elicits traveling waves (Prechtl et al., 1997; Ferezou et al., 2007; Xu et al., 2007; Wu et al., 2008; Sato et al., 2012; Muller et al., 2014) and that these waves can improve stimulus discriminability (Agarwal et al., 2014). Few, however, have attempted to relate stimulus or behavioral features to the structure of waves [as exceptions, the amplitude of traveling waves has been linked to saccade direction in macaque area V4 (Zanos et al., 2015) and arm movement in macaque motor cortex (Rubino et al., 2006)].
We observed complex propagating patterns in the phase of single-trial LFP oscillations. These phase patterns very likely reflect true propagating waves instead of stimulus-driven latency effects because they are highly variable from trial to trial instead of following static response patterns and display coherent dynamics even in the absence of stimulus. We only rarely observed global propagating waves that swept uniformly across the recording array, as has been reported previously (Rubino et al., 2006; Muller et al., 2014). This difference likely arises because we used a high-contrast stimulus presented over a wide region of visual space. Such stimuli are known to suppress planar traveling waves in visual cortex (Sato et al., 2012). In addition, sleep and anesthesia are associated with slow fluctuations in cortical excitability (Vanhatalo et al., 2004; Cheong et al., 2011; Pietersen et al., 2017), which may influence the trial-to-trial wave timings and dynamics. In future studies, it would be interesting to determine whether propagating waves generated by small, low-contrast stimuli have different sensitivity to stimulus features than the varied propagating activity patterns reported here and to examine the relationship between stimuli and cortical propagating patterns in different behavioral states.
Neurons in area MT have very large and overlapping receptive fields, so a coordinated population response may be important in stimulus encoding (Chen et al., 2015). One way to achieve this coordination is through synchronous feature binding (Singer and Gray, 1995), in which neurons excited by features of the same object fire together. Consistent with previous studies (Thiele and Stoner, 2003; Palanca and DeAngelis, 2005), we did not find stimulus dependence of LFP synchrony measures. Instead, we observed complex wave patterns that could help to discriminate motion direction for a few hundred milliseconds after stimulus onset. We speculate that these transient propagating patterns may serve to coordinate neural responses and communicate stimulus-related information to different brain regions (Gong and van Leeuwen, 2009; Sato et al., 2012). Propagating wave patterns are widely observed, occurring over several spatial scales throughout multiple brain regions (Rubino et al., 2006; Wu et al., 2008; Sato et al., 2012; Muller et al., 2016), so the pattern-based representation revealed by our study could be of general applicability to understanding neural representations.
Footnotes
This work was supported by the Australian Research Council (Grants DP160104316 and CE140100007). We thank S.C. Chen, N. Zeater, S.K. Cheong, and A. Pietersen for help with data acquisition.
The authors declare no competing financial interests.
- Correspondence should be addressed to Dr. Pulin Gong, School of Physics, The University of Sydney, NSW 2006, Australia. pulin.gong{at}sydney.edu.au
