2.1.

Participants

Eleven healthy adults (nine males and two females, age range 20 to 49 years) were recruited from the local community. Participants could not have a past or current psychiatric disorder, history of a significant medical disorder, a presently unstable medical condition, caffeinism (i.e., withdrawal symptoms with three days of abstinence), nicotine use, be pregnant or breast feeding, or currently taking any medication. The participants were screened with the Structured Clinical Interview for DSM-IV Axis I Disorders (SCID-I),²⁷ Transcranial Magnetic Stimulation Adult Safety Screen (TASS) form,²⁸ medical history review, and physical exam. A urine sample was obtained for a drug screen and a urine pregnancy test (if the participant was a woman with childbearing potential). Only eligible participants after the screening underwent simultaneous rTMS/fNIRS. The study protocol was approved by the University of Texas Southwestern Medical Center (Dallas, Texas) Institutional Review Board. Written informed consent was obtained from every participant prior to any study evaluations or procedures.

2.2.

Experimental Procedures

The study involved two visits that were two to three days apart to reduce the chance of carryover effects from Visit 1 to Visit 2. In each visit, 1-Hz rTMS were applied over the participant’s motor cortex and prefrontal cortex in two separate sessions using a Neuronetics Model 2100 CRS system (Malvern, Pennsylvania).⁹ A continuous-wave fNIRS system (CW5, TechEn Inc., Massachusetts)²⁹ was used to measure the cortical activities during rTMS. The CW5 system had magnetically compatible fibers and probe holder, which were inspected and confirmed in MRI right before the study.

Visit 1: Participant was positioned in the TMS device chair and adjusted to ensure comfort throughout the experiment. The primary motor cortex was investigated first because this region produced an observable and quantifiable behavioral change, which was used to adjust the stimulation intensity for each participant. An fNIRS probe was attached against the scalp with Velcro straps and elastic bandages. Inside the probe, the tips of the fibers were bent into an $L$ shape so that the fibers could be pressed on the scalp. The fNIRS probe consisted of eight sources and 16 detectors, as shown in Fig. 1(a), which covered the primary motor cortices on both hemispheres. It provided a total of 28 channels (measurements) at a source-detector distance of 3.2 cm, 14 channels on each hemisphere. Once the fNIRS probe was in place, a TMS coil was placed over the left hemisphere and on top of the fNIRS probe, as shown in Fig. 1(b). The coil consisted of two loops of wire in a figure-of-eight arrangement; each loop was approximately 7 cm in diameter. The two loops were left-right oriented along a coronal line and adjusted to have approximately equal distance from the scalp. The location of maximal stimulation of the right abductor pollicis brevis (left motor cortex) was determined using the visual method.³⁰ The participant’s individual RMT was determined using the T.M.S. Motor Threshold Assessment Tool³¹ three times and averaged. The stimulation intensity during the experiment was 120% RMT or 100% of TMS machine output if it could not reach 120% RMT.

Fig. 1

Schematics of simultaneous rTMS/fNIRS measurement on participant’s head: (a) the geometry of the fNIRS probe that consisted of eight sources (filled circles) and 16 detectors (filled squares). The probe provided 28 nearest source-detector pairs (lines) in a distance of 3.2 cm, 14 pairs on each hemisphere; (b) the real placement of the fNIRS probe and TMS coil on a participant’s motor cortex.

Repetitive TMS was performed while the participant sat stably and awake. The stimulation epoch consisted of 10 s of 1-Hz magnetic pulses (i.e., 10 pulses) and 80 s of rest, which was repeated 15 times. The fNIRS data acquisition was initiated 1 min prior to the first epoch and ceased 1 min after the last epoch was completed, which resulted in a total recording time of 24 min and 30 s. After stimulation, the coordinates of coil position on the TMS device chair were recorded that would be used in Visit 2. The position of the fNIRS probe on the participant’s head was determined by referring to the distance from the nasion to the center of the probe in a para-sagittal line.

The participant was given a short break of approximately 10 min after the stimulation session on the primary motor cortex. Then the fNIRS probe and the TMS coil were moved to the prefrontal cortex that was defined as 5 cm anterior to the motor cortex in a para-sagittal line. The stimulation was performed with the same intensity, frequency, and times as it was on the motor cortex, while the brain activity was recorded by fNIRS.

