Functional Connectivity Analysis of Transcutaneous Vagus Nerve Stimulation (tVNS) Using Magnetoencephalography (MEG)

Altered brain functional connectivity has been observed in conditions such as schizophrenia, dementia and depression and may represent a target for treatment. Transcutaneous vagus nerve stimulation (tVNS) is a form of non-invasive brain stimulation that is increasingly used in the treatment of a variety of health conditions. We previously combined tVNS with magnetoencephalography (MEG) and observed that various stimulation frequencies affected different brain areas in healthy individuals. We further investigated whether tVNS had an effect on functional connectivity with a focus on brain regions associated with mood. We compared functional connectivity (whole-head and region of interest) in response to four stimulation frequencies of tVNS using data collected from concurrent MEG and tVNS in 17 healthy participants using Weighted Phase Lag Index (WPLI) to calculate correlation between brain areas. Different frequencies of stimulation lead to changes in functional connectivity across multiple regions, notably areas linked to the default mode network (DMN), salience network (SN) and the central executive network (CEN). It was observed that tVNS delivered at a frequency of 24 Hz was the most effective in increasing functional connectivity between these areas and sub-networks in healthy participants. Our results indicate that tVNS can alter functional connectivity in regions that have been associated with mood and memory disorders. Varying the stimulation frequency led to alterations in different brain areas, which may suggest that personalized stimulation protocols can be developed for the targeted treatment of different medical conditions using tVNS.


I. INTRODUCTION
I N THE brain, information is processed and transmitted via networks of neurons that often involve many different areas.Connections between brain areas may be functionally linked but spatially separate and form the basis of an important part of brain network analysis called functional connectivity [1].Functional brain connectivity can be estimated by determining whether signals recorded in different areas of the brain are correlated.
Magnetoencephalography (MEG) is a non-invasive method of neuroimaging that detects the magnetic fields that are induced by neuronal currents using sensors placed outside the head.The data collected by the sensors can be projected into source space to estimate the neural activity in different regions within the brain.MEG has excellent temporal resolution and has full coverage of both cortical and deep brain areas which allows the real-time analysis of brain response to various stimuli [2], [3], [4], [5].There are multiple methods to measure functional connectivity using MEG data including a range of methods which measure the phase-synchronization between neural sources [6].
Abnormal connectivity within or between brain networks has been suggested to play a role in medical conditions such as schizophrenia, dementia, autism and depression [7], [8], [9], [10].In recent years, invasive brain stimulation methods such as deep brain stimulation (DBS) or vagus nerve stimulation (VNS) have been investigated for their potential to regulate dysfunctional connectivity in the brain [11], [12], [13], [14], [15], [16], [17], [18].However, due to the invasive nature of these techniques and the risks associated with implantation, they are not routinely offered.Non-invasive stimulation techniques that offer similar therapeutic effects whilst overcoming some of the challenges of invasive stimulation are therefore of great interest.
Transcutaneous vagus nerve stimulation (tVNS) is a noninvasive method of stimulating the vagus nerve via an afferent branch in the ear and has been found to have similar therapeutic effects to the invasive method in a number of clinical studies investigating various medical conditions [19], [20], [21], [22].In particular, tVNS has been investigated recently for its use as a therapy option for major depressive disorder (MDD) [23], [24], [25], [26], an increasingly prevalent condition that is associated with high morbidity and mortality.A number of neuroimaging studies investigating the mechanism of depression have suggested that the symptoms arise from a dysregulation of brain networks encompassing both cortical and limbic brain regions [27], [28].It is still unclear how exactly tVNS modulates brain activity, but one theory is that the anti-depressive effects may be due to tVNS altering the functional connectivity between these same networks which may potentially affect mood [29], [30].
There are three main networks thought to play a part in the pathology of depression [31], [32], [33], [34].In addition to depression, these sub-networks are also linked to chronic pain, mood, memory and reward centers in the brain and dysfunction within and between these sub-networks are linked to a number of medical conditions [35], [36], [37], [38], [39].
The default mode network (DMN) which is involved in emotional regulation and self-referential processing, incorporates the medial prefrontal cortex, posterior cingulate cortex, precuneus and angular gyrus [32].The DMN has also been associated functionally with the inferior parietal lobe, the parahippocampus, the hippocampus and the lateral temporal cortex [32].Increased connectivity has been observed in parts of the DMN in MDD compared with healthy controls [40], [41].This increased connectivity is thought to decrease after treatment with anti-depressant pharmaceuticals [42].
The central executive network (CEN) is often also called the frontal-parietal network and involves the lateral prefrontal cortex, the posterior parietal cortex, the frontal eye fields and dorsolateral prefrontal cortex.This network is linked to goal-directed behavior [43].A number of studies have observed decreased connectivity between the dorsolateral prefrontal cortex and other regions within the CEN in depression [32].
The salience network (SN) is thought to act as a switch between the DMN and the CEN and consists of the insula, the dorsal anterior cingulate cortex, the amygdala and the temporal poles, the putamen, thalamus and ventral striatum [44].It has been observed that decreased connectivity in the right insula is inversely correlated with the severity of depression symptoms, and also aligns with decreased connectivity between the DMN and the CEN, which supports the theory that the SN works between the other two networks [45].Decreased connectivity between the amygdala and a number of other areas both within and outside the SN were also observed in a number of studies of depression [32].One study suggested that decreased connectivity between the dorsolateral prefrontal cortex and several regions within the SN could predict non-response to antidepressant treatment [46].
There are currently very few studies that have combined tVNS with a neuroimaging technique [19], [47], [48], [49].Simultaneous imaging and stimulation is challenging due to the presence of large stimulation artifacts that arise from the electrical current being delivered by tVNS, which result in magnetic fields many orders of magnitude larger than the underlying brain activity.In our previous paper [50], we showed that stimulation artifacts can be removed from the data, which allowed reconstruction of underlying brain dynamics, and demonstrated a successful example of combined MEG and auricular tVNS.
Different stimulation frequencies of tVNS have also been suggested to lead to brain responses in different areas [50], [51].Understanding this effect may be key to unlocking the potential for tVNS to be used to target specific brain areas, with various stimulation frequencies being utilized to alter connectivity in the desired brain region or network.
In this paper, we investigated whether functional connectivity changes can be determined from simultaneous MEG and tVNS, and whether various stimulation frequencies lead to changes in the functional connectivity in different areas of the brain.In addition to considering all regions across the brain (whole-head), analysis was also restricted to regions of interest (ROI) to investigate whether changes in connectivity as a result of different stimulation frequencies could be observed in networks that have been associated with mood, memory and pain.The findings from this study suggest that tVNS has the potential to be a viable therapy option for various medical conditions and mood disorders.The results suggest that tVNS can be personalized using particular stimulation frequencies to lead to responses in specific brain areas or networks which may lead to the targeted treatment of particular medical conditions.

