MEG sensor and source measures of visually induced gamma-band oscillations are highly reliable

High frequency brain oscillations are associated with numerous cognitive and behavioral processes. Non-invasive measurements using electro-/magnetoencephalography (EEG/MEG) have revealed that high frequency neural signals are heritable and manifest changes with age as well as in neuropsychiatric illnesses. Despite the extensive use of EEG/MEG-measured neural oscillations in basic and clinical research, studies demonstrating test–retest reliability of power and frequency measures of neural signals remain scarce. Here, we evaluated the test–retest reliability of visually induced gamma (30–100 Hz) oscillations derived from sensor and source signals acquired over two MEG sessions. The study required participants (N = 13) to detect the randomly occurring stimulus acceleration while viewing a moving concentric grating. Sensor and source MEG measures of gamma-band activity yielded comparably strong reliability (average intraclass correlation, ICC = 0.861). Peak stimulus-induced gamma frequency (53–72 Hz) yielded the highest measures of stability (ICCsensor = 0.940; ICCsource = 0.966) followed by spectral signal change (ICCsensor = 0.890; ICCsource = 0.893) and peak frequency bandwidth (ICCsensor = 0.856; ICCsource = 0.622). Furthermore, source-reconstruction significantly improved signal-to-noise for spectral amplitude of gamma activity compared to sensor estimates. Our assessments highlight that both sensor and source derived estimates of visually induced gamma-band oscillations from MEG signals are characterized by high test–retest reliability, with source derived oscillatory measures conferring an improvement in the stability of peak-frequency estimates. Importantly, our finding of high test–retest reliability supports the feasibility of pharma-MEG studies and longitudinal aging or clinical studies.

Given its potential role in routing information within brain networks (Fries, 2015) rhythmic activity at gamma frequencies in response to visual stimuli has received considerable interest in both invasive and non-invasive electrophysiology. Specifically, visually induced brain responses at gamma frequencies have been shown to vary with stimulus properties (Jia et al., 2013;Perry et al., 2013; see also Box 1 in Tan et al., 2013) such as contrast (Ray and Maunsell, 2010), stimulus type (Hermes et al., 2014;Hermes et al., 2015), stimulus velocity (e.g. Friedman-Hill et al., 2000;Gray et al., 1990;Lima et al., 2011;Muthukumaraswamy, 2013), and the temporal expectation of reward (Lima et al., 2011). Moreover, the peak frequency of visually induced gamma brain oscillations shows high heritability in monozygotic twins (van Pelt et al., 2012), with concordances comparatively lower for heterozygotic twins and lowest between non-related individuals.
In addition to their role during normal brain functioning, visuallyinduced oscillations have been investigated in several neuropsychiatric disorders, such as schizophrenia and autism spectrum disorders, as a means of deriving insights into the underlying circuit dysfunctions (e.g. Sun et al., 2013;Tan et al., 2013). Collective experimental and theoretical studies provide evidence that cortical gamma-band activity is predominantly generated through rhythmic synaptic inhibition, which temporally coordinates windows of excitability in principal cells (Bartos et al., 2007;Buzsáki and Wang, 2012;Wang, 2010;Whittington et al., 2000). These recurring frames of excitability provide an efficient and elegant means of organizing and coordinating functional cell ensembles for neural communication (Akam and Kullmann, 2014;Buzsáki and Watson, 2012;Buzsáki, 2010). Along with mutually connected (GABA-mediated) inhibitory neurons (Traub et al., 1996;Wang and Buzsáki, 1996), networks of recurrent (AMPA-mediated) excitatory-inhibitory neurons also contribute to the genesis of gamma oscillations (Brunel and Wang, 2003;Tiesinga and Sejnowski, 2009;Wang, 2010). Importantly, changes in cellular parameters have been found to affect the coordination of excitatory and inhibitory processes (Bernard et al., 2000;Hashimoto et al., 2008;Lewis et al., 2005) and may lead to impairments in generating high-frequency oscillations (Ben-Ari et al., 2004;Ben-Ari et al., 2012;Gonzalez-Burgos and Lewis, 2008;Schnitzler and Gross, 2005;Uhlhaas and Singer, 2012).
Despite prevalent use of visually-induced gamma-band responses during normal brain functioning and in neuropsychiatric populations, evidence for the reliability of EEG/MEG-derived oscillatory measures remains scarce. However, for measures of high-frequency neural activity to be useful neurophysiological "spectral fingerprints" (Siegel et al., 2012) and potentially serve as biomarkers or endophenotypes in clinical research, it is essential that these parameters are highly reliable and robust. To the best of our knowledge, there are only two studies that systematically assessed the test-retest reliability of visually induced neural oscillations from EEG/MEG signals. The study by Muthukumaraswamy et al. (2010) revealed high intraclass correlation (ICC) values (0.8-0.98) for spectral measures of visually elicited gamma (40-60 Hz) band responses to static visual gratings across multiple assessments for source-derived MEG-signals. Specifically, they found highest repeatability for peak frequency, followed by its corresponding spectral signal change and bandwidth. In addition, Fründ et al. (2007) observed that the magnitude and frequency of participants' EEG responses to large (vs. small) foveally presented static visual gratings were strongly (Pearson's) correlated between sessions, highlighting the dependency of visually induced gamma activity on stimulus properties.
To further examine the test-retest reliability of visually-induced responses in neuromagnetic data, we employed a modified moving visual stimulus protocol (van Pelt et al., 2012) developed by Hoogenboom et al. (2006) that has been shown to yield gamma-band responses with robust signal-to-noise-ratio (SNR). Specifically, the current study assessed estimates of peak frequency, spectral modulation and spectral bandwidth across measurements to comprehensively assess the reliability of visually-induced responses from both sensor and source derived MEG signals. Additionally, we evaluated the SNR for sensor vs. source estimates of visually-induced high-frequency activity, which is a question relevant to both basic and applied MEG research.