Visit 2: Participants returned for the second visit two to three days after Visit 1. Each of the participants was interviewed by a physician to ensure there were no changes to medical condition and to assess for any adverse effects related to the Visit 1 study procedures. Using the coordinates on the TMS device chair, the participant was placed back in the same position as in the first visit. The fNIRS probe was placed on the scalp at the same position of Visit 1 by referring to the distance from the nasion to the probe center in a para-sagittal line. The location of maximal stimulation of the abductor pollicis brevis was confirmed, and the RMT was determined again using the T.M.S. Motor Threshold Assessment Tool. This was to confirm that no significant change in RMT (greater than 10% change) occurred since Visit 1. The rTMS was performed with same intensity, frequency, and times as it was in Visit 1, on the primary motor cortex and the prefrontal cortex, respectively.

2.3.

Data Screening and Processing

While some of the experimental data were reported in our previous publication,²⁰ there exist two major updates or differences in data analysis between the last study and this paper: (a) data screening for qualified channels and subjects and (b) selection of approximate values of differential path-length factor (DPF) at two wavelengths. Our earlier study focused on channel-wise data analysis.²⁰ For this study, more rigorous selection criteria to screen the data were needed because a topographic image can be formed only when all measurement channels provide high-quality data.

In this study, the fNIRS data from each stimulation session were processed in the following steps:

First, the raw data were inspected for the entire time course to exclude: (a) some channels that had high-frequency instrumental noise larger than 15% of the baseline intensities and (b) some epochs that had significant discontinuities caused by facial muscle twitching.

Second, in order to derive a reliable image of the cortex based on all possibly good-quality channels available, we excluded stimulation sessions that had: (a) bad channels beneath the center of TMS coil as well as the contralateral counterparts; or (b) more than two bad channels on each hemisphere. After such data screening, the qualified subjects to be analyzed in this study were fewer than those in our previous report.²⁰

Third, a low-pass filter with a cut-off frequency at 0.2 Hz was applied to remove the high-frequency physiological noises in the qualified data, such as the arterial and respiratory waves.

Fourth, the channel-wise $\mathrm{\Delta }[{\mathrm{HbO}}_2]$ and $\mathrm{\Delta }[\mathrm{Hb}]$ were calculated based on the Modified Beer-Lambert Law,¹² in which we used $\mathrm{DPF}=6.8$ for laser wavelength at 690 nm and $\mathrm{DPF}=5.8$ for wavelength at 830 nm, according to the published data on the adult heads.³²

At last, the channel-wise $\mathrm{\Delta }[{\mathrm{HbO}}_2]$ and $\mathrm{\Delta }[\mathrm{Hb}]$ were averaged across all of the qualified epochs in the session to get epoch-averaged hemodynamic responses to rTMS.

2.4.

Topographic Imaging

For each stimulation session, hemodynamic images of the cortical region under the fNIRS probe and the TMS coil were generated based on the epoch-averaged responses across all channels. To do so, first the fNIRS measurements were co-registered to a standard brain atlas in order to determine the cortical regions under the probe. Then a rapid 2-D topography approach with a spherical approximation of the brain surface was developed to form the hemodynamic images of the cortex.

Co-registration: A fundamental limitation of fNIRS is that it measures the hemodynamic changes from the scalp without any accurate knowledge of the underlying cortical structures. In other words, fNIRS as a stand-alone technology can provide only functional information without structural information. The ideal solution to this issue is acquiring each participant’s brain MRI with fNIRS optodes attached, which, however, significantly reduces the convenience of fNIRS. As an alternative approach, probabilistic registration has been widely used in this field,²⁶ which registers the fNIRS measurements onto a standard brain template, instead of on each participant’s own brain. In this approach, the optode positions on the participant’s head are measured using a 3-D digitizer along with several cranial landmarks. The cranial landmarks are used to estimate the affine transformation from the real-world coordinate system to a target brain coordinate system.