A. Data Acquisition
Seventeen healthy volunteers (10 female, 7 male, age range 19-61, mean age 37.8, all right-handed) were recruited and selected based on predetermined inclusion and exclusion criteria based on handedness and MEG and magnetic resonance imaging (MRI) suitability.All participants gave written informed consent.The studies were approved by the Ethics Committee of Swinburne University of Technology and the experimental sessions were carried out at the Swinburne Neuroimaging Facility, Swinburne University of Technology.
MEG data was collected in a magnetically shielded room using a whole-head 306-channel (102 magnetometers, 204 planar gradiometers) ElektaNeuromag TRIUX (Elekta Instrument AB Stockholm) located at Neuroimaging Facilities, Swinburne University.
The recording was performed in a dimly lit room and participants were asked to keep their eyes open for the study and focus on a small cross projected on a screen in front of them.MEG data was collected in continuous mode and sampled at 5kHz to allow capture of the fine structure of the short stimulation artifact.A band-pass filter was applied between 0.1 -1660 Hz.Structural T1-weighted MRI scans for co-registration with the MEG data were obtained for each participant with a Siemens Tim Trio 3T.

B. Transcutaneous Vagus Nerve Stimulation
The stimulation was applied simultaneously with MEG recording using a Digitimer DS5 bipolar current stimulator which was placed outside the magnetically shielded room and isolated from the mains to ensure the safety of the participant at all times.The output of the stimulator was connected to three disposable stick-on electrodes (Ag/AgCl, Ambu Neuroline, area 18 mm2) that were placed on three locations of the left ear (Figure 1).The stimulation electrode was placed on the cymba concha, and the return electrode was placed on the tragus, which are both regions that are believed to contain afferent vagus nerve projections [52].The sham electrode was placed on the earlobe, which is a location that has been shown to be free from vagus nerve projections [53], [54].The left ear was selected to avoid adverse cardiac effects, as the right vagus nerve has efferent fibers leading to the heart [55], [56].
Four stimulation protocols were used in this study: 1) 24 Hz Sham delivered to the ear lobe; 2) 24 Hz delivered to the cymba concha; 3) Pulse Frequency Modulation (PFM) delivered to cymba concha; 4) 1 Hz delivered to cymba concha.In Protocol 1 the stimulation was delivered to the ear lobe whilst in the other protocols (Protocol 2,3,4) the stimulation was delivered to the cymba concha.In Protocol 3 the stimulation was pulse-frequency modulated based on a modulation frequency of 6 Hz which delivered the same number of pulses and energy as in Protocols 1 and 2 but pulses were delivered irregularly based on the modulation frequency.This was included to investigate whether the brain would respond to the modulation frequency of this stimulation, or whether this stimulation leads to less neural adaptation [57].
Protocol 4 used a lower frequency of 1 Hz and was included based on a study by Straube et al. who noted this frequency was more successful than higher frequencies in the treatment of migraine [51].
The order in which the stimulation was delivered was randomized for each participant.Each stimulation protocol was delivered for 10 minutes whilst simultaneously recording using MEG, and in between each 10-minute block the sensors were heated to clear any trapped flux that may have been caused by the electrical stimulation.The stimulation was delivered using biphasic pulses (cathodic phase first, no interphase gap, 1 ms duration) and the current was personalized to each participant based on their individual pain threshold (range 0.2 -2 mA).