Methods
The University of Glasgow College of Science and Engineering Ethics Committee approved the experimental protocol, in which the present experiment was part of a battery of sensory processing tasks performed during each MEG session.

Participants
Fourteen healthy participants (4 Females; mean age (±SD) = 25 (±4) years) took part in our study to assess the test-retest reliability of visually induced high frequency neural oscillations over two MEG sessions (range 1-11 days; mean (±SD) = 4 (±3) days apart). Participants were recruited from the University of Glasgow School of Psychology participant pool, provided informed consent prior to the experiment and were compensated (at the standard rate of £6/h) for their time. All participants were right handed (Edinburgh Handedness Test; Oldfield, 1971), characterized by normal or corrected vision and had no known neurological disorders.
Prior to each MEG session, scheduled at the same time of the day, each participant filled in a brief questionnaire which assessed differences in caffeine intake, smoking habits, alcohol consumption, and hours of sleep as well as general well-being prior to each measurement. As previously reported (Tan et al., 2015) participants' responses did not differ across sessions. Additionally, female participants took part in the study within the first 5-10 days during the follicular phase of their menstrual cycle in both MEG sessions to control for potential influence of hormonal fluctuations (Epperson et al., 2002).

Stimuli and task
We employed a foveally presented moving visual grating stimulus (Supplementary SFig. 1) that has been observed to induce robust MEG gamma band response in the human visual cortex (Hoogenboom et al., 2006). The visual grating was presented at a viewing distance of 186 cm in front of the seated participants. Each trial began with the presentation of a central fixation spot (Gaussian diameter: 0.5°) for 500 ms. The contrast of the fixation spot is subsequently reduced by 40% for 1500 ms, indicating the upcoming presentation of the moving circular sine wave gratings. The ensuing concentric moving grating (2.7 cycles/°; contrast: 100%; 5°visual angle; velocity: 0.75°/s) contracted towards the fixation spot and accelerated (velocity: 1.2°/s) randomly between 750-3000 ms post grating presentation onset. Participants were tasked to indicate the detection of this acceleration with a button press within 700 ms of its occurrence. Each trial lasted for~4-6 s long and during the inter-trial interval (1000 ms) participants were provided with feedback as to whether the speed of their response was adequate, too fast or too slow. Rare incidences (10%) in which no acceleration occurred were interspersed within a sequence of 80 trials that made up a block of the visual task. We provided performance accuracy feedback during the short break after each of the 3 task blocks.

