A novel concurrent TMS‐fMRI method to reveal propagation patterns of prefrontal magnetic brain stimulation

Abstract Major depressive disorder (MDD) is a severe mental disorder associated with high morbidity and mortality rates, which remains difficult to treat, as both resistance and recurrence rates are high. Repetitive transcranial magnetic stimulation (TMS) of the left dorsolateral prefrontal cortex (DLPFC) provides a safe and effective treatment for selected patients with treatment‐resistant MDD. Little is known about the mechanisms of action of TMS provided to the left DLPFC in MDD and we can currently not predict who will respond to this type of treatment, precluding effective patient selection. In order to shed some light on the mechanism of action, we applied single pulse TMS to the left DLPFC in 10 healthy participants using a unique TMS‐fMRI set‐up, in which we could record the direct effects of TMS. Stimulation of the DLPFC triggered activity in a number of connected brain regions, including the subgenual anterior cingulate cortex (sgACC) in four out of nine participants. The sgACC is of particular interest, because normalization of activity in this region has been associated with relief of depressive symptoms in MDD patients. This is the first direct evidence that TMS pulses delivered to the DLPFC can propagate to the sgACC. The propagation of TMS‐induced activity from the DLPFC to sgACC may be an accurate biomarker for rTMS efficacy. Further research is required to determine whether this method can contribute to the selection of patients with treatment resistant MDD who will respond to rTMS treatment.