The probabilistic registration was implemented in this study. To estimate the locations of fNIRS optodes over the motor cortex and the prefrontal cortex, three adult volunteers (two males and one female) were recruited afterwards. The fNIRS probe was placed on their heads in the same way as that above-mentioned during simultaneous rTMS/fNIRS. The positions of optodes and five cranial landmarks (the nasion, inion, left and right earholes, and the vertex) were measured with a 3-D digitizer (PatriotTM, Polhemus, Colchester, Vermont, U.S.A.) and then loaded into a publically available toolbox (NIR-SPM, http://bisp.kaist.ac.kr/NIRS-SPM.html)³³ to elicit the positions of optodes and channels (defined as the midpoint of the source-detector pair) in the standard Montreal Neurological Institute (MNI) space.³⁴ The co-registration procedure was performed for each individual volunteer and then the results were averaged. The brain atlas we had used was ICBM 152 nonlinear, asymmetric template (ICBM 2009a Nonlinear Asymmetric $1\times 1\times 1\mathrm{mm}$ template, http://www.bic.mni.mcgill.ca/ServicesAtlases/ICBM152NLin2009).^35,36 This template covered a broad range of age (18.5 to 43.5 years), and the anatomical structures were segmented using $\text{ANIMAL}+\text{INSECT}$ algorithm.³⁷ Figure 2 shows the cortical surface of the brain extracted from the ICBM 152 template.

Fig. 2

Brain surface for optical topography with two coordinates systems. The brain surface was extracted from the standard ICBM 152 template in stereotaxic coordinates $(x,y,z)$. In optical topography, the stereotaxic coordinates were converted into spherical coordinates $(r,\theta ,\phi )$ based on a hemispherical approximation of the brain surface, where r is radial distance, $\theta$ is azimuthal angle, and $\phi$ is polar angle.

Topography: Given the sparse probe geometry and the single source-detector distance used in this study, the current fNIRS data was expected to have low specificity in depth. Therefore, topographic images of the cortical surface in the stimulated regions were generated by interpolating the epoch-averaged, channel-wise data. With limitations in depth specificity, the topography approach provides us with a fast and efficient approximation to obtain hemodynamic images of the cortical region in response to the combined rTMS/fNIRS. It avoids the complex computation to simulate the light propagation in heterogeneous tissues and to solve the inverse problem that would be encountered in optical tomography.

2.5.

Formation of “EasyTopo”

Brain atlas-based topography has been used in a few published reports.^33,38 In general, the data is interpolated on a stereotaxic convex defined by a cohort of 2-D surfaces, which indeed does not match the brain surface perfectly. Sometimes the data have to be further extrapolated. In this study, we developed a rapid topography for brain imaging (referred to as “EasyTopo”) by implementing angular interpolation in a spherical coordinate system. EasyTopo was derived based on the fact that the cortical layer of the brain is approximately a hemispherical surface. Therefore the stereotaxic coordinates of the brain surface in $(x,y,z)$ were converted into spherical coordinates in $(r,\theta ,\phi )$ (refer to Fig. 2). As schematized in Fig. 3(a) to 3(c), the brain activation in the spherical coordinate system represented a distribution mainly in 2-D $(\theta ,\phi )$ space. Therefore, 2-D angular interpolation was conducted in the $(\theta ,\phi )$ space to form an image of the brain activation, which was then projected back onto the brain surface in original 3-D stereotaxic coordinates.

Fig. 3

Schematics of 2-D optical topography in a spherical coordinate system: (a) the $(\theta ,\phi )$ lattices on the brain surface defined in the spherical coordinate system. As the brain surface is approximately in a hemispherical shape, the cortical activation is a distribution mainly on the $(\theta ,\phi )$ lattices regardless the radius; (b) the topographic image of motor cortex generated through data interpolation in $(\theta ,\phi )$ space. The solid curves in the image outline the ROIs identified by $k$-means clustering; (c) the topographic image projected back from $(\theta ,\phi )$ space to the original stereotaxic coordinate system and overlaid on the brain surface.

It is noted that EasyTopo does not choose a constant radius. Each surface pixel is marked by its own radius, $r_i$, as well as its azimuthal and polar angles (${\theta }_i$ and ${\phi }_i$). However, $r_i$ is not reflected in 2-D $(\theta ,\phi )$ map. Analogy to this procedure is how the earth is drawn in a latitude-longitude map: each local place (i.e., a surface pixel) has its elevation (i.e., radius) as well as its latitude and longitude degrees, but only the latitude and longitude degrees are used for mapping. In other words, we can compress all of the surface pixels of the brain onto a 2-D $(\theta ,\phi )$ map even though their radii or elevations differ. In EasyTopo, the hemodynamic data are interpolated on such a map [see Fig. 3(b)], which contains all of the pixels from the brain surface with respective radius values, not a uniform or equal radius. After all the surface pixels with interpolated values are projected back onto the brain according to their radii, the eventual topographic image [see Fig. 3(c)] is continuous. In addition, we wish to point out a special case: it is possible for a pixel on the gyrus and a pixel on the suci to get the same projected value in 2-D $(\theta ,\phi )$ space if they share the same latitude and longitude degrees (i.e., the suci pixel is beneath the gyrus pixel, and they are located along the same radical direction). In this case, we should ignore the suci pixel since fNIRS is much more sensitive to the shallower tissues.

Compared with the previous topography methods,^33,38 EasyTopo was more computationally efficient and did not require any data interpolation/extrapolation in 3-D stereotaxic space. Another advantage of EasyTopo is that the data between two adjacent channels is interpolated along their included angle in the spherical coordinate system, rather than along a straight line going under the brain surface. The former interpolation is more physiologically meaningful, agrees better with the reality, and serves more accurate in image representation.

2.6.

Identifying Regions of Interest

To characterize the spatial pattern of rTMS-evoked response, a $k$-means clustering algorithm was used to identify regions of interest (ROIs) from the background. This method had been used in our previous study³⁹ to identify positive motor response from the background. In this study, the rTMS-evoked response could be either positive or negative. Therefore the algorithm was modified, and it segmented the topographic image into two to three regions through one of the following cases:

For convenience, the clustering algorithm was applied on the interpolated images in $(\theta ,\phi )$ space rather than the eventual images in stereotaxic coordinates, as schematized in Fig. 3(b). Once the ROIs were identified, the solid angle of each pixel in the ROIs was computed as:

In this study, the $A_{\mathrm{ROIs}}$ was estimated in the hemisphere ipsilateral to the rTMS (i.e., left hemisphere) and in the one contralateral to the rTMS (i.e., right hemisphere) separately. To evaluate the laterality of $A_{\mathrm{ROIs}}$ between two hemispheres, we introduced a lateralization factor, $L$, which was computed as:⁴⁰

2.7.

Assessing Test-Retest Reliability

To comprehensively evaluate the test-retest reliability of the rTMS-evoked hemodynamic response, we used three types of indices that have been introduced in Ref. 25: (a) the intraclass correlation coefficients, ICCs; (b) reproducibility of the identified ROIs in sizes, $R_{\text{size}}$, in two repeated sessions; and (c) degree of overlap, $R_{\text{overlap}}$, between the identified ROIs in two repeated sessions.

The first index, ICC, assesses the reliability of the quantified hemodynamic changes in two repeated sessions based on a one-way random effect model for consistency measurements.^22,24,25,41 In this study, the hemodynamic changes within the identified ROIs of each session were averaged as a measure of rTMS-evoked response in amplitude. Two ICCs were defined to evaluate the reliability of a single measurement, ${\mathrm{ICC}}_{\text{SINGLE}}$, and the reliability of the mean of multiple measurements, ${\mathrm{ICC}}_{\text{AVERAGE}}$:

The second index, $R_{\text{size}}$, assesses the agreement of the identified ROIs in size between two repeated sessions⁴² and was defined as:

The third index, $R_{\text{overlap}}$, assesses the degree of overlap of the identified ROIs between two repeated sessions and was defined as:

For all of the indices described above (i.e., ${\mathrm{ICC}}_{\text{SINGLE}}$, ${\mathrm{ICC}}_{\text{AVERAGE}}$, $R_{\text{size}}$, and ${\mathrm{ROI}}_{\text{overlap}}$), values of $\ge 0.80$ were considered as highly reliable, $\ge 0.60$ as moderately reliable, and $<0.60$ as weakly reliable, according to the criteria given by Charter.⁴³

3.1.

Quantification and Topography of rTMS-Evoked Hemodynamic Response

Among all of the 11 participants, three participants had consistently noisy fNIRS data across all the channels in motor stimulation sessions. It was presumably due to dense and/or dark-colored hair, which could lead to loose contact of fNIRS optodes on the head or high light absorption⁴⁴ if some hair was left between the optode and scalp. Thus these three participants were excluded from data analysis (i.e., $N=8$ for motor stimulation sessions). For the same reason, we excluded one participant when rTMS was on the prefrontal cortex (i.e., $N=10$ for prefrontal stimulation sessions).

For all of the remaining participants with qualified data, we further excluded a few channels with high-frequency instrumental noise above the threshold. Those excluded channels were approximately 5.8% of total channels in motor stimulation sessions and 7.3% of total channels in prefrontal stimulation sessions. Also, the excluded channels were not directly beneath the center of TMS coil or its contralateral counterparts, which likely had the maximum response to the stimulation. Given the diffusion nature of fNIRS, the response signal at every excluded channel was approximately replaced with an averaged signal from its nearest neighboring channels. This procedure ensured the topographic image to be relatively continuous across the region of a bad channel. Note that the location of a channel was defined as the midpoint between its nearest source and detector. The nearest neighboring channels of an excluded channel was those that have the minimal Euclidean distance to the excluded channel.

For the qualified data after screening, previously we have reported a deactivation pattern at the group level, indicted by robust decreases in $\mathrm{\Delta }[{\mathrm{HbO}}_2]$ and slight increases in $\mathrm{\Delta }[\mathrm{Hb}]$, during and shortly after rTMS.²⁰ As an example, Fig. 4(a) and 4(b) shows the grand-averaged temporal responses seen in Visit 1 by selecting two to four channels under the center of TMS coil (i.e., on the ipsilateral side), at the motor cortex and the prefrontal cortex, respectively. At both the motor cortex and the prefrontal cortex, the most significant changes are seen during the first 20 s of the epoch. Therefore, the channel-wise data in every session was averaged across the first 20 s to form the topographic images.

Fig. 4

Grand-averaged temporal brain responses to rTMS (solid curve: mean; shaded region: standard deviation) at: (a) the primary motor cortex ($N=8$, Visit 1); (b) the prefrontal cortex ($N=10$, Visit 1). At each cortical region, the hemodynamic responses (red: $\mathrm{\Delta }[{\mathrm{HbO}}_2]$; blue: $\mathrm{\Delta }[\mathrm{Hb}]$) were attained by selecting two to four channels under the center of TMS coil for each individual participant. The vertical lines (dark green) in each plot indicate the rTMS pulses, which were applied in the first 10s of the epoch.

It is noted that the curves in both Fig. 4(a) and 4(b) look slightly different from those seen in our previous report (Fig. 2 in Ref. 20). Two reasons cause the difference: (a) In this study, eight subjects in the motor stimulations and 10 subjects in the prefrontal stimulation sessions had been selected to form topographic images of respective cortices. The number of qualified subjects in Ref. 20 was 11 for both motor and prefrontal stimulation studies, based on channel-wise analysis. (b) In this study, we used $\mathrm{DPF}=6.8$ for 690 nm and $\mathrm{DPF}=5.8$ for 830 nm in computation of hemoglobin concentrations. These values were estimated based on the previous publication¹² to approximate light scattering effects in the human brain. In Ref. 20, we had used $\mathrm{DPF}=1$ for both 690 and 830 nm, which significantly deviated from the real DPF values in human brain tissues. The current method corrects or improves the accuracy of quantification greatly with respect to the previous calculation.

Figure 5(a) to 5(d) shows the grand-averaged ($N=8$) topographic images at the motor cortex, obtained in Visit 1. In the ipsilateral hemisphere, the deactivation pattern is localized under the TMS coil, which is indicated by distinct $\mathrm{\Delta }[{\mathrm{HbO}}_2]$ decrease around the motor cortex strip. The deactivation pattern in the contralateral hemisphere is more diffused. Figure 6(a) to 6(d) shows the grand-averaged ($N=10$) topographic images at the prefrontal cortex, again from Visit 1. Similarly, the deactivation pattern in the ipsilateral hemisphere is more localized.

Fig. 5

Grand-averaged topographic images ($N=8$, Visit 1) at the motor cortex: (a) the $\mathrm{\Delta }[{\mathrm{HbO}}_2]$ image in the ipsilateral hemisphere; (b) the $\mathrm{\Delta }[\mathrm{Hb}]$ image in the ipsilateral hemisphere; (c) the $\mathrm{\Delta }[{\mathrm{HbO}}_2]$ image in the contralateral hemisphere; (d) the $\mathrm{\Delta }[\mathrm{Hb}]$ image in the contralateral hemisphere. Each image was normalized by its minima [(a) and (c)] or maxima [(b) and (d)]. The location and orientation of TMS coil over the motor cortex is indicated in [(a) and (b)].

Fig. 6

Grand-averaged topographic images ($N=10$, Visit 1) at the prefrontal cortex: (a) the $\mathrm{\Delta }[{\mathrm{HbO}}_2]$ image in the ipsilateral hemisphere; (b) the $\mathrm{\Delta }[\mathrm{Hb}]$ image in the ipsilateral hemisphere; (c) the $\mathrm{\Delta }[{\mathrm{HbO}}_2]$ image in the contralateral hemisphere; (d) the $\mathrm{\Delta }[\mathrm{Hb}]$ image in the contralateral hemisphere. Each image was normalized by its minima [(a) and (c)] or maxima [(b) and (d)]. The location and orientation of TMS coil over the prefrontal cortex is indicated in [(a) and (b)].

3.2.

Assessment of Test-Retest Reliability for rTMS/fNIRS

The test-retest reproducibility was evaluated based on the $\mathrm{\Delta }[{\mathrm{HbO}}_2]$ images only because it demonstrated robust changes associated with rTMS. We neglected the $\mathrm{\Delta }[\mathrm{Hb}]$ images since the $\mathrm{\Delta }[\mathrm{Hb}]$ changes were insignificant at the group level. Most participants had mainly decreased or negative $\mathrm{\Delta }[{\mathrm{HbO}}_2]$ images. Consequently, in $k$-means clustering, the cluster initialized by the minimum value of the image represented the ROIs that were significantly deactivated by rTMS. Interestingly, we did observe one participant in two repeated motor stimulation sessions and another participant in two repeated prefrontal stimulation sessions with consistently increased or positive $\mathrm{\Delta }[{\mathrm{HbO}}_2]$ images. In these special cases, the cluster initialized by the maximum value of the image represented the ROIs. Table 1 summarizes the resultant areas of ROIs in each hemisphere and the laterality between two hemispheres after $k$-mean clustering. At the primary motor cortex, the rTMS evoked a grand-averaged area of 6.1 to $11.5{\mathrm{cm}}^2$ in each hemisphere. At the prefrontal cortex, the rTMS evoked a grand-averaged area of 12.1 to $16.6{\mathrm{cm}}^2$ in each hemisphere, which is significantly larger than that at the primary motor cortex (t-test, $p<0.01$ for the ipsilateral hemisphere and $p<0.02$ for the contralateral hemisphere). At both the primary motor cortex and the prefrontal cortex, the lateralization factor indicates that the ROIs in the ipsilateral hemisphere are slightly larger than those in the contralateral hemisphere.

Table 1

Area of ROIs in the ipsilateral hemisphere, AROIs_ipsi (cm2), and in the contralateral hemisphere, AROIs_contra (cm2), as well as laterality between two hemispheres, L, quantified after k-means clustering.

At the primary motor cortex, high reliability for single measure of $\mathrm{\Delta }[{\mathrm{HbO}}_2]$ (mean value within the identified ROIs) with ${\mathrm{ICC}}_{\text{SINGLE}}=0.86$ and high reliability for average measure of $\mathrm{\Delta }[{\mathrm{HbO}}_2]$ with ${\mathrm{ICC}}_{\text{AVERAGE}}=0.92$ were attained at the group level ($N=8$). However, it is noteworthy that the positive $\mathrm{\Delta }[{\mathrm{HbO}}_2]$ response from one participant greatly enlarged the between-subject mean square (BMS) when calculating ${\mathrm{ICC}}_{\text{SINGLE}}$ and ${\mathrm{ICC}}_{\text{AVERAGE}}$. If we excluded this participant from data analysis, moderate reliability for single measure of $\mathrm{\Delta }[{\mathrm{HbO}}_2]$ with ${\mathrm{ICC}}_{\text{SINGLE}}=0.68$ and high reliability for average measure of $\mathrm{\Delta }[{\mathrm{HbO}}_2]$ with ${\mathrm{ICC}}_{\text{AVERAGE}}=0.81$ were attained ($N=7$). At the individual level, the identified ROIs showed averagely high agreement in size ($R_{\text{size}}=0.81\pm 0.11$) and nearly moderate degree of overlap ($R_{\text{overlap}}=0.56\pm 0.16$). The individual results of mean $\mathrm{\Delta }[{\mathrm{HbO}}_2]$ within the identified ROIs, $R_{\text{size}}$ and $R_{\text{overlap}}$ are summarized in Table 2.

Table 2

Individual results of mean Δ[HbO2] (μM) within the identified ROIs, Rsize and Roverlap quantified in the motor stimulation sessions.

Similarly, at the prefrontal cortex, high reliability for single measure of $\mathrm{\Delta }[{\mathrm{HbO}}_2]$ with ${\mathrm{ICC}}_{\text{SINGLE}}=0.81$ and high reliability for average measure of $\mathrm{\Delta }[{\mathrm{HbO}}_2]$ with ${\mathrm{ICC}}_{\text{AVERAGE}}=0.89$ were attained at the group level ($N=10$). Here it is also noted that the positive $\mathrm{\Delta }[{\mathrm{HbO}}_2]$ response from one participant greatly enlarged BMS. If we excluded this participant from data analysis, moderate reliability for single measure of $\mathrm{\Delta }[{\mathrm{HbO}}_2]$ with ${\mathrm{ICC}}_{\text{SINGLE}}=0.60$ and moderate reliability for average measure of $\mathrm{\Delta }[{\mathrm{HbO}}_2]$ with ${\mathrm{ICC}}_{\text{AVERAGE}}=0.75$ were attained ($N=9$). At the individual level, the identified ROIs showed averagely high agreement in size ($R_{\text{size}}=0.88\pm 0.09$) and moderate degree of overlap ($R_{\text{overlap}}=0.67\pm 0.15$). The individual results of mean $\mathrm{\Delta }[{\mathrm{HbO}}_2]$ within the identified ROIs, $R_{\text{size}}$ and $R_{\text{overlap}}$ are summarized in Table 3.

Table 3

Individual results of mean Δ[HbO2] (μM) within the identified ROIs, Rsize and Roverlap quantified in the prefrontal stimulation sessions.

4.1.

Spatial Characterization of rTMS-Evoked Response

A rapid optical topography on the basis of a standard brain atlas was developed to spatially characterize the rTMS-evoked responses on both ipsilateral and contralateral hemispheres at the primary motor cortex and the prefrontal cortex, respectively. Spatial pattern of the topographic images, represented by ROIs, was identified with a $k$-means clustering algorithm. These methods enabled objective determination of the rTMS-evoked response measured with fNIRS, and quantitative evaluations of the test-retest reliability of this combined technology at the individual level. Interestingly, we also found that the rTMS-evoked responses had greater (nearly two times) hemodynamic changes (see Fig. 4) and larger area of ROIs (see Table 1) at the prefrontal cortex as compared with those at the primary motor cortex. The difference in the areas of ROIs was statistically significant. These differences probably result from the differences in neurophysiology or neuro-anatomy of the different brain regions in response to rTMS. One possible explanation is as follows: the intensity of magnetic pulses is attenuated by the distance from the coil; the optical signals are also attenuated greatly by the scalp and the skull before reaching the brain. Therefore, a shorter scalp-to-brain distance could result in (a) a larger area of magnetic stimulation that generates the electrical potential above the firing threshold of the neurons as well as (b) a higher measurement sensitivity of fNIRS. Unfortunately, the sparsity of our current fNIRS probe did not support reliable 3-D tomographic imaging to verify this explanation; it needs to be further investigated or confirmed in future. Nevertheless the findings that the rTMS-evoked responses had greater [${\mathrm{HbO}}_2$] changes and larger extended ROIs in the prefrontal cortex than in the primary motor cortex are new, providing valuable evidence and observation to neuroscientists and psychiatrists for a better understanding of neurological effects of rTMS on the human brain.

4.2.

Test-Retest Assessment of rTMS/fNIRS

For the second aspect of this paper, namely, the test-retest assessment, it is noteworthy that several measurement conditions in this study are essential: first, the participant was placed in the same position by referring to the coordinates on the TMS device chair; second, the fNIRS probe was placed in the same position by referring to the distance from the nasion to the center of the probe; third, the individual RMT was determined in both visits and confirmed unchanged. The intensity of rTMS was identical between the two visits. The overall lessons learned in test-retest reliability of simultaneous rTMS/fNIRS are as follows: At the primary motor cortex, moderate-to-high test-retest reliability between two repeated sessions could be achieved in quantity of $\mathrm{\Delta }[{\mathrm{HbO}}_2]$, represented by ${\mathrm{ICC}}_{\text{SINGLE}}$ and ${\mathrm{ICC}}_{\text{AVERAGE}}$. The identified ROIs showed overall high reproducibility in size and nearly moderate degree of overlap between two repeated sessions (see Table 2). At the prefrontal cortex, the quantity of $\mathrm{\Delta }[{\mathrm{HbO}}_2]$ showed moderate test-retest reliability between two repeated sessions. The identified ROIs showed overall high reproducibility in size and moderate degree of overlap between two repeated sessions (see Table 3).

Several earlier studies given by Refs. 2223. to 24 focused on the fNIRS reliability with the test-retest interval between three weeks and one year. The purpose of this study was to test the stability of the rTMS/fNIRS measurement method, not to test the stability of the brain. The acute effects of a single session of rTMS are generally thought to last about 20 to 30 min with standard pulse sequences as used in this protocol.⁴⁶ The time interval of two to three days used in this study was chosen for the exact purpose of reducing any chance of a carryover effect while minimizing possible physiologic changes in the brain that could account for a difference. Even if there was a carryover effect from the rTMS, it would be expected to reduce the reliability of the technique and not improve it. The limitation of such a protocol design is that it does not test the variability over a wide range of time intervals; further work is planned in this direction.

Further exploration of the test-retest reliability of combined rTMS/fNIRS can be neurological, technical, or both. Neurologically, any cumulative brain effects across the two visits could certainly affect the test-retest reliability, which had been a concern in long-term daily rTMS treatments.⁴⁶ However, cumulative effects were very unlikely to occur in this study because the stimulation sessions were short (150 s of 1-Hz stimulation per session), and the two visits were two to three days apart. Technically, the quality of fNIRS data can be affected by the individual brain anatomy, spatial geometry, and placement-variability of the fNIRS probes.^22–25 In particular, two factors need to be considered when referring to the results in the present paper: first, although the location of motor cortex was determined based on the maximal stimulation of the right abductor pollicis brevis, the location of prefrontal cortex was defined as 5 cm anterior to the motor cortex in a para-sagittal line. The exact same prefrontal region could not be guaranteed to have been stimulated across all subjects, depending on the cranial size. Future studies should consider using other orientation points, such as markers from the 10 to 20 international EEG system, to determine the location of prefrontal regions more precisely. Second, because of the initial experimental design, the locations of the fNIRS probe on each participant’s head were not measured for co-registration with a standard brain atlas during simultaneous rTMS/fNIRS. To estimate the cortical regions under the fNIRS probe, digitizer readings were obtained from other three volunteers who did not undergo simultaneous rTMS/fNIRS. Their co-registration results were averaged to provide a rough, probabilistic estimation about the cortical location interrogated by the probe in order to implement the brain atlas-based topography. However, the head shape and anatomical structures vary across participants, so nonsubject-specific co-registration may introduce location errors. For future studies, brain topography should be performed according to each participant’s own optode-brain relationship in order to minimize the location-variability.

The practice of replacing excluded channels with an averaged signal from its nearest neighboring channels was an approximate step to avoid a discontinuous or void area in topographic images. Since each channel in the probe array (Fig. 1) corresponded to a certain spatial region in the tomographic image, the averaged signal from nearest neighboring channels (corresponding to other regions) was only an estimation to reflect or replace the missing signal at the excluded channel. Such a practice would introduce spatial errors at the excluded channels for the identified spatial distribution of rTMS-evoked hemodynamic response and ROIs, thus leading to errors in the quantified $R_{\text{size}}$ and $R_{\text{overlap}}$. The above statement implies that the more measurement channels were excluded for data analysis, the less reliable the test-retest results could be. Therefore it is critical and important in the experimental setup to improve optode-scalp contact and thus to reduce the number of excluded channels for data analysis. In this study, while we did have several measurement channels excluded due to highly noisy signals, they were only 5.8% of total channels in motor stimulation sessions and 7.3% of total channels in prefrontal stimulation sessions, also not directly beneath the center of TMS coil or its contralateral counterparts. Therefore the test-retest assessment of combined rTMS/fNIRS on the motor and prefrontal cortices should be reasonably reliable, as reported in Sec. 3.

Finally, the hemodynamic responses presented in this paper are only valid for a specific rTMS paradigm. TMS activates a mixture of neuron populations with different electrical thresholds; its eventual outcomes depend on stimulation intensity, frequency, shape and orientation of the coil, and other factors. For example, several previous studies have suggested that stimulation at a higher frequency than 1 Hz tends to increase cortical excitability, whereas stimulation at a lower frequency ($\le 1\mathrm{Hz}$) tends to decrease cortical excitability.^47,48 In this study, cortical deactivations associated with 1-Hz rTMS were identified among most participants. However, we did observe that one participant during motor stimulation and another participant during prefrontal stimulation had reproducible cortical activations. These special cases reveal and imply the complexity of rTMS-evoked brain responses, especially when the participant pool becomes large. In our future studies, we plan to use other stimulation frequencies besides 1 Hz to examine frequency dependence of the rTMS-evoked hemodynamic responses and thus to gain a better understanding of the brain under rTMS treatment.