Study of complex hemodynamic fluctuations in the human brain by simultaneous near-infrared spectro-imaging and functional magnetic resonance imaging

In this paper we discuss temporal and spatial patterns of brain hemodynamics under rest and motor stimulation conditions obtained by functional magnetic resonance imaging and simultaneous fast multi-channel near-infrared spectro-imaging in the human motor cortex. Our data indicate that the main difference between the brain hemodynamics under the repetitive stimulation and the rest conditions is not in the appearance of hemoglobin concentration changes during the stimulations (since fluctuations occur at rest as well), but in their more regular, i.e. phase-synchronous with the stimulation behavior.


INTRODUCTION
The main advantages of functional magnetic resonance imaging (fMRI) in brain mapping studies are high spatial resolution and the absence of any penetration limits. However, the blood-oxygen-level-dependent (BOLD) effect is a complex biophysical phenomenon 1,2 . The advantages of near-infrared spectroscopy (NIRS), such as non-invasiveness, high temporal resolution and relatively low cost, make it an effective method for studying the dynamics of physiological processes, particularly brain hemodynamics 3 . The aim of our study was to compare functional cerebral hemodynamic signals obtained simultaneously by near infrared spectroscopy (NIRS) and by functional magnetic resonance imaging (fMRI), and to explore the possibility of detecting hemodynamic changes in the brain using near-infrared optical signals. In this paper we discuss temporal and spatial patterns of brain hemodynamics under rest and motor stimulation conditions obtained by fast multi-channel NIRS measurements in the motor cortex area. Particularly, we study the relationship between underlying hemodynamic fluctuations and responses to stimulation. The exercise epochs consist of the repetition of stimulation/relaxation cycles. The relative duration of stimulation and relaxation is different for different epochs. Recently a new time-domain method of data analysis was developed, which allows detection and quantification of phase synchronization between two signals 4 . The phase synchronization analysis is important because, as shown in Ref. 4 , phase synchronization is not equivalent to coherence or frequency synchronization, being an independent characteristic of interrelationship between two processes. We compare dynamic features and optical maps of hemodynamic signals obtained using different data analysis techniques, including both the traditional folding average and Fourier transform methods, as well as the new quantitative phase synchronization analysis.