C. Stimulation Artifact Removal
A major source of artifacts in the raw data occurred from the electrical stimulation which was performed concurrently with the scan.These artifacts had amplitudes that were an order of magnitude larger than the underlying brain activity and lasted approximately 4 ms from the peak of each biphasic stimulus pulse (Figure 2).
To remove stimulation artifacts, firstly the raw data was passed through a temporal signal space separation filter (Max-Filter 2.1) to remove general noise and artifacts and allow for correction due to any head movements that had occurred during the scan [58].The part of the signal that was contaminated with stimulation artifact was then removed which corresponded to 4 ms from the onset of the stimulus pulse and 1 ms before the onset.The stimulus delivery was synced to the MEG data recording using an output from the DS5 stimulator.The signal was then interpolated using an autoregressive model which was implemented using the MATLAB fillgaps function and used a prediction-sequence length of 80 samples and a model order of 40 [59], [60].This interpolation method was applied to the maxfiltered data prior to down-sampling and filtering, for further details please see [50].
For consistency across all datasets, the same sequence of gap interpolations was applied to all data sets, which was based on the 24 Hz regular and PFM spacing of the pulses, in total removing 5 ms of data at 36 time points per second, corresponding to a loss of 180 ms out of each second of data.This was applied to all four protocols.
After the major sources of artifacts were removed, the data was down-sampled from 5 kHz to 500 Hz and then processed and filtered for further analysis.

D. Data Analysis in Source Space
The sensor level data was projected into source-space by applying a spatial filter, in the form of a linearly constrained minimum variance (LCMV) beamformer [61] to the interpolated and down-sampled data which was split into 4 second long epochs and filtered into six frequency bands of interest; Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.delta (2 -4 Hz), theta (4 -8 Hz), alpha (8 -13 Hz), beta (13 -30 Hz) and gamma (30 -80 Hz) along with a lower frequency band; sub-delta (0.5 -2 Hz) which was selected to observe any direct effects of the low frequency 1 Hz stimulation protocol.
Co-registration of individual structural MRI scans to MEG coordinates was performed by aligning the fiducial points to construct source volumes.A grid of 12,000 points with 5mm spacing was constructed across the whole brain by warping the Montreal Neurological Institute (MNI152) template brain to individual head spaces.The singe shell model was derived from individual head shapes and the lead field model was used to create the forward model along with a covariance matrix (regularization 0.05) to determine beamformer filter weights.These weights were applied to the sensor data to obtain source estimations.Time series were obtained for each source point using the single orientation synthetic aperture magnetometry method of beamforming, which calculates a single direction with the maximum power from the three orientations [62].