Unlike ECT, rTMS is a targeted noninvasive brain stimulation method with only mild side effects and has proven to be effective in the treatment of TR-MDD (George, Lisanby, Avery, Mcdonald, & Durkalski, 2014). TMS is a means of using electromagnetic induction in order to stimulate a brain region. Repetitive delivery of TMS pulses (rTMS) to a brain region modulates the excitability of the stimulated area, inducing changes in neural plasticity (Allen, Pasley, Duong, & Freeman, 2007). These neuroplastic changes outlast the duration of stimulation and are believed to be induced through long-term potentiation/depression mechanisms (Esser et al., 2006;Ishikawa et al., 2007). High (>5 Hz) or low (<5 Hz) frequency stimulation results in a lasting increase or decrease in excitability, respectively (Esser et al., 2006;Valero-Cabré, Payne, & Pascual-Leone, 2007). Stimulation is applied to a focal region in the brain, but the effects of TMS are not limited to the stimulated brain region, but spread to other cortical areas (Bestmann, Baudewig, Siebner, Rothwell, & Frahm, 2005;De Weijer et al., 2014;Rahman et al., 2013). Repetitive stimulation has been shown to induce a clinically relevant effect in MDD and has recently obtained FDA approval for its application in MDD (George et al., 2014;O'Reardon et al., 2007). Patients with MDD are treated through high frequency stimulation of the dorsolateral prefrontal cortex (DLPFC), which causes an antidepressant effect in around 25% of patients with TR-MDD (O'Reardon et al., 2007).
Although already applied clinically, little is known about the mechanism of action of high frequency stimulation of the DLPFC. Over the last decades neuroimaging studies have investigated changes in the brain related to MDD intensively and several neuroanatomical regions have been found to exhibit abnormal activity in patients with MDD, with the subgenual anterior cingulate cortex (sgACC) attracting most attention. Neuroimaging studies have shown that baseline metabolic activity in the sgACC is increased in patients with MDD and that normalization of the activity in the sgACC correlates with relief of depressive symptoms (Kennedy et al., 2001;Mayberg et al., 2000Mayberg et al., , 2005Videbech, 2000). It is hypothesized that rTMS of the DLPFC induces an antidepressant effect through direct or indirect neuromodulation of the abnormal activity in the sgACC (Fox, Buckner, White, Greicius, & Pascual-Leone, 2012). MRI functional connectivity studies show that treatment outcome positively correlates with functional connectivity strength between the sgACC and the DLPFC, providing some evidence for this hypothesis (Baeken et al., 2014). However, the evidence is limited as only a small number of patients received rTMS treatment. Furthermore, Fox and Baeken assume that a functional connection based on resting state fMRI data provides the pathway for TMS effects evoked in the DLPFC, while there is no evidence that shows that TMSinduced activity can actually propagate from the DLPFC to the sgACC.
The antidepressant effect of rTMS is restricted to a quarter of the patients with TR-MDD and is known to vary substantially between patients (O'Reardon et al., 2007). This stresses the need for a better understanding of the effects of rTMS treatment of the DLPFC in order to identify the patients who will respond to this type of rTMS treatment. Identification of the individual propagation pattern of TMSevoked activity in response to stimulation of the left DLPFC can reveal the mechanisms of action of rTMS treatment and provide clues for further improvement of such treatment using individualized treatment methods. This can be done by using a patient's individual response to TMS to optimize the effects after subsequent treatment with repetitive stimulation. More specifically, identification of the propagation patterns can reveal whether TMS-induced activity evoked at the DLPFC has the ability to propagate to the sgACC, and potentially modulate the activity in the sgACC (Fox et al., 2012). Therefore, we investigated the propagation pattern of TMSinduced activity after stimulation of the left DLPFC using single pulse TMS. Because TMS effects are strongly affected by TMS coil placement (with respect to individual brain morphology), we also investigated the effect of TMS coil placement on propagation patterns of TMS-induced activity (Fitzgerald et al., 2009;Sack et al., 2009). In order to achieve these goals, we applied single pulses of TMS to the left DLPFC during a functional MRI recording in 10 healthy participants, using a novel concurrent TMS-fMRI technique (De Weijer et al., 2014).

| MATERIALS AND METHODS
The experimental procedure was approved by the medical ethical committee of the University Medical Center Utrecht (UMCU), Utrecht, The Netherlands. All participants provided written informed consent and were screened for MRI and TMS exclusion criteria. MRI data were acquired in 10 right-handed participants (Table 1). One participant had to be excluded due to unavailability during the follow-up session. During the experimental procedure, we strictly adhered to the guidelines and recommendations for TMS endorsed by the International Federation for Clinical Neurophysiology (Rossi et al., 2009).
Our concurrent TMS-fMRI setup has previously been used to successfully detect TMS-induced activity in the motor network in response to TMS pulses delivered to the primary motor cortex (M1) of healthy participants (De Weijer et al., 2014). We used this setup to stimulate M1 again, but primarily investigated TMS-induced activity in response to stimulation of the left DLPFC of healthy participants. Both experiments were done using a biphasic Magstim Rapid 2 magnetic stimulator, an MR-compatible TMS coil and a custom designed TMS filter box (all manufactured by Magstim Inc., Whitland, The United Kingdom, www.magstim.com;De Weijer et al., 2014;Mandija, Petrov, Neggers, Luijten, & Berg, 2016). All MR sequences (MRI-only and concurrent TMS-MRI) were performed in a 3T MR scanner (Achieva, Philips Healthcare, Best, The Netherlands).

| Data acquisition
The experiment was divided into two parts: an intake session (MRIonly) and a TMS session (concurrent TMS-MRI) (Figure 1). These sessions are described below.

| TMS session
For each participant, the T1 weighted image acquired during the intake session was segmented with SPM12 to obtain skin, skull, cerebrospinal fluid (CSF), white matter, and gray matter (GM) masks (Penny, Friston, Ashburner, Kiebel, & Nichols, 2011). The segmentations were used to visualize the 3D brain and skin surface in the Neural Navigator (Brain Science Tools, The Netherlands, www. neuralnavigator.com). The location of M1 was obtained from the statistical map acquired during the intake session and marked on the 3D   Next, the TMS coil was placed over the left M1 guided by neuronavigation to determine the resting motor threshold (RMT). This was done by applying single pulses of TMS (with an inter-stimulus interval of 7 s) to the primary motor cortex while increasing the TMS stimulator output until a response in the APB muscle was visible in 5 out of 10 TMS pulses (Schutter & van Honk, 2006).
The concurrent TMS-MRI session was divided in two parts: In the first part TMS was applied to M1 and in the other part TMS was applied to the DLPFC. For the concurrent TMS-MRI experiments, a custom made setup was used. The head was positioned in a custom designed setup between two MR receive coils ( Figure 2c). The TMS coil was attached to a custom made mount which was positioned over the participant's cranium, eliminating potential TMS coil movement.
Additionally, to minimize Lorentz forces on the TMS coil wings, the angle between the TMS coil plane and the MRI static magnetic field was limited to 25 . The head was tilted backward and rotated slightly to match the coil position with the markings on the bathing cap. The head and neck of the participant were supported to increase comfort and to minimize head movement during the scan.
After TMS coil positioning, two sequences were acquired. First, a T2-weighted scan with a TR/TE of 13,609.0/80.0 ms, flip angle of 90 , voxel size of 2 × 2 × 2 mm 3 , scan duration of 3.6 min. The purpose of this scan was to visualize the coil location and orientation with respect to the head. This was done by attaching six custom made markers (small capsules filled with water) to the back of the TMS coil ( Figure 2b), which appear hyper intense on the T2-weighted scan ( Figure 1). Second, a single-shot EPI sequence was acquired with 500 dynamics, a TR/TE of 2,000.0/23.0 ms, flip angle of 70 , FOV of 256 × 119.6 × 208 mm 3 , matrix of 64 × 63, voxel size of 4 × 4 × 4 mm 3 , scan duration of 17 min, 30 slices with a slice thickness of 3.6 mm, and a slice gap of 0.4 mm. We acquired 500 dynamics to make sure that we had sufficient power to detect the effects of single pulses of TMS. During the EPI sequence, single pulses of TMS with an intensity of 115% RMT were interleaved with pulses with an intensity of 60% RMT. TMS pulses were delivered with a random interval of 5-8 dynamics (10-16 s) to avoid habituation.

| Data analysis
The data obtained for stimulation of M1 and the DLPFC were analyzed in the same way. Analysis of the structural and fMRI data was performed with custom scripts and SPM12 (Penny et al., 2011) in the MATLAB R2014a environment (MathWorks Inc., Natik, MA, The United States).
The T1-weighted image was segmented with SPM through unified segmentation with six tissue type priors to obtain a gray matter, white matter, and CSF mask (Ashburner & Friston, 2005).
All EPI volumes were inspected to determine image quality and to identify the presence of potential artifacts. This revealed small random These artifacts are discussed in more detail in the section "Image artifact." All EPI volumes were realigned and normalized using the nonlinear normalization parameters obtained from segmentation of the T1-weighted scan. The EPI volumes were subsequently resliced at a resolution of 4 × 4 × 4 mm 3 and smoothed with a FWHM of 8 mm.
The location of the TMS coil isocenter was reconstructed with respect to the brain based on the location of the fluid markers on the TMS coil. In order to reconstruct the TMS coil isocenter with respect to the brain, the T2-weighted scan was first co-registered to the anatomical T1-weighted scan using SPM12. Then, based on the location of the TMS coil markers in the T2-weighted scan, we were able to reconstruct the TMS coil position and isocenter, as described in De Weijer et al. (2014). Thereafter, the EPI volumes were realigned using SPM12 and the mean EPI scan was co-registered to the T1-weighted scan. The inverse of the EPI to T1 affine transformation and the inverse of the EPI realignment affine transformations were used to create a head movement-corrected reconstruction of the location of the TMS coil isocenter. Thereafter, the center of gravity (COG) of the We performed subject-level event-related GLM analyses in SPM12. A group-level analysis was not performed due to the limited sample size (Desmond & Glover, 2002;Thirion et al., 2007). The generalized linear model (GLM) included two events: single pulses of 115% RMT and 60% RMT. The BOLD response was modeled with the canonical HRF and its first-order derivative. The mechanism through which TMS induces brain activity is different from conventional MRI tasks which investigate voluntary brain activity. Inclusion of the first-order derivative allows more variability in the hemodynamic response, which allows more accurate modeling of the BOLD response during TMS pulse delivery. Two nuisance regressors were included in the analysis: the average BOLD signal in the white matter and the CSF. BOLD signals were filtered with a high pass filter of 80 Hz before construction of the GLM. Statistical images were constructed based on the contrast between TMS pulses of 115% RMT and TMS pulses of 60% RMT. An F-statistic was used to test for variance explained by either the canonical HRF or the first-order derivative as the sum of both weighted regressors in a first-order Taylor expansion, using a threshold at p < .05, FWE corrected (Penny et al., 2011). TMS pulses of 115% RMT were contrasted with TMS pulses of 60% RMT to minimize the neural response to the sensation or sound that is accompanied by TMS pulse delivery.

| Voluntary versus TMS-induced motor activity (M1 session)
We compared motor network activity in response to voluntary hand movements with TMS-induced activity after stimulation of M1. The full findings on TMS-induced activity are summarized in Table 3   the right hemisphere of the cerebellum and the bilateral putamen and thalamus (p < .05, FWE corrected; Figure 3b). Because the participants were instructed to move their thumb based on auditory cues, activity can also be observed in the primary auditory cortex (A1).
During the concurrent TMS-fMRI session, the TMS coil isocenter was located slightly medio-anterior to the maximum BOLD response in the precentral gyrus, with limited displacement of the TMS coil isocenter due to head movement during the session (Figure 3a). TMS-induced activity was defined as the contrast between TMS pulses of 115% RMT and 60% RMT. The participant reported thumb movement in response to high intensity TMS pulses throughout the majority of the session and reported no thumb movements for low intensity TMS pulses. TMS-induced activity can be observed in the bilateral M1 and thalamus, the left SMA and putamen and the right hemisphere of the cerebellum (p < .05, FWE corrected; Figure 3c). Because TMS pulse delivery is accompanied by an auditory "click," activity can also be observed in A1. Additionally, TMS-induced activity can be observed in the bilateral primary visual cortices.

| TMS target (DLPFC session)
The locations of the COG of the TMS coil isocenters were located well within the anatomical landmarks of the DLPFC in all participants ( Figure 4a). The TMS coil isocenter remained within the DLPFC for the majority of the session in all participants, despite small head movements.
Head movement resulted in a maximum displacement of the TMS coil isocenter from the COG of 2.1-6.1 mm (mean: 4 mm) depending on the participant (Table 1), which shows that TMS coil placement was accurate and head movement was limited. In the majority of the participants, the

| TMS-induced activity (DLPFC + M1 session)
An overview of neuroanatomical regions that show TMS-induced activity in response to TMS pulses delivered to the left DLPFC (Table 2) and the left M1 (Table 3) (Fox et al., 2012). That work showed that the strength of the functional connection between the DLPFC and the sgACC correlated positively with treatment outcome after rTMS, implying that rTMS utilizes this functional connection (Baeken et al., 2014). However, evidence for a direct structural connection between these regions was lacking.
Our observations provide the first direct evidence that TMS-induced activity can propagate to the sgACC (at least in some individuals).

| Variability in propagation patterns
Although TMS coil placement was well-controlled and head movement was limited, we observe that TMS-induced brain activity is variable between participants ( Figure 5). More specifically, we observe that stimulation of the same neuroanatomical region results in substantially different propagation patterns of TMS-evoked activity. An explanation of the variability in the propagation patterns are differences in structural brain connectivity, i.e. differences in propagation pathways, between participants (Bürgel et al., 2006;Thiebaut de Schotten et al., 2011). Another explanation comes from differences in (functional) neuroanatomy (Amunts et al., 1999;Der Malsburg, Phillips, & Singer, 2010;Watson et al., 1993). Stimulation of the same neuroanatomical region does not necessarily mean stimulation of the same functional region.
Stimulation of different functional regions can lead to differences in propagation pathways, depending on individual functional connectivity.
Finally, functional brain connectivity shows substantial variability between participants, with strong variability observed in the prefrontal cortex (Mueller et al., 2013). Moreover, resting state functional connectivity studies show that brain functional connectivity shows substantial temporal variability, especially in the sgACC (among other regions) (Allen et al., 2014;Handwerker, Roopchansingh, Gonzalez-Castillo, & Bandettini, 2012;Mueller et al., 2013). To summarize, the structural and functional organization of the brain together with the dynamic nature of functional connectivity could explain the variability in propagation patterns of TMS-induced activity.
Finally, the complex interaction of the TMS E-field with the neuronal populations in the cortex has been shown to be brain-state dependent and therefore not consistent over time (Romei et al., 2008;Sauseng, Klimesch, Gerloff, & Hummel, 2009). EEG recordings have revealed that the pre-existing neuronal oscillatory activity in the TMS target region affects properties of TMS-evoked activity. For example, the amplitude of TMS-induced motor-evoked potentials depends on pre-existing oscillatory activity of the primary motor cortex while phosphene thresholds depend on pre-existing activity in the visual cortex (Romei et al., 2008;Sauseng et al., 2009). Since local neuronal oscillatory activity shows temporal and spatial variability, the TMS effects will also show temporal and spatial variability, which is likely to contribute to the observed variability of propagation patterns of TMSinduced activity (Allen et al., 2014;Arieli, Sterkin, Grinvald, & Aertsen, 1996;Fox, Corbetta, Snyder, Vincent, & Raichle, 2006).
Our findings are in agreement with prior literature, which indicates that propagation of evoked activity is a complex process that, similar to functional connectivity, varies significantly between individuals (Fox et al., 2006;Handwerker et al., 2012;Mueller et al., 2013;Sauseng et al., 2009). Consequently, the ability of rTMS to modulate specific brain regions, and thus inducing an antidepressant effect, is likely to depend on the state of the individual brain network. Concurrent TMS-fMRI allows identification of the present individual structural and functional brain organization and connectivity, related to stimulation of a specific area of interest, such as the DLPFC. This concept can potentially be used to predict treatment outcome or to increase treatment efficacy of repetitive stimulation of the DLPFC in MDD. For example, by targeting treatment at the functional region near the DLPFC that leads to individual sgACC activation, assuming the pathway from prefrontal cortex to sgACC is the mechanism of action of the antidepressant effect of rTMS, as has been suggested (Baeken et al., 2014;Fox et al., 2012). In this way, future studies can investigate whether propagation of TMS-induced activity to sgACC is an accurate predictor of a clinical response to rTMS treatment of the DLPFC in a clinical population of patients with MDD.
We discussed possible (functional) neuroanatomical and state dependent explanations for the variability in TMS-induced activity between participants. However, the observed variability might originate from other sources. Variability can arise from the physiological interaction between the MRI static magnetic field and the electrophysiology of the individual brain or from the physical interaction between the TMS coil position with respect to the MRI static magnetic field.
A possible explanation is that in some cases the TMS-induced currents depolarize the descending white matter tracts, rather than the cell bodies, bypassing synaptic transmission in the TMS target area. As synaptic transmission makes up the majority of the hemodynamic signal measured using fMRI, concurrent TMS-fMRI is predominantly sensitive to synaptic transmission evoked by TMS (Tagamets & Horwitz, 2001).

| Image artifact
The application of concurrent TMS-fMRI is challenged by numerous technical difficulties, a few of which have already been addressed in other works (Bestmann et al., 2004;Ruff et al., 2008). A technical issue which was not previously described is that we observed short deflections (one sample) in baseline activity in a single slice in the vicinity of the TMS coil in EPI volumes during inspection of the BOLD signal (section "Data analysis"). These artifacts were only observed during sessions in which TMS intensity was changed during the experiment and were absent when the machine output was kept constant.
This suggests that currents leaked in the TMS coil, creating a local magnetic field around the TMS coil. This local magnetic field perturbed the MRI static magnetic field during image acquisition, resulting in image distortions. Unfortunately, the cause of the artifact could not be identified with complete certainty.
Fortunately, the effect on the results is negligible, as only a few slices were affected per participant (in <0.5% of all acquired slices per participant) and the artifact in these slices could effectively be removed through interpolation. However, TMS-induced activity in the vicinity of the TMS coil should be interpreted with caution.

| CONCLUSIONS
In conclusion, TMS pulses delivered to the left DLPFC induce activity in a number of connected brain regions, including the sgACC in some participants (4 out of 9). This indicates that the effects of stimulation of the left DLPFC stimulation have the ability to propagate to the sgACC, depending on individual structural connectivity, providing a potential mechanism for the antidepressant effect of rTMS delivered to the left DLPFC. This individual propensity could potentially be used as a predictor of the response to rTMS treatment in patients with MDD. Additionally, the propagation patterns of TMS-evoked activity show substantial variability between participants while TMS coil placement was well-controlled during image acquisition. Combining our observations with prior literature, this implies that the propagation pattern of TMS-induced activity, and thus the ability of rTMS to modulate specific brain areas, could depend on the current state of the individual brain network.

ACKNOWLEDGMENTS
This work was supported by the DeNeCor project as part of the ENIAC Joint Undertaking, (ENIAC131003; http://eniac.eu). The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript. The authors would like to acknowledge Sanne Schuite-Koops for her support in the organization of this study.