2.1.

Protocols and Instrumentation

Concurrent NIRS and fMRI resting-state studies were conducted in six healthy volunteers participating in an ongoing protocol (three men and three women, $\text{average age}\pm \mathrm{SD}$, $30\pm 7\text{years}$). Participants were asked to lie quietly in the scanner and view a gray screen with a fixation point in the center. The Institutional Review Board at McLean Hospital approved the protocol, and volunteers were compensated for their participation.

All MR data were acquired on a Siemens TIM Trio 3T scanner (Siemens Medical Systems, Malvern, PA) using a 32-channel phased-array head matrix coil (200 time points, $\mathrm{TR}/\mathrm{TE}=3000/30\mathrm{ms}$, flip angle 90 deg, $\text{matrix}=64$ by 64 on a 224- by 224-mm FOV, 50 2.5-mm slices with 0.625-mm gap parallel to the AC-PC line extending down from the top of the brain). Physiological waveforms (pulse oximetry and respiratory depth) were recorded using the scanner’s built-in wireless fingertip pulse oximeter and respiratory belt.

The optical device used in this experiment is a near-infrared tissue imager (Imagent, ISS, Champaign, IL) comprising 16 laser sources at 690 and 830 nm (average power of about 1 mW) and four optical detectors (photomultiplier tubes, Hamamatsu Photonics R928). The imager was placed in the control room, and 10-m-long optical fibers (illumination fibers with 400 μm in core diameter) and fiber bundles (a collection fibers with 3 mm in diameter) were used to connect the laser diodes and detectors with optical probes in the 3T scanner. Two (and sometimes three) NIRS optical probes were made from MRI-compatible materials (ABSplus and neoprene) using a three-dimensional (3-D) printer (uPrint SE 3D, Stratasys, Eden Prairie, MN). Each probe included one collection fiber and one pair of illumination fibers (690 and 830 nm). The separation between the collection and illumination fibers was 1.5 cm. The probes were placed over the tip of each participant’s left middle finger and left big toe (and also the right big toe in three participants), as shown in Fig. 1(a) and 1(b), respectively, with a black elastic band. The sampling rate of the NIRS data acquisition was 12.5 Hz to fully sample the cardiac waveform. fMRI data were collected for 10 min; NIRS data were recorded continuously during this time, and for several minutes before and after the fMRI acquisition.

Fig. 1

Experimental setup and example of near-infrared spectroscopy (NIRS) temporal traces. Schematic drawing of the placement of the NIRS probes on the middle fingertip (a) and the big toe (b). Raw temporal traces of changes in total hemoglobin concentration Δ([tHb]) obtained by NIRS at the fingertip and left toe (c). Enlarged section of (c) (indicated by gray block) shown in (d). The low-frequency oscillation (LFO) signal of the finger and toe from (c) permit calculating that the signal measured at the toe is 2.72 s later than that of the finger (e).

2.2.

NIRS Regressor Preparation

Each pair of raw NIRS time courses (690- and 830-nm data) was converted into three time courses representing Δ[HbO], Δ[Hb] and total hemoglobin concentration change (Δ[tHb]) using the modified Beer-Lambert law^32,33 in Matlab (The Mathworks, Natick, MA). For the differential pathlength factors, we used published values of 6.51 at 690 nm and 5.86 at 830 nm.³⁴ The low-frequency component ($0.01\sim 0.15\mathrm{Hz}$) of Δ[HbO], Δ[Hb], and Δ[tHb] was obtained by using a zero-delay Fourier-domain bandpass filter in Matlab. Although BOLD is sensitive to the Δ[Hb] signal, it has been shown to be sensitive to changes in cerebral blood volume (CBV) and cerebral blood flow (CBF) as well.³⁵ In particular, changes in CBV are more closely associated with physiological LFOs;⁹ CBV is best reflected by the Δ[tHb] NIRS waveform.³⁶ Moreover, the Δ[tHb] obtained by NIRS has a higher signal-to-noise ratio than Δ[Hb] and captures BOLD variations in the low-frequency range with greater statistical significance than the Δ[Hb] time course. Therefore, the temporal trace of Δ[tHb] from NIRS was selected for correlation with BOLD signals.

2.3.

Correlations