Neuroimaging acquisition
MEG data were acquired using a 248-channel magnetometer system (MAGNES® 3600 WH, 4D-Neuroimaging, San Diego) while participants engaged in the task, sitting upright within an electromagnetically shielded room. For each participant, a suitable MEG seat position was determined and marked during the first session. Every attempt was taken to keep this seat position and the MEG system's helmet (housing the SQUID sensors) in the same configuration prior to each acquisition so as to minimize the variance of participants' head and sensors' positioning across runs and sessions. Head position stability was assessed before and after each acquisition run via five indicator coils attached relative to the (left, right preauricular and nasion) fiducials, and were codigitized with participants' head-shape (FASTRAK®, Polhemus Inc., VT, USA) for subsequent co-registration with individual MRI (1 mm 3 T1weighted; 3D MPRAGE). The MEG, touch-pad response (LUMItouch™, Photon Control Inc., BC, Canada) and eye-tracker (EyeLink 1000; SR Research Ltd., Ontario, Canada) signals were sampled synchronously at 1017.25 Hz, with online 0.1 Hz high-pass filtering.

MEG data processing
All data processing and analyses were performed using Fieldtrip Toolbox functions (http://fieldtrip.fcdonders.nl; Oostenveld et al., 2011) and additional scripts developed within MATLAB® (The MathWorks, Natick, MA). Faulty sensors (mean (± SEM) = 10 ± 2 per session, visually identified) with large signal variance or whose signals were flat were removed and interpolated using nearest-neighbor averaging procedure. One MEG-measurement was corrupted by global noise and technical issues during one of the two acquisition sessions. Accordingly, this participant was excluded from the analyses reported here (i.e. N = 13).
Raw MEG signals from correctly responded trials were epoched from −1500 to +2000 ms relative to grating stimulus onset (0 ms), with linear trends removed, power-line (50 Hz) notch-filtered, and 'de-noised' relative to reference MEG channel signals. Raw trials were visually inspected and trials with obvious artifacts (muscle, squid jumps etc.) were excluded. Subsequently, Independent Component Analysis was used to isolate and to reject ocular-movement and cardiac components from the MEG signals, yielding on average (±SEM) 180 (±7) artifactfree trials for each participant and session.

Time-frequency analysis on sensor and source signals
For sensor level analysis, artifact-free neuromagnetic time series were transformed to planar gradient signals (Bastiaansen and Knosche, 2000) prior to time-frequency analyses and subsequently recombined. Similar to previous work (e.g. Hoogenboom et al., 2006;van Pelt et al., 2012) we focused our sensor-level analysis on the spectral power time-series derived from the 23 parieto-occipital sensors At the source level, prior work (Hoogenboom et al., 2006;Muthukumaraswamy et al., 2010) and preliminary assessment of source-level data indicated that strongest signals were generated within the calcarine, visual lingual, and occipital areas. Given that prior findings have consistently reported visually-induced cortical sources significantly associated with the high frequency oscillations, average signals from visual cortical regions (bilateral calcarine, cuneus, lingual, superior, mid and inferior occipital cortical areas) were initially Fig. 2. Summary of all participants' maximally modulated voxels in response to moving visual grating stimulus for both sessions. Maximally modulated voxels shown in (i) axial, (ii) coronal and (iii) sagittal views within 'glass' brain volume. Each participant's maximally modulated voxel is color-coded as in Figs. 4-6 with square and diamond markers denoting non-parametric T-stats derived maximum in sessions 1 and 2, respectively. (iv) 3-D brain volume view of maximally modulated voxels surrounded by their neighboring 26 voxels whose signals were incorporated in the time-frequency analyses. Refer to Methods for further details.   Fig. 1 (ii) ). However, this led to a loss of signal power, which suggested source-specificity in the stimulus-induced modulation (Supplementary SFig. 3 (iii) , see also SFigs. 5-7). We therefore determined the maximally induced voxel for each individual for each MEG session.
To this end, a 30-100 Hz broadband time-frequency decomposition (frequency centered at 65 Hz; ±35 Hz taper smoothing; 10 ms temporal resolution; 50% overlapping with 500 ms time window) was performed on the artifact-free epochs prior to the derivation of common source spatial filters using the DICS inverse-solution algorithm (Gross et al., 2001). Subsequently, bootstrap resampling source statistics was performed (with 1000 Monte Carlo repetitions) between stimulus duration of interest (StimDur) and baseline (500-2000 ms and −1500-0 ms relative to moving grating onset, respectively) to determine individual participants' maximal source statistic (see Fig. 1 (iii) ; Supplementary SFig. 2). For each participant and session, virtual sensors' signals were extracted from the maximally modulated source (FDR-corrected, alpha = 0.05) as well as its corresponding 26 surrounding voxels using individual MNI-normalized source model grid (6 mm resolution).