INSTRUMENTATION AND METHODS
Magnetic resonance imaging was performed using a 1.5 Tesla whole body MR scanner (Signa, General Electric Medical Systems, Milwaukee, WI) equipped with echospeed gradients and a standard circularly polarized birdcage head coil. Sagittal T1-weighted localizer scans were used to determine the correct plane for the functional scans. Gradientecho echo-planar images were acquired using a data matrix of 64 x 64 complex points, TR=1280 ms, TE = 40 ms, FOV = 240 mm, slice thickness = 7 mm, no inter-slice gap, receiver bandwidth 62.5 kHz, and flip angle 60o. The voxel size was 3.75◊3.75◊7 mm3. Multi-modality radiological markers (IZI Medical Products Corp, Baltimore, MD) were embedded into the optical sensor to facilitate correct orientation of the MRI slices with respect to the sensor and to enable recovery of the sensor orientation for data analysis.
We use a two-wavelength instrument for near-infrared probing of tissues in which light emitted by laser diodes (758 and 830 nm) is guided to the tissue through multi-mode silica optical fibers (400 µm core diameter). Twelve laser diodes (six per each wavelength) operate in a sequential multiplexing mode with 10 ms on-time for each diode. Two glass fiber bundles (3.2 mm internal diameter) collect the scattered light and conduct it to the photomultiplier tube (PMT) detectors. In order to obtain cerebral optical maps, we designed a probe shown in Fig.1. Two detector fibers were fixed in the central part of the probe. The paired (758 and 830 nm wavelength) source optodes were attached to the probe pad at 6 positions. In addition to the optical signals probing the brain hemodynamics, we acquire the heart rate, and the arterial saturation by means of a pulse oximeter N-200 (Nellcor) with the sensor attached to the left hand index finger, and the respiratory signal with the monitoring system Resp-EZ (Sleepmate/Newlife Technologies). All these physiological signals are acquired by the PC computer simultaneously with the near-infrared signals.
A written informed consent was obtained from each subject before measurements. The probe was positioned on the left side of the head in such a way that the point shown in Fig. 1 with the cross mark coincided with the measured middle point (C3-position) of primary motor cortex. This point was found as the one situated on the line connecting the vertex and the left earlobe, 4 cm down from the vertex. During measurements the subjects were comfortably supine and instructed not to speak or to make unnecessary movements. Each measurement consisted of the rest epoch (10 min) and three exercise epochs. During the rest epoch the baseline data were acquired. During the exercise epoch, subjects were asked to begin or stop performing a finger-motion (a light palm squeezing) exercise by the right hand. The squeezing rhythm (1.5 Hz) was maintained by means of a metronome. Exercise epochs E1, E2, and E3 differed by the duration of the stimulation/relaxation period: 20/60 s, 20/20 s, 10/17 s, respectively, and consisted of 5, 10 and 10 periods, respectively.
Mathematically each exercise epoch can be associated with a rectangular stimulation wave, whose magnitude is zero during relaxation phases and one during stimulation phases. We used such stimulation rectangular waves in the analysis of coherence and phase synchronization between the stimulation and hemoglobin signals. We converted optical intensity data into hemodynamic concentration changes assuming a model for light transport in a strongly scattering medium based on the diffusion approximation to the Boltzmann transport equation, and that the absorption variations in the tissue are only due to oxy-and deoxy-hemoglobin 5 .
We analyzed the autocorrelation power spectra of the oxy-and deoxy-hemoglobin signals. The autocorrelation power of the data series (for example, of [Hb] k or [HbO 2 ] k series, where k is the light channel) at frequency ν can be obtained by averaging the squared magnitudes of complex Fourier transforms ) ( f k ν of successive (possibly overlapping) data subsets . The bandwidth and the frequency resolution of the autocorrelation and coherence spectra are determined by the sampling rate and by the length of the data subset, respectively. For the measurements reported in this work the acquisition time per point is 160 ms, which corresponds to a frequency bandwidth of 3.125 Hz. The duration of each subset data trace to produce a Fourier-transform is 164 s, which corresponds to a spectral resolution of 0.0061 Hz. We estimated the phase synchronization strength between the stimulation rectangular wave and the oxy-and deoxyhemoglobin signals. According to the definition 19 , two periodic quasi-harmonic signals having phases φ 1 (t) and φ 2 (t) (defined on the infinite interval) are phase-locked, or phase-synchronized if there are two integer numbers n and m such that |n φ 1 (t) -m φ 2 (t)| < const. For periodic, quasi-harmonic signals the phase synchronization condition is equivalent to the frequency locking condition n ω 1 = m ω 2 , since for periodic oscillators ) t ( (where the bar denotes the time average). Recently, phase synchronization was understood as a specific relationship between two signals of arbitrary nature, including non-periodic and noisy signals 4 . This interpretation required a generalization of the concept of phase and of the mathematical condition restricting the relative phase change of two signals. For an arbitrary real function of time, F(t), the phase may be defined as the phase of the complex analytical continuation of the given function into the complex plane. The real part of this analytical function is given by the original function F(t), and the imaginary part Q(t) is given by the Hilbert transform of F(t): The condition |n φ 1 (t) -m φ 2 (t)|< const, where φ 1,2 (t) are the phases of analytical continuations of two signals taking their values on the interval ) , ( ∞ −∞ , may still be considered as a phase-locking condition. However, as it was shown in 4 , in the case of noisy signals this condition may not be always fulfilled. Instead, the authors of 4 proposed a characterization of the phase synchronization strength by means of the statistically computed phase synchronization index (PSI) nm η . It is important to note that in the case of a strong (n:m)-phase synchronization characterized by the integer numbers n and m, there is also a strong (n:m)-frequency locking in terms of the sharpness of the distribution for the function Ω nm = nω 1 are the instantaneous frequencies and ) t ( k φ are the phases of the complex analytical signals. However, this statement may not be generally inverted. In other words, frequency synchronization may indicate phase synchronization, but is not equivalent to it. Thus, two signals may be frequency-locked, but not phase-synchronized. Furthermore, as it was shown in Ref. 4 , the (1:1) -phase synchronization between two signals is not equivalent to their coherence: the signals may be coherent, but not phase-locked.
In practice, in order to detect and measure the degree of phase synchronization one should try different pairs of integer numbers n and m to find the pair which gives the maximal value of nm η . Note that since nm η measures the deviation of a given distribution of phase difference values from the uniform distribution, any two finite data sets exhibit some non-zero PSI values. Therefore, to avoid spurious detection of synchronization between two signals represented by the time series of the length L, one should disregard nm η values that are below the significance level determined by the statistical distribution of PSI values corresponding to all possible L -length sets of random data having homogeneous distribution on the interval (-π,π) in the limit ∞ → L . In order to calculate PSI, we filtered the time series using a non-recursive digital filter having a pass-band centered at the repetition frequency of stimulations, heart rate, or breathing rate. Then we performed the numerical Hilbert transform 4 of filtered signals to obtain imaginary parts of the corresponding complex analytical signals. After that, the relative phase ) t ( nm Ψ calculated from the phases of analytical signals and reduced to (-π,π) interval was analyzed statistically to obtain PSI values. Only significant values were considered. We assumed the significance level to be equal to the mean plus one standard deviation of the PSI value distribution for all possible equiprobable random data sets of the same length as the analyzed one.
We used values resulting from the phase synchronization analysis for mapping brain motor activity. The optical maps were made of six pixels corresponding to the squares labeled as 1 through 6 in Fig. 1. The value of each pixel was obtained using data from the corresponding 2.8 cm source-detector pair. The six pixels of an optical map correspond to six zones in the motor cortex area. . Spatial maps of brain oxy-and deoxy-hemoglobin signals during epochs E1, E2, and E3 in terms of phase synchronization with the stimulation wave. The plotted value is PSI for the 2:1 type of synchronization (epoch E1), and 1:1 type of synchronization (epochs E2 and E3). Figure 2 shows power spectra of fluctuations in the cerebral [HbO 2 ], arterial pulse, respiration and heart rate in the frequency range between 0.01 Hz and 10 Hz as obtained by NIRS on a healthy right-handed male subject at rest conditions. A typical [HbO 2 ] power spectrum, as the one shown in Fig. 2, includes a heartbeat peak and the structures towering around the average respiration frequency and some frequency lower than the respiratory one, which for the given subject are at about 1.1 Hz, 0.3 Hz, and 0.1 Hz, respectively. One can identify the respiratory structure in the [HbO 2 ] spectrum by comparison with the power spectrum of the respiratory monitor signal. Note that the power spectrum of the pulse oximeter heart-rate signal shows peaks both at the respiratory frequency (due to the sinusarrhythmia) and the lower-frequency peak of the [HbO 2 ] spectrum (~0.1 Hz). Therefore, one can believe that the heartrate variability is the factor contributing to the [HbO 2 ] spectrum below 0.5 Hz at rest. In different subjects the average breathing frequency varies approximately between 0.1 Hz and 0.4 Hz. Besides the sinus-arrhythmia, the heart-rate variability may exhibit slow changes in the range of 0.01 -0.1 Hz. The width of the corresponding spectral peaks varies depending on the regularity of the heart rate variations and respiration. Hemodynamic fluctuations in different zones of motor cortex at rest conditions were not identical. However, their coherence magnitude in the frequency band below 0.3 Hz is high and their fluctuation power spectra in this band have similar shapes and magnitudes.

RESULTS
In order to determine the location of tissues contributing to the functional response in optical signals during stimulation, we performed a correlation analysis of the BOLD signals using the negative of the [Hb] signal as a predictor (the negative of the [Hb] signal was used since an increase in BOLD signal should correspond to a decrease in [Hb].) We have found a good spatial collocation between the optical and fMRI signals.
The power spectrum analysis of signals acquired during motor stimulation reveals frequency locking of the [Hb] and [HbO 2 ] time series with stimulation. Since the frequency locking may be an indication of phase synchronization, we calculated the corresponding PSI's. The maps of PSI spatial distribution for the oxy -and deoxy-hemoglobin signals during epochs E1-E3 are presented in Figs.3. The zones showing no statistically significant phase synchronization are labeled with clear bars. Note that the PSI values for the epoch E1 correspond to the (2:1)-phase synchronization. During E1 no statistically significant synchronization of any (n:m)-type was detected in any of [Hb] signals. Our analysis showed that the phase synchronization maps always show a localized response, while relatively large areas in the brain may exhibit fluctuations that are coherent with the stimulation in terms of temporal correlation.

CONCLUSION
In all subjects we found a good collocation of the brain activity centers revealed by both fMRI and NIRS. We also found a high temporal correlation between the BOLD fMRI signal and the deoxy-hemoglobin concentration (NIRS) in the subjects who exhibited low fluctuations in superficial head tissues.
Using multichannel near-infrared cerebral spectroscopy we have found that the hemoglobin changes in the motor cortex under periodic motor stimulation can be highly regular or irregular depending on the duration of stimulation/relaxation period. We have shown that in the case of complex hemodynamics, the response to stimulations can be detected combining statistical, folding average, power spectrum, coherence, and quantitative phase synchronization analyses. Using these techniques, we found that the oxy-hemoglobin concentration changes are significant and phasesynchronous with the stimulation at most of the exercise conditions. The deoxy-hemoglobin response to stimulations depends on the stimulation and relaxation timing, partially due to the interference with the background fluctuations, and partially due to the possible stimulation-response phase synchronism. Using the power spectrum, coherence and phase synchronization analysis, we have shown that functional stimulation can cause local frequency-and phasesynchronization of cerebral hemodynamic fluctuations. We have shown that while a relatively large area of the brain may exhibit fluctuations that are coherent with the stimulation wave, when phase synchronization occurs, the phase synchronization maps always show a localized response both in the oxy-and deoxy-hemoglobin signals.
Significant hemodynamic fluctuations in the brain occur not only during exercises, but also at rest. Our results indicate that the main difference between the brain hemodynamics under the repetitive stimulation and the rest conditions is not in the appearance of hemoglobin concentration changes during the stimulations (since fluctuations occur at rest as well), but in their more regular, i.e. coherent or phase-synchronous with the stimulation behavior.