For each participant, the correlation between fingertip and toe data was calculated. The LFO portion of Δ[tHb] at the fingertip and toes were first extracted using a zero-delay Fourier-domain bandpass filter ($0.01\mathrm{Hz}<\mathrm{LFO}<0.15\mathrm{Hz}$) in Matlab. The two derived LFO signals were then cross-correlated in Matlab to calculate the maximum correlation and corresponding temporal shifts between these two signals. For the concurrent analyses, the fingertip LFOs were used because the blood vessels that supply blood to the arm and the brain are much closer on the vascular tree than those supplying the toes. Moreover, the fingertip is a commonly used blood oxygenation and pulse recording site.

We then applied RIPTiDe^30,31 to both the fingertip NIRS and fMRI datasets for each participant. This analysis was carried out using a combination of in-house Matlab tools and FEAT, part of the FSL analysis package (FMRIB Expert Analysis Tool, v5.98, http://www.fmrib.ox.ac.uk/fsl, Oxford University, UK).³⁷ First, the LFO component of the fingertip Δ[tHb] NIRS waveform was extracted using the same bandpass filter. Because the LFOs detected at the finger are presumed to travel with the blood and arrive in each voxel with variable latency, we applied a range of time shifts to the NIRS signal prior to filtering to capture all the voxels that correlated with the NIRS signal at any latency. Then the temporal signal was downsampled to the 3-s TR of the fMRI acquisition. The resulting time traces were used as regressors in a GLM analysis of the fMRI data. A positive time-shift value corresponded to events that occurred before the time the NIRS data was recorded at the fingertip, meaning the signal arrived earlier in the brain than in the fingertip. We considered 181 time shift values covering the range of $-14.4$ to $+14.4\mathrm{s}$ with respect to the unshifted signal, in steps of 0.16 s (for Video 1, a longer covering range was chosen for the purpose of dynamic display). The range of the time shifts was determined on the basis of a preliminary analysis to find the range over which significant correlations were observed.

Video 1

The dynamic passage of the LFO wave is shown in red-yellow. The green dot in the middle panel indicates the time lags NIRS regressor (from fingertip) was shifted. The total time shift is from $-16\mathrm{s}$ to $+16\mathrm{s}$ at steps of 0.16 s. The video was played repeatedly, each time with shifted coronal and axial views marked by the cross. The color bar indicates the colors corresponding to the $z$ value from 2.3 to 8 (QuickTime, 1.9 MB) [URL: http://dx.doi.org/10.1117/1.JBO.17.10.106004.1].

For each participant, standard fMRI preprocessing steps [including motion correction, high-pass filter ($>0.01\mathrm{Hz}$), and slice timing correction] were applied to the original BOLD data before further analysis. The RIPTiDe procedure determined two independent quantities: the time delay value between the NIRS and fMRI data, and the contribution of the delayed NIRS signal to the BOLD variance. The time delay was determined by fitting a Gaussian to the cross-correlation calculated in the time domain³⁸ (which is more robust for determining subTR time shifts in the presence of aliasing than the direct frequency domain method). This cross-correlation was implemented using 181 invocations of film_gls (part of FSL) with differently time-shifted NIRS Δ[tHb] time courses used as a regressor. The resulting thresholded $z$-statistic maps of the LFO signals were concatenated over all time lags in sequence, and a voxelwise Gaussian fit was performed of the peak of signal in the time lag dimension to obtain the voxelwise delay values. The $z$-statistic of the correlation between the NIRS data and the BOLD at this optimum time lag was then calculated and used to form a new 3-D thresholded $z$-statistic map (here referred to as a “max $z$-statistic map”). The max $z$-statistic map represents all the voxels that highly correlate with the NIRS Δ[tHb] timecourse (fingertip) at any time.

2.4.

Spectral Analyses

To understand the spectral content of the signals being studied, we performed spectral analysis on the peripheral NIRS data and coherence analysis between peripheral NIRS and cerebral BOLD. For each participant, we calculated the power spectral density (PSD) of the Δ[tHb] signal collected at the finger by NIRS using Matlab to understand the spectral components of NIRS data. Furthermore, coherence analyses for each participant between the NIRS data (Δ[tHb] from the finger) and the BOLD from all voxels were performed using the Matlab (function “mscohere”), to understand how well the peripheral LFO (from NIRS) corresponds to BOLD fMRI at each frequency. The magnitude squared coherence (mscohere) is a function of frequency with values between 0 and 1 that indicates how well one signal corresponds to the other at each frequency. The acquisition rate of 0.33 Hz in BOLD ($\mathrm{TR}=3\mathrm{s}$) limits the maximum frequency of the result to be 0.16 Hz. To have meaningful coherence analysis, the BOLD signals were preprocessed using the standard “prestats” analysis in FSL, which includes motion correction, slice time correction, spatial filtering (5 mm), and high-pass filtering ($>0.01\mathrm{Hz}$). The NIRS data (Δ[tHb]) was also preprocessed using spike removal³⁹ (Matlab), bandpass filtering ($0.01\sim 0.15\mathrm{Hz}$) and downsampling to 0.33 Hz. For each voxel, the magnitude squared coherence was calculated from 0 to 0.16 Hz. The coherences from all voxels were averaged. The final steps involved (1) finding the frequency of maximum averaged coherence in each participant and (2) generating the corresponding coherence map at that frequency. Each coherence map was then transformed into standard space for the purpose of comparison.

3.1.

Power Spectral Density of NIRS Data

Figure 2 shows the PSD of the NIRS data collected at the fingertip. The figure confirms the presence of the cardiac pulsation in all the participants and clearly shows the respiration signal in three participants (3, 4, and 6). Participants 1, 2, and 5 had no clear signal at the respiratory frequency. This observation is consistent with our previous work,³¹ in which respiration signals were not always observed by NIRS, while the cardiac pulsation signals were always present. Moreover, the graph indicates the presence of high energy in the LFO band in the NIRS data (as in BOLD). This shows the importance the LFO in both signals, since it is the dominant component. Last, we show that by choosing the low-frequency band ($0.01\sim 0.15\mathrm{Hz}$), we not only preserve the high energy in the LFO, but also reduce the contamination of LFO with other physiological noise (respiration). This is critical in the study of the LFO, which may be independent of other physiological processes.

Fig. 2

Power spectral densities (PSDs) of the temporal traces of Δ[tHb] obtained by NIRS at the fingertip from six participants. Vertical line in each graph at 0.15 Hz indicates the upper limit of what is considered the LFO band.

3.2.

Temporal Correlations between LFOs from the Peripheral Vasculature

Figure 1(c) shows the temporal trace of Δ[tHb] collected by NIRS at the left middle fingertip and on the left big toe of one representative participant. Both signals show visible LFO and cardiac signals [visible in Fig. 1(d), an enlarged part of Fig. 1(c)]. Figure 1(c) also shows the LFOs of these two signals. The maximum correlation between the signals was 0.80 at a time shift of 2.72 s (the data from the toe lags by 2.72 s) which is depicted more clearly in Fig. 1(e), in which the bandpass filter (0.01–0.15 Hz) was applied to both signals to present only the LFO signal. Our previous studies³¹ along with others using laser Doppler^24,25,40 suggested that these LFOs travel with the blood, and thus are closely associated with blood circulation. Consequently, the LFOs will, in general, arrive at different sites (finger and toe) at different times reflecting the pathlength, vessel size, and flow rate of the blood to the target site. The longer pathlength from the heart to the foot (relative to the hand) would explain the time delays between the LFOs detected at the finger and toe. Table 1 summarizes the maximum cross-correlations and time delays between the LFOs measured at different peripheral sites with NIRS. The LFO signal reaches the toe 2.16–4 s later than the finger (time delays: $T_{\text{delay}}=3.07\pm 0.81\mathrm{s}$). For three participants, NIRS data was also collected at the right big toe; the LFOs from the two toes had maximal correlations ($r_{\max }=0.85\pm 0.09$) with small time shifts between sides ($T_{\text{delay}}=-0.02\pm 0.57\mathrm{s}$). Relative to the results obtained comparing finger to toes ($r_{\max }=0.70\pm 0.08$), the correlation is higher between the data from the two toes and exhibits less temporal shift. We also calculated the correlation between the time shifts (fingertip and toe) and the participants’ heights (data shown in Table 1). No linear relationship was found, which suggests that the time delays may be affected not only by height, but also by other conditions, such as age, sex, and general health.

Table 1

A summary of maximum cross-correlations and time delays between the LFOs of left finger and left toe, or left toe and right toe, as well as the participants’ heights. The group means ±SD are shown.

3.3.

Dynamic Propagation of LFOs throughout the Brain

We can similarly correlate the NIRS fingertip Δ[tHb] with the cerebral BOLD signal to determine whether this signal is present in the brain. Video 1depicts the spatiotemporal pattern of the LFO correlations over a range of time lags relative to the NIRS data, demonstrating the propagation of the maximum correlations throughout the brain revealed by RIPTiDe and illustrating several important points. First, the voxels, which are highly correlated with NIRS data, are widely and symmetrically distributed throughout the brain, with the highest correlation appearing in the draining veins, although there is also significant correlation throughout the gray matter. This global signal confirms that a significant portion of the LFO signal in the brain is related to systemic blood circulation variations. Second, the dynamic pattern reflects the variable arrival times of the LFOs at different parts of the brain, just as it arrives at the finger and the toe with different time delays. This latter observation supports the contention that the LFO signal directly reflects bulk blood flow and confirms our previous, brain-only measurements.

To explore the relationship between the LFOs in the brain and the finger, two example points were selected in the brain of one participant, indicated by the green cross and the blue cross in Video 1. These regions were chosen to represent areas with large positive and negative delay values relative to the finger signal. Figure 3(a) and 3(b) illustrates the corresponding BOLD fMRI temporal traces (dashed gray line) together with the NIRS regressors derived from the finger (solid lines): one without any temporal shift (solid gray line) and one with the optimal temporal shift (solid black line), as estimated by RIPTiDe. The temporal shifts of the NIRS LFOs at the two locations differ substantially. They represent the arrival time of the LFOs at these two locations relative to that of the finger. The LFO signal arrived at the indicated voxels 2.16 s earlier and 1.98 s later than it reached the fingertip, respectively, for this participant, suggesting that as in the peripheral NIRS data, the wide temporal shifts of the LFOs exist in the BOLD fMRI. Clearly it is necessary to determine the optimal delay of the LFO signal for each voxel to account for the dynamic evolution of this signal.

Fig. 3

Temporal traces of blood oxygen level-dependent functional MRI (BOLD fMRI) and comparison with that of NIRS at the fingertip. In one participant, two points in the brain were chosen as shown in the green cross and the blue cross in Video 1; the corresponding BOLD signals (dashed line) are depicted in (a) and (b) together with the temporal traces of Δ[tHb] obtained by NIRS at the fingertip (solid gray line) and its temporal-shifted version (solid black line). A.U., arbitrary unit.

Figure 4 shows the maximum $z$-statistic map for each participant, registered to the MNI152 brain⁴¹ in axial, sagittal, and coronal slices. The activated voxels represent areas where global circulatory LFOs accounted for a significant portion of the signal variance at some time lag within the tested range. The supra-threshold area is large and symmetric. The maximum $z$-statistic map summarizes the correlation maps at 181 different time-shift values. This map simply combines and displays the highest significance values found without conducting any statistical analysis.

Fig. 4

Max $z$-statistics maps of each participant overlaid on the standard brain showing areas highly correlated with LFOs collected at the fingertip. The color bar indicates the colors corresponding to the $z$ value from 3 to 6.

3.4.

Coherence between Peripheral LFO and BOLD fMRI

Figure 5 shows, for each participant, the averaged coherences over all the voxels [Fig. 5(a)] with a blue cross indicating the maximum and the corresponding coherence map at the frequency of the maximum averaged coherence [Fig. 5(b)]. From Fig. 5(a), the frequency of highest coherence between peripheral LFO and BOLD fMRI were participant specific. Only three participants (1, 5, and 6) showed clear peaks, at 0.029, 0.097, and 0.082 Hz, respectively. For the rest of the participants, there was no clearly dominant frequency. The patterns of high coherence at the dominant frequency, shown in Fig. 5(b), are similar to the patterns shown in the maximum $z$-statistic maps (Fig. 4).

Fig. 5

Averaged coherence of each participant calculated between the Δ[tHb] at the fingertip and the BOLD from each voxel (a), blue cross marks the maximum coherence point. Corresponding map of coherence of each participant obtained at the frequency of max coherence (blue cross). The color bar indicates the colors corresponding to the coherence values from 0.3 to 0.7.