Derivation of oscillatory parameters and reliability analysis
For both sensor and source derived signals multi-taper fast-fourier time-frequency decomposition (± 2 Hz taper smoothing and 10 ms temporal resolution; 2 Hz resolution from 30 Hz to 100 Hz) was performed with 50% overlapping 500 ms time window on the artifactfree epochs. All spectral power time-series were expressed as relative change to baseline (from −1500 ms to 0 ms prior to visual grating stimulus onset; Supplementary SFig. 4; Fig. 3). Induced sustained spectral modulations were averaged over the period of 500 ms-2000 ms post grating stimulus onset (see Supplementary SFig. 4; SFig. 3) for each frequency interval within the 30-100 Hz range; avoiding transient visual response onset and any preparatory behavioral responses.
A 1st order Gaussian fit was performed on these time-averaged spectral time-series (Campbell et al., 2014;Haegens et al., 2014) to determine (a) peak response frequency, and the corresponding (b) signal change modulation and (c) bandwidth (i.e. by deriving the full-width-at-halfmaximum; FWHM) at this peak frequency (Supplementary SFig. 4). As in previous research (e.g. Hoogenboom et al., 2006), in cases where a participant manifested double gamma-band peaks (e.g. S02), the higher gamma-band peak was selected for subsequent analysis.
Adopting a similar approach to previous reliability assessments (Muthukumaraswamy et al., 2010;Tan et al., 2015) we calculated the intraclass correlation (ICC; Shrout and Fleiss, 1979) using Matlab Central file-exchange ICC.m function (A. Salarian 2008; implemented with statistical testing based on McGraw and Wong, 1996a;McGraw and Wong, 1996b) to assess the degree of consistency of these spectral variables. Defined as the ratio of between-subject variance and the total variance, ICC assesses the reliability of the repeated measures of an individual's oscillatory parameters by comparing the between-measures variability of each individual to the total variation across all measures and participants. An ICC value of 1 indicates perfect within subject reliability of neural oscillatory measures derived on differing occasions from the same participants, while ICC of 0 indicates no reliability. ICCs were assessed for both sensor and source-derived neuromagnetic parameters. The distributions of parameters of interest were similar across sessions (i.e. insignificant 2-sample Kolmogorov-Smirnov tests), although those of signal change at sensor level and spectral bandwidth at source level from session 1 were marginally skewed (as determined by Lilliefors test of normality). For appropriate application, these, together with their corresponding distributions from session 2, were square root transformed prior to ICC assessments.