E. Connectivity Methods
Connectivity can be estimated by considering the relationship between time series from two brain regions.Using the spectral connectivity package on MNE Python [63], we selected a phase-based metric, weighted phase-lag index (WPLI) which is robust to spurious zero-phase connections potentially arising due to remaining stimulation artifact [64].
The change in connectivity between edges across the whole brain, based on different stimulation protocols was then calculated.Three protocol comparisons were constructed by subtracting the absolute values in the connectivity matrices for each comparison: 24 Hz -Sham, 24 Hz -1 Hz and 24 Hz -PFM.
In addition to calculating total connectivity changes across the whole brain, ROI single-seed analysis was performed separately on three sub-networks: the DMN, the CEN and the SN (Figure 3).For the DMN, seeds were placed in the medial prefrontal cortex, the posterior cingulate cortex, precuneus, angular gyrus, the parietal lobe, parahippocampus and hippocampus and the connections from these seed areas to other areas in the brain were calculated which allowed for increased power in the statistical analysis.For the CEN, seeds were placed in the prefrontal cortex.In the SN, seeds were placed in the anterior cingulate cortex, insula, amygdala, temporal pole, putamen, thalamus and accumbens.These areas were selected based on their implication in the pathology of mood and memory disorders [32].

F. Statistical Analysis
At the individual level, the values of connectivity difference between protocols were converted using a Fisher z-transformation to ensure that values were approximately normally distributed.
A non-parametric permutation procedure was performed at the group level using maximum statistics to control for the family-wise error rate (α = 0.05, two-tailed) to find threshold values for edges showing a statistically significant difference in connectivity.Visual inspection of p-values and the nulldistribution generated by the permutation procedure was also performed to look for outliers which may be indicative of an effect that is too weak to appear with the current statistical power.

III. RESULTS
On inspection of p-values and the null-distribution from maximum statistics, a number of slightly suprathreshold edges were observed with t-scores just over the threshold value.It was decided to add any outliers with a p-value of less than 0.05 that had previously not passed the significance test (Figure 4).Further outliers that had p < 0.05 were not included, although these edges may be considered significant if the power in the statistics was increased, or if a slightly less conservative statistical method was used.Particular suprathreshold areas of interest were noted in deep brain regions such as the insula, the hippocampus and brain-stem.

A. Whole Head
Initially, connections across the whole head were considered, without restriction to areas of interest.Connections from this analysis indicate particularly statistically strong results.In the 24 Hz -Sham contrast, decreased connectivity was noted in both the delta and theta bands, between areas of the DMN and other brain areas.Increased connectivity was observed in the higher frequency beta band between the SN and areas of the temporal lobe (Figure 5).
Increased connectivity between SN and DMN was also observed in the beta band of the 24 Hz -1 Hz contrast (Figure 6).Finally, in the 24 Hz -PFM contrast increased connectivity was observed in the gamma band between the DMN and temporal areas (Figure 7).

B. Sub-Networks
The analysis was restricted further to focus on areas included in the three sub-networks associated with depression to increase the power in the statistics.I for a summary of all significant edges.
In the 24 Hz -Sham contrast (Figure 5), altered connectivity was observed across all six frequency bands between and within various sub-networks, with altered connectivity being observed in the DMN in all of the frequencies apart from sub-delta band.In the delta band, connectivity within the DMN was observed to be both increased and decreased in multiple regions, with increased connectivity being noted between the hippocampus and the anterior parahippocampal gyrus, whilst decreased connectivity was observed between the hippocampus and the posterior cingulate cortex.In this  I for a summary of all significant edges.Fig. 7.All significant connectivity changes (p<0.05) in the 24 Hz -PFM comparison (both whole-head and ROI combined), *denotes edges with p values < 0.025.Please see Table I for a summary of all significant edges.frequency band, connectivity was also increased between the DMN and the SN.Decreased connectivity between the frontal lobe (CEN) and the parietal lobe (DMN) was also observed in the alpha band.In the gamma band, connectivity was increased both within the SN and within the DMN.The results from this contrast suggest that tVNS delivered at a frequency of 24 Hz to the cymba concha alters connectivity within and between all three of the sub-networks that have been implicated in depression when compared with the sham placement on the ear lobe at the same frequency.
In the 24 Hz -1 Hz contrast (Figure 6), the majority of connections were altered in the highest frequency band, gamma (30 -80 Hz).Here, increased connectivity was observed between the DMN and the CEN, notably between frontal areas (CEN) and the hippocampus (DMN), as well as increased connectivity within the DMN.Increased connectivity between the amygdala (SN) and both frontal areas (CEN) and temporal areas was also observed.This contrast suggests that the brain reacts differently to the different stimulation frequencies of tVNS, particularly in the higher frequency bands.These results may indicate that certain stimulation frequencies are more suitable for altering functional connectivity in the three subnetworks.
Finally, in the 24 Hz -PFM contrast (Figure 7), connections were altered in mainly the sub-delta and gamma bands.Connectivity was observed to both increase and decrease Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.from various regions of the DMN in the sub-delta band, with decreased connectivity being noted between the posterior cingulate cortex (DMN) and the insula (SN) and increased connectivity being observed between the DMN and other brain areas.In the gamma band, connections were increased in areas involved in the SN, DMN and temporal lobe.This result suggests that the brain responded differently to the various stimulation types.It appears that the brain responded more strongly to the PFM than 24 Hz in the sub-delta band (indicated by the decreased connectivity between and within the DMN), whereas the regular spacing 24 Hz increased connectivity across the three networks in the gamma band.This result suggests that different frequencies or pulse-spacing of stimulation should be considered for altering connectivity in particular brain areas or frequency bands.