Distributions of peak frequency and signal change within visual cortical regions
The derivation of oscillatory parameters within each participant's set of 12 AAL-parceled (calcarine, cuneus, lingual and occipital gyri;~800 voxels including both hemispheres) visual cortical regions were repeated to further assess the distribution of peak frequency and corresponding signal change within visual cortical regions for each session.

Results
Participants demonstrated high response accuracy; mean (±SD): 88 (± 5) % and 93 (± 4) % for sessions 1 and 2, respectively. On average (±SD), their response times were 554 ± 71 ms and 536 ± 60 ms, for sessions 1 and 2, respectively. While response accuracy improved on the second session (t 12 = − 4.451; p b 0.001), mean response times did not differ between sessions (t 12 = 1.417; p = 0.182). Maximally modulated visual voxels from all participants were distributed predominantly within the primary visual lingual, calcarine, and occipital areas (Fig. 2). Table 1 lists the SPM coordinates and anatomical labels of each participant's maximum voxels for both sessions. The locations of participants' maximally modulated voxels are not identical between sessions, but are mostly clustered within neighboring voxels (6 mm resolution voxels). The mean (± SD) intra-participant spatial variability of maximally modulated voxel location in our participant dataset is 21 (±15) mm.
Spectral changes during stimulus presentation for both sessions are shown in Fig. 3 for sensor (A) and source (B) derived MEG signals for all participants. For most participants, the induced gamma frequency response is sustained from about 350 ms post moving grating stimulus  onset until end of stimulus presentation and confined to a frequency range of 50-80 Hz. The induced visual gamma responses for both sessions from either sensor or source derived signals show good resemblance. Spectral responses derived from maximally modulated visual ROIs were generally much stronger, reaching a maximal of~450% signal change compared to sensor estimates (~200%). For participants with stronger spectral change relative to baseline (N 40%) at the sensor level, the source-derived frequency response revealed stronger modulations (100-400%). Participants for whom spectral change was weak at sensor level (b25%), the induced spectral modulations were similarly recovered at source level. For a quantitative measure of repeatability of visually induced gamma response, we assessed the retest reliability of peak frequency, the corresponding spectral signal change and response frequency bandwidth, as determined by the full-width-at-half-maximum (FWHM) of the 1st order Gaussian fit (Supplementary SFig. 4; see Methods). These spectral parameters are summarized in Table 2 and Fig. 4 for both MEG sessions and their respective grand-average time-frequency plots for each session (A, B) are shown in Supplementary SFig. 3. Individual peak frequencies ranged approximately 53-72 Hz for both sensor and source derived spectral measures (Fig. 4A (i) , B (i) ) with an average (± SD) of~61 (± 5) Hz. The signal change ranged between 8-200% (mean (±SD): 74 (±51) %) and 12-450% (mean (±SD): 175 (±114) %) for sensor and source derived peak frequencies (Supplementary SFig. 4A (ii) , B (ii) ), respectively. The bandwidth of the peak frequency spanned the range between 12 and 46 Hz for sensor and source derived peak frequencies (Fig. 4A (iii) , B (iii) ) were 23 Hz (±7) Hz on average (±SD).
Sensor as well as source-derived peak frequency and spectral resolution did not differ significantly within participants between sessions (Fig. 5A,B (i; iii) ). While within participant power modulation did not differ significantly between sessions at the sensor level (Fig. 5A (ii) ), corresponding source spectral modulation was significantly larger in the 2nd compared to the 1st session (t 12 = − 3.391, p b 0.001), and this  .g. Fig. 4A, B). difference appeared to be driven by a subset of participants (Fig. 5B (ii) ). We observed significant within participant sensor vs. source differences (t 12 = − 3.370; t 12 = − 3.306, for sessions 1 and 2 respectively; p b 0.01) for peak frequency of induced oscillations despite small average differences. In line with observed spectral responses (Fig. 3), power modulations corresponding to peak frequencies were significantly lower (p b 0.01) at sensor compared to source level (t 12 = − 3.391; t 12 = − 3.858, for sessions 1 and 2 respectively; Table 2). Finally, there were no significant differences in the corresponding peak frequency bandwidth derived from sensor and source for each session.
ICC values for peak frequency, corresponding spectral signal change and peak frequency bandwidth measures indicated overall strong reliability (mean ICC = 0.861; ICC range: 0.622-0.966; p b 0.001; refer to   Fig. 6) and comparable sensor-(mean ICC = 0.895) and source-derived (mean ICC = 0.827) oscillatory measures. ICC values were highest for peak frequency, followed by corresponding spectral signal change and frequency bandwidth, for both sensor and source derived measures. Assessment of peak frequency yielded strong reliability (ICC N 0.90) for both sensor (ICC PEAK-FREQ = 0.940; p b 0.0001) and source (ICC PEAK-FREQ = 0.966; p b 0.0001) derived peak frequencies ( Fig. 6A (i) , B (i) ). Similarly, we observed high ICC values for both sensor (ICC SPECTRAL = 0.890; p b 0.0001) and source (ICC SPECTRAL = 0.893; p b 0.0001) derived spectral modulations at peak frequencies ( Fig. 6A (ii) , B (ii) ). The bandwidth of peak frequency (Fig. 6A (iii) , B (iii) ) yielded higher reliability for sensor (ICC FREQ-RES = 0.856; p b 0.0001) compared to source derived (ICC FREQ-RES = 0.622; p b 0.01) signals.
Beyond the mean cluster signals derived from our maximally modulated visual voxels, we observed (Supplementary Figs. SFig. 5 (i) and SFig. 5 (ii) ) narrow ranges of peak frequencies (SFig. 5 (i)(ii) , A) across visual cortical regions in the range comparable to those found for the maximally task-induced voxel clusters (Fig. 4A, B (i) ). We also noted a consistent trend that the majority of visual voxels have lower spectral signal changes, with a smaller set of voxels yielding the largest signal changes (SFig. 5 (i)-(ii) , B). Supplementary Figs. SFig. 6 (i) and SFig. 6 (ii) provide a 3D overview of the variability of peak frequencies derived for each participant's~800 AAL-labeled visual cortical voxels in both sessions (SFig. 6 (i)-(ii) , A). We noted the variability in signals of correspondent AAL-parceled voxels between sessions and the overlap of maximally task induced voxel clusters with visual cortical voxels that generally manifested high signal changes (SFig. 6 (i) , SFig. 6 (ii) , B, C). Furthermore, assessment of peak frequencies and spectral signal changes indicated similar distributions within each AAL-parceled visual cortical region (SFig. 7 (i)-(xiii) ). Additionally, voxels manifesting the largest spectral signal changes tended to be found within the calcarine, lingual gyri, and on occasion, the occipital lobules.

Discussion
The present study assessed the test-retest reliability of visually induced high-frequency (30-100 Hz) oscillations derived from MEG sensor and estimated source signals. Overall, estimates of individual peak gamma frequencies, spectral modulation and peak frequency bandwidth were remarkably stable across measurements. Although the spectral signal change was stronger in the repeated session for sourcederived peak frequency spectra, both sensor and source spectral modulations exhibited comparably high repeatability. Peak oscillatory frequency yielded highest measure of reliability followed by its corresponding spectral signal change and peak frequency bandwidth. These stability measures of induced gamma oscillatory activity corroborate those repeatability assessments at source level reported by Muthukumaraswamy et al. (2010) for static visual stimulus. Importantly, the current study suggests that spectral estimates of both sensor and source derived parameters are both highly reliable. In addition, the current study suggests that the source-space approaches significantly improve SNRs of high frequency oscillations.
Consistent with previous findings (Hoogenboom et al., 2006;Muthukumaraswamy et al., 2010;Schwarzkopf et al., 2012), measures of visually-induced high frequency oscillations manifested individual variability, particularly for peak frequency and its corresponding spectral modulation. Previous study by Schwarzkopf et al. (2012) revealed higher peak gamma frequency being associated with larger primary visual cortical surface area, and suggested that differences in visually induced peak gamma frequency could be attributed to individual differences in the structural and functional architecture of visual cortex. Interestingly, cortical environments that are highly similar in e.g. receptor density, cytoarchitecture and/or coupling strength are thought to enable greater consistency in oscillatory activity (e.g. Breakspear et al., 2010). These observations led Schwarzkopf et al. (2012) to further suggest the role of lateral intra-areal inhibitory processes (Alitto and Dan, 2010;Edden et al., 2009) in sharpening sensory responses. Relevantly, the narrow range of individual peak gamma frequency observed across the extent of visual cortical areas assessed may arise from synchronyenhancing mechanisms of interneuronal dendritic gap junctions in spatially extended interneuron networks (Traub et al., 2001). On the other hand, theoretical studies (Cannon et al., 2014;Serenevy and Kopell, 2013) have indicated that heterogeneous cell properties and their connections (e.g. those that allow inputs to arrive at target networks with a range of phases covering a large part of the gamma cycle) may facilitate network entrainment (e.g. through differential recruitment of fastspiking interneurons) that might otherwise not be feasible if driving phases were highly similar. It is therefore reasonable to infer that neural oscillations at any given frequency (and brain region) may arise from various mechanisms (e.g. Ainsworth et al., 2012).
As with previous study by Muthukumaraswamy et al. (2010) we observed that maximally modulated voxels do not necessarily manifest within the exact location in the repeated session. However, the 3dimensional Euclidean distance between maximal loci suggested close proximities given an imaging resolution of 6 mm. Furthermore, spectral modulation within the visual areas reflected some inhomogeneity, as averaging across all (a priori) visual voxels led to the loss of signal. This inference corroborated with findings from our assessment of peak frequencies of individuals' voxels within the visual cortical regions, revealing that a majority of voxels manifested lower spectral modulations relative to a small proportion of voxels exhibiting the strongest (within participants) spectral signal change. Additionally, corresponding AAL-parceled voxels revealed between session variability in spectral signal change. These observations highlight the challenges in extracting signals from source as well as the likely mixing of signals of differing strengths that are picked up by the MEG sensors. In particular, the strength of MEG source signals is dependent on its orientation, increasing as the source orientation deviates from the radial towards the tangential axis (Hämäläinen et al., 1993). These are likely to contribute to the variability in spectral signal change and bandwidth corresponding to the induced peak gamma frequency both across participants and between sessions. Additionally, varying beamforming algorithms and approaches employed in previous (e.g. Muthukumaraswamy et al., 2010;Schwarzkopf et al., 2012) and current research may have also contributed to the differences in the reported spectral modulation. Employing a broad-band time-frequency decomposition with multi-tapers to extract the maximally modulated voxels and subsequent finer-resolution source signals spectral decomposition, we were able to recover the spectral modulations observed in sensor-derived signals or exceeded them at source level.
Interestingly, we observed that gamma spectral modulation was relatively higher in the repeated session for many participants, particularly at the source level. It has previously been shown that gamma band activity is enhanced in neurons driven by attended stimuli (Fries et al., 2001) and associated with corresponding improvements in perceptual performance (e.g. Taylor et al., 2005;Womelsdorf et al., 2006). With repeated performance of the same visual task, our participants became faster (although not significantly so) and more accurate in detecting stimulus acceleration. More recently, neurophysiological studies by Ray et al. (2013) and Chalk et al. (2010) highlighted that where attention is called upon in a task, the modulation in gamma spectral power may be a reflection of the changes in the underlying excitationinhibition activity, which could be accounted for by normalization. The normalization model of attention posits that a neuron's (gamma) response is suppressed by the overall response of its neighboring neurons and predicts that attention increases its excitatory drive, which in turn increases normalization. Crucially, the findings of Ray et al. (2013) demonstrated that even with attentional load fixed, increased normalization, e.g. when a cell's receptive field processes its nonpreferred (vs. preferred) motion direction, led to an increase in gamma spectral power. It is conceivable that as the visual task becomes more familiar it calls upon less directed attention towards the center fixation spot for adequate performance. From this perspective, and the variable maximally modulated voxel location, subtle modulation in participants' focal attention and corresponding normalization could underlie the differing excitability that yielded varying induced gamma signal change between MEG sessions. Although beyond the scope and resolution of the present analysis, changes in signal strength between testing sessions might be induced through stimulus repetition. Recent study in primate visual cortex (Brunet et al., 2014) reported that repeated presentation of a pair of orthogonal visual gratings during a change detection task was associated with increased local and inter-areal gamma activity, and suggested the role of lateral inhibitory activity in sharpening the underlying stimulus representation. Interestingly, the authors further noted a corresponding increase in peak gamma frequency with stimulus repetition and observed enhanced synchrony, which is not observed in the present study.
Last but not least, it has been reasoned that source analysis may yield more reliable estimates of MEG-activity compared to sensor-derived measures because exact positioning of participants under the MEG sensors across repeated recordings is not a prerequisite. Furthermore, the beamforming approach in source analysis acts as spatial filters in suppressing background activity that may lower the reliability of sensor signals and thereby improve the signal estimates. Here, we replicated the findings observed in our previous reliability assessment of auditory steady state responses from sensor and source derived signals (Tan et al., 2015), which highlighted that even without available continuous head position information (e.g. Deuker et al., 2009) careful monitoring of head and sensors' positioning can yield highly reliable estimates of oscillatory measures from MEG sensor signals. We note that source estimated spectral modulation in some participants exceeded those derived from sensor signals. However, relative to sensor-derived estimates source-derived measures yielded less consistency in its corresponding spectral bandwidth. Importantly, our findings affirm that peak frequency measures of gamma-band brain oscillations can be very reliably estimated from sensor (van Pelt et al., 2012) and source (Muthukumaraswamy et al., 2010) derived MEG signals. Additionally, our assessments indicated that while spectral signal changes are more variable across visual cortical areas, individual peak gamma frequency manifested a narrow range across large areas of the visual cortex.

Conclusion
The present study further substantiates the view that MEG-derived oscillatory signatures of visual (and other sensory) responses are highly reproducible. This finding is important as MEG is increasingly used as a tool for the identification of biomarkers in clinical research (e.g. Georgopoulos et al., 2010;Sun et al., 2013) and for investigation of rhythmic activity during normal brain functioning. Specifically, our research underscores peak frequency of visually induced brain oscillations as particularly reliable. Given that gamma oscillation is generated by well-coordinated inhibitory and/or excitatory neuronal networks (Buzsáki and Wang, 2012) and peak frequency of induced visual gamma has been shown to decrease with age (Gaetz et al., 2012;Muthukumaraswamy et al., 2010), estimates of peak gamma frequency could be useful in tracking individual's underlying neural excitability (Ray and Maunsell, 2015) over time. No doubt, further studies are needed to better link macroscopic measures of gamma activity e.g. M/ EEG with quantifiable proxies of microscopic (e.g. molecular; network-level) processes. Nonetheless, these insights are particularly encouraging for larger-scale, longitudinal, and/or clinical studies that require repeated MEG measurements, and for which study outcomes are not contingent on source-derived oscillatory readouts.