C. Regions Showing Most Alterations in Connectivity
Over all three contrasts, connectivity was compared between 24 Hz and another stimulation frequency.Therefore, we can consider the effect of vagus nerve stimulation generally, relative to the 24 Hz.Several areas within these networks were noted to be highly connected and had a large number of altered connections to other areas in response to tVNS (Figure 8).In summary, these nodes show altered connectivity both within and between the three sub-networks, notably the SN and the DMN, and represent regions that appear to be highly modulated by vagus nerve stimulation.

IV. DISCUSSION
The aim of this study was to use MEG to measure brain functional connectivity in response to different stimulation frequencies of tVNS in healthy individuals.Using concurrent imaging allowed the direct visualisation of brain response to stimulation, and four protocols of stimulation were compared to investigate the effect of stimulation frequency on functional connectivity.
In addition to evaluating functional connectivity changes across the whole brain, the analysis was also focused on regions that have been linked to mood.The results from this study suggest that tVNS can be used to alter functional connectivity in regions involved in three sub-networks that have been linked to mood, memory and pain: the DMN, CEN and SN.It was observed that different frequencies of stimulation led to changes in functional connectivity in various brain areas.
Connectivity changes are reported as increases or decreases relative to the 24 Hz stimulation frequency.The findings suggest that the 24 Hz leads to more overall changes in connectivity within and between the 3 sub-networks than both the lower frequency 1 Hz and the irregularly spaced PFM protocol.This suggests that tVNS delivered at a stimulation frequency of 24 Hz is the most suitable of those considered in this study for the treatment of mood disorders such as MDD.
In MDD, dysfunction both within and between the subnetworks has been observed [65].The DMN is commonly split into two sub-networks: the anterior sub-network which involves the medial prefrontal cortex, and the posterior subnetworks which centers on the posterior cingulate cortex, precuneus and the angular gyrus [32].Zhu et al. noted both increased and decreased connectivity within the DMN for first episode, treatment naive patients, which appeared to be split into the anterior and posterior parts of the DMN, with the anterior part showing increased connectivity and the posterior part showing decreased connectivity [40].We also identified variation between the anterior and posterior parts of the DMN, with these areas showing opposite effects to the stimulation, which supports the notion that both parts of the DMN act separately in MDD.
The amygdala was also shown to have altered connectivity with various regions in MDD, with Ramsubbu et al. observing reduced connectivity between the amygdala and the prefrontal cortex, insula and middle superior temporal lobe, whilst also observing increased connectivity between the amygdala and the temporal pole, which was negatively correlated with symptom severity and anxiety scores [66].Overall, they noted that the left amygdala showed more reduced connectivity with cortical regions than the right amygdala.Our results showed increased connectivity in response to tVNS delivered at 24 Hz between the left amygdala and frontal and subcortical regions that may counteract the reduced connectivity observed by Ramsubbu et al. and produce an anti-depressive effect.
In a study using on patients with MDD, Fang et al. noted that functional connectivity between the DMN and both the insula and parahippocampus decreased after four weeks of Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.tVNS delivered twice daily at 20 Hz [29].We observed a similar decrease between the DMN and the insula in the sub-delta band of the 24 Hz -PFM contrast; however, we noticed an opposite effect with the parahippocampus showing increased connectivity with many areas of the DMN across numerous frequency bands and contrasts.Fang et al. observed an increase in connectivity between the DMN, precuneus and orbital prefrontal cortex compared with sham stimulation, which we did not observe in our study on healthy participants, although we did observe increased connectivity between all these regions and other areas of the brain.The difference in results may be due to the prolonged period of stimulation that Fang et al. used prior to fMRI after four weeks of daily tVNS compared with the concurrent stimulation and MEG utilized in our study.
In addition to MDD, altered connectivity within and between these sub-networks has been linked to memory conditions such as Alzheimer's Disease and dementia.Grieder et al. noted that functional connectivity between the hippocampus and the posterior cingulate cortex was decreased in Alzheimer's disease and observed a general decrease in connectivity with the DMN [67].Schwab et al. also observed reduced connectivity between the right temporal lobe and the right orbitofrontal cortex in both dementia and Alzheimer's disease [68].
In a study using tVNS on patients with mild cognitive impairment, Williamson et al. identified increased connectivity between multiple nodes of the SN [69].Our results suggest that tVNS at a stimulation frequency of 24 Hz increased connectivity in all of these same areas.Significantly, the right hippocampus has been linked to memory, so altered connectivity from the right hippocampus to other areas that we observed in our results suggests that tVNS has the potential to be used to target these areas which may have a therapeutic effect on both memory and cognition.
In addition to clinical contexts, functional connectivity has also been investigated for it's role as a biomarker to predict individual response to VNS.Mithani et al. used resting state functional connectivity to identify differences between positive responders and non-responders to pediatric VNS for epilepsy [70].Responders were observed to have greater resting state functional connectivity in areas of the limbic and sensorimotor networks compared to non-responders.Further research may uncover whether functional connectivity could be utilized as a screening method to identify those who would respond best to tVNS as a therapy tool.Further research into individual variability may uncover whether functional connectivity could be utilized as a screening method to identify those who would respond best to tVNS as a therapy tool (Figure 9).
In our study, stimulation was delivered for a period of 10 minutes, and MEG was recorded during this time frame.Brain response was not measured either before or after the tVNS was delivered.Our study was focused on whether stimulation and MEG could be performed concurrently.The question of how long tVNS effects last is left for future investigations.
The data was down-sampled from 5kHz to 500 Hz, however since the brain activity of interest occurred between low delta (0.5 Hz) and high gamma bands (125 Hz), this satisfied the Nyquist theorem.
In this study, healthy participants were used.In the future, it is important to investigate whether the changes in connectivity observed here are similar in a cohort of patients with mood disorders.Resting state data was also not collected in this study to limit the participant's time in the scanner, therefore future studies should contrast the participant's response to tVNS with resting state in order to quantify alterations in connectivity which may be compared with the literature to determine the effects of tVNS on MDD and other mood disorders.

V. CONCLUSION
tVNS was observed to alter functional connectivity in healthy participants in a number of brain areas, notably regions linked to the DMN, SN and CEN.Altered resting state connectivity in and between these sub-networks has been observed in MDD and other mood disorders.Altered connectivity was observed in many more regions as a result of higher frequency 24 Hz stimulation when contrasted with 1 Hz stimulation.This study was performed on healthy participants, and thus any clinical implications remain to be elucidated.The results in this manuscript suggest that tVNS delivered at 24 Hz should be investigated further for it's potential as a treatment option for MDD and other mood disorders.

Fig. 1 .
Fig. 1.The positioning of the stimulation electrodes on the left ear.The return electrode was located on the tragus for both the sham site and cymba concha site stimulation.Insets: Schematic of stimulation protocols investigated in this study.Green bars indicate timing of stimulation pulses.PFM: pulse frequency modulation, dashed blue line indicates modulation frequency (6 Hz).

Fig. 2 .
Fig. 2. Duration and timing of stimulation artifacts in the raw data.Blue line indicates the magnetometer sensor readings and black lines indicate timing of the peak of the stimulation pulses.Red lines bound the region to be interpolated.

Fig. 3 .
Fig. 3.Regions involved in three sub-networks implicated in the pathology of depression.

Fig. 4 .
Fig. 4. Inspection of p-values for individual edges.Grey histogram indicates the null-distribution from maximum statistics, green dots indicate non-significant p-values, blue dots indicate p-values considered significant and red line indicates two-tailed significance threshold from permutation test.

Fig. 5 .
Fig. 5.All significant connectivity changes (p<0.05) in the 24 Hz -Sham comparison (both whole-head and ROI combined), *denotes edges with p values < 0.025.Please see TableIfor a summary of all significant edges.

Fig. 6 .
Fig. 6.All significant connectivity changes (p<0.05) in the 24 Hz -1 Hz comparison (both whole-head and ROI combined), *denotes edges with p values < 0.025.Please see TableIfor a summary of all significant edges.

Fig. 9 .
Fig. 9. Boxplot of individual variability for all significant edges (abbreviated, please see Table I for a full description of all edge names).

TABLE I SUMMARY
Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.