Diffuse optical imaging (DOI) is a spectroscopic method capable of noninvasively measuring changes in the concentrations of oxy- and deoxyhemoglobin in the human brain during functional activity.1, 2 As a tool for neurology, this technique offers several additional attributes that complement other existing imaging methods, such as functional MRI (fMRI).3 MRI and optical imaging complement each other by their spatial, temporal, and spectroscopic resolutions. In recent years, this has led to an increased interest in multimodality functional human imaging4 and image reconstruction techniques.5, 6, 7 However, as these multimodality methods are developed, careful studies must be performed to explore the underlying biophysical relationships between the different imaging modalities. Properly interpreting data together from both modalities requires an understanding of the relationships between the underlying physiological hemodynamic states and the physical measurements obtained from the instruments and, with that, experimental validation of the biophysical models supporting each modality.
In recent years, a number of studies have examined these relationships through investigations of the temporal correspondence between DOI and blood oxygen level-dependent (BOLD)-fMRI, 8, 9, 10, 11, 12 DOI and arterial spin labeling (ASL)-fMRI,12 or DOI and cerebral blood volume (CBV)-fMRI.11 These and many other collected works have attempted to reconcile theoretical models describing the relationship between these fMRI contrast mechanisms and the physiological deoxyhemoglobin, blood volume, and blood flow parameters through the measurements recorded by optical imaging. While it is clear that the differing vascular sensitivities of the modalities may play a role in precise details of this relationship, the conclusion of many of these studies has been that the temporal dynamics of measurements by fMRI and optical imaging are consistent with their respective theoretical biophysical underpinnings. However, although there have been many such investigations, which have explored the temporal relationships between these two modalities, little has been published exploring the spatial correspondence between fMRI and DOI in a quantitative way. Although a few qualitative reports suggest agreement between the spatial profiles of fMRI and optical methods,13 no detailed and quantitative comparisons have been published. This lack of quantitative comparisons might be attributed to the ill-posed problem of image reconstruction by optical imaging. Optical imaging is based on the topographic reconstruction of a discrete set of absorption change measurements between pairs of optical source-detector positions on the head. This reconstruction is generally underdetermined, with far more unknown parameters than actual measurements. Tomographic (three-dimensional) reconstruction is even more ill-posed, since the added depth degree of freedom adds more unknowns, as well as, an exponentially decaying sensitivity to absorption changes. Over the years, several groups have continued to make progress with better reconstruction algorithms and methodologies, which exploit anatomical 6, 7, 14, 15 or functional5 MRI information. By improving the spatial localization of reconstructed optical images, spatial priors may help to improve the quantitative accuracy of DOI. However, for these reconstructions, the method and amount of regularization used to constrain the solutions will have a quantitative effect on the resulting image. As a result, it has been difficult to quantitatively compare the spatial profiles of fMRI and DOI independently of reconstruction technique. Without this independent assessment, we cannot be confident in the use of fMRI information to guide the optical inverse problem.
In this work, we perform the first quantitative assessment of the spatial correlation between optical and fMR imaging. To avoid the optical inverse problem, rather than reconstructing images from the optical data using these regularization techniques, we directly tested the fMRI data as a solution to the optical inverse problem. Here, given the anatomical structure of the head from MRI, we use photon-migration theory and Monte Carlo techniques to calculate the Green’s function describing the optical sensitivity profile to absorption changes in the underlying structures of the head (i.e., the optical “forward” model) as described in Boas 16 The overlap of this optical sensitivity profile with the brain activation identified by fMRI provides a means of predicting the response amplitudes, which would be measured optically between each source and detector pair. This approach allows us to account for the decaying depth sensitivity and partial volume errors of the optical measurements and avoid the ill-posed inverse problem by using the optical forward equation in the forward direction. Rather than increasing the dimensionality of the optical data by image reconstruction through regularization, knowledge of the optical measurement model is used to reduce dimensionality of the fMRI data into the source-detector–based (measurement) space of the optical probe, thus allowing quantitative comparisons between DOI and fMRI. Performing comparisons of the spatially and temporally varying hemodynamic response across the probe layout allows us to explore the spatiotemporal relationships between optical and fMRI while avoiding the need for regularized inversion techniques.
In this experimental study, we performed DOI simultaneously with BOLD- and ASL-based fMRI during an event-related finger-tapping task in order to quantitatively test the spatial correlation between the optical and fMRI methods. This data was published in Huppert, 12 who investigated the temporal correlation between the region-of-interest averaged signals. Here, we provide a significant follow-up to that analysis to investigate the spatial correlation between these modalities.
DOI of brain activation is generally based on the measurement of spatiotemporal changes in the absorption of light. Near-infrared light is introduced into the head and propagates through the dense scattering layers of the scalp and skull into the brain. Hemodynamic changes in oxy- and deoxyhemoglobin concentrations affect the absorption properties of the brain and result in a change in the intensity of light as it migrates back out of the head. Using multiple measurements taken between an array of light source and detector positions spaced several centimeters apart, DOI attempts to spatially resolve these changes in hemoglobin. However, the reconstruction of images requires knowledge of the spatial profile of the light propagation through the head, which determines the spatial sensitivity of these measurements. To derive the distribution of photons in a medium with a complex distribution of absorption and scattering properties [ and , respectively], such as the human head, the photon density must be modeled by empirical means using computerized simulations. In Monte Carlo–based modeling, such as is performed in this study, photons are “launched” into the medium. The distribution of these photons is statistically modeled based on the probability of an absorption or scattering event at each region of space as described by and .16, 17, 18
For brain activation, the changes in the optical properties are small, and thus a linear approximation is reasonably accurate for predicting the changes in optical measurements produced by localized changes in the optical properties. The forward model for such optical measurements is of the formis the vector of measured optical signal changes for each source-detector pair; is a vector of the absorption changes in different discrete volume elements; and is a three-point Green’s function matrix describing the linear transformation from absorption changes within the volume to optical signal changes between each measurement pair. Each row of the matrix describes the light propagation between a particular optical source and detector pair, which is often characterized as a “banana-shaped” profile, which describes the spatial sensitivity for that measurement. The matrix is thus a projection operator, which integrates absorption changes over the volume to predict the measurements between sources and detectors. It has been shown that this linear operator is approximated by the adjunct product between the photon density distributions for each given source and each given detector involved in each optical measurement.19 Because the propagation of light depends on the optical properties of the medium, each of these parameters ( , , and ) are dependent on the optical wavelength .
The forward propagation of fMRI
The forward model, shown in Eq. 1, describes the linear transformation of absorption changes, which occur at particular volume locations, to the measurement of these changes between the optodes of the DOI probe. Rather than inverting this equation to achieve an image of these changes as is typical of optical imaging, the forward model can be used to predict how the changes shown by fMRI images should look when measured optically. Thus, Eq. 1 can be used to “forward model” the optical measurements based on the changes measured at high spatial resolution by fMRI. This provides us with a direct way of testing the hypothesis that the optical and fMRI spatial profiles of the response amplitudes are consistent with one another.
However, the direct application of the optical forward operator to project the fMRI image requires information about the relationship between the fMRI signal and optical absorption changes. A change in the fMRI-BOLD signal is proportional to a change in reduced hemoglobin (HbR) and not necessarily proportional to changes in optical absorption at a single given wavelength since the optical absorption is a function of both HbR and oxygenated hemoglobin changes. Thus, Eq. 1 cannot be directly applied without considering how light propagation differs for each optical wavelength. Rather then using the optical forward operator determined at a single wavelength, we must consider a multispectral forward equation.20, 21
In order to derive a forward operator to directly project hemoglobin changes, we assume the absorption is dominated by hemoglobin. Optical absorption changes are linear combinations of the oxy- and deoxyhemoglobin changes, such thatare the extinction coefficients of and HbR, which are wavelength dependent, and indicates the concentration of . The absorption and hemoglobin concentration variables in Eq. 2 are vectors of these values at each discrete spatial location (volume element). Making use of Eq. 2 allows us to write the multispectral forward operator for measurements at two different wavelengths in terms of the perturbations in the oxy- and deoxyhemoglobin concentrations (following the notation of Li 20): is the Kronecker product of and , where the identity matrix has dimension equal to the number of columns of . and are column vectors of the measurement at wavelengths and , where each element in the vector represents a different source-detector pair.
Rather than using Eq. 1 to forward project the fMRI image, the multispectral forward equation given in Eq. 3 should be considered. The BOLD-fMRI signal is hypothesized to represent the spatio-temporal variation in . Equation 3 can thus be used to predict the optical signals and due to a change in deoxy-hemoglobin22, 23
In the modified Beer-Lambert equation [Eq. 5], is the linear distance between the position of a source and detector and is the differential path-length factor, which is a unitless coefficient defined as the effective path length through the head divided by the source-detector seperation.22, 23, 24 Using Eqs. 4, 5, the BOLD signal is changed to a spectrally resolved optical absorption change, projected through the optical forward equation at each measured wavelength, and finally converted back to an estimated concentration change using the modified Beer-Lambert law. Thus, we obtain an estimate of the deoxyhemoglobin changes for each source-detector pair, predicted from the spatial profile of the BOLD signal. This overall procedure is similar to the method detailed by Strangman 24 for examining the errors in the estimate of the hemoglobin changes arising from the use of the modified Beer-Lambert equation. These errors have been found to be negligible for the wavelengths used in our study.24, 25, 26 Similarly, the fMRI-ASL signal depicts changes in blood flow, which we have previously noted has a strong temporal correlation with total hemoglobin changes .12 Again, because this represents absorption changes at multiple wavelengths and must be propagated through a multispectral forward model, we approximate the ASL signal as a spatiotemporal variation in . We predict the optical signals measured for each source-detector pair using Eqs. 3, 4, 5, but now projecting both oxy- and deoxyhemoglobin changes assuming an oxygen saturation of 65%.24, 27
The data used in this analysis has been previously reported in Huppert 12 for comparison of the temporal characteristics of optical and fMR imaging. The methods are briefly summarized here.
In this study, we enrolled 11 healthy, right-handed subjects (8 male, 3 female). The subjects were 20 to 60 years old (mean ). The Massachusetts General Hospital Institutional Review Board approved the study and all subjects gave written informed consent. For BOLD:DOI comparisons (study I), a total of six subjects were scanned. One subject failed to show any significant activation in either modality and was excluded from further data analysis. All of the remaining five subjects (four male, one female) showed significant localized activation patterns in both modalities. In the ASL-DOI experiments (study II), again a total of six subjects were scanned. One of these subjects had to be excluded because of a low signal-to-noise ratio in the DOI measurement, due to poor coupling of the probe to the head. The remaining five subjects (four male, one female) were included in the full analysis. One subject was repeated for both studies (subject B/F).
Prior to placement in the bore of the MRI scanner, subjects were briefed on the details of the finger-tapping task. During data acquisition, they were instructed to watch a visual display from a laptop computer, which was projected from the back of the scanner room onto a screen in the bore of the magnet and was visible from a mirror mounted on the top of the MR head coil above the subjects’ eyes. On this display was shown the text “STOP” or “GO.” For the duration of the GO visual cue, subjects were instructed to tap their thumb and each of the fingers on their right hand sequentially at a self-paced rate (approximately 2 to ). After , the GO text was changed to the STOP command; at which time, subjects were instructed to stop tapping and keep their hand relaxed. Aside from the finger-tapping, subjects were otherwise instructed to remain motionless for the duration of the scanning session.
During the scan, the interstimulus interval (ISI) between finger-tapping periods was pseudorandomly chosen and optimized to provide the uniform temporal coverage necessary for deconvolution with a 500-ms time step.28 The length of the ISI ranged from 4 to with an average ISI period of . The timing of the stimulus presentation was synchronized with the MR image acquisitions and generated with a custom-written MATLAB script (Mathworks, Sherborn, Massachusetts). Each run lasted and consisted of between 27 and 32 stimulus periods. This was repeated 4 to 6 times for each subject during the course of 1 scan session (approximately ) for both the ASL and BOLD measurements. The experimental procedure, subject instructions, and stimulus parameters were identical for both the BOLD and ASL functional scans.
DOI Acquisition and Processing
In these experiments, we used a multichannel continuous-wave optical imager (CW4, TechEn Inc., Milford, Masschusetts) to obtain the measurements as previously described.29 The DOI imager has 18 lasers—9 lasers at and 9 at —and 16 detectors. Of these, only four source positions (two wavelengths each) and eight detectors were used here. Laser wavelengths were chosen to minimize cross talk between the two hemoglobin species at and ( and , respectively).24, 25, 26
For detectors, the imager employs 16 avalanche photodiodes [(APDs), Hamamatsu C5460-01), each of which is digitized at approximately . A bandpass filter follows each APD module with a cut-on frequency of approximately to reduce noise and the 60-Hz room light signal, and a cutoff frequency of approximately to reduce the third harmonics of the square-wave signals. The signal is then passed to a programmable gain stage, which is used to match the signal levels with the acquisition level on the analog-to-digital converter within the computer. A digital demodulation and low-pass filter is used off-line to separate the individual source signals on each detector channel, which were frequency encoded as described in Franceschini 29 This low-pass filtering was done via a seven-pole infinite-impulse response filter Butterworth design with a 10-Hz passband.
The DOI probe is made of flexible plastic strips with plastic caps inserted in it to hold the ends of the 10-m source-detector fiber optic bundles. The probe consists of two rows of four detector fibers and one row of four source fibers arranged in a rectangular grid pattern and spaced between nearest neighbor source-detector pairs as demonstrated in Fig. 1 . This separation distance was chosen to give enough depth penetration to provide sensitivity to the brain.30 The plastic probe was secured to the subject’s head via Velcro straps and foam padding and centered over the contralateral primary motor cortex (M1). The 10-m fibers were run through the magnet bore to the back of the scanner and through a port into the control room where they were connected to the DOI instrument. The CW4 instrument was operated from the MR control room.
After collection, demodulation, and low-pass filtering at a 10-Hz bandwidth to separate individual source contributions, the data was further processed using a custom MATLAB data analysis program (HOMER).31 Signals were further low-pass filtered at using a zero-phase forward and reverse digital filter to remove the heart signal. Changes in optical density for each source-detector pair were then high-pass filtered with a 1/30-Hz drift correction and finally converted to change in hemoglobin concentrations using the modified Beer-Lambert relationship22, 23, 32 with a differential path-length correction calculated from the Monte Carlo simulations (Sec. 2.5).
From the concentration data, the individual subject hemodynamic responses were calculated using an ordinary least-squares (OLS) linear deconvolution and implemented within the HOMER program (e.g., Ref. 33). No nuisance variable polynomial trends were fit in the linear regression, since the low- and high-pass filtering steps had already removed these. In contrast to the generalized linear modeling (GLM) methods commonly used in fMRI analysis, which make use of assumed functional forms of the hemodynamic response [i.e., SPM (Frackowiak 34)], OLS deconvolution estimates the entire response function and thus does not assume any fixed canonic responses (i.e., Gaussian or gamma variants). Because both the fMRI and DOI time-series analyses are preformed in this way, this allows comparison of the response functions from both modalities without concerns of introducing errors that might have arisen from such assumed response functions. Epoch timing was synchronized to the MRI images so that subject response times were kept constant between the fMRI and DOI. Regions of time showing significant motion artifacts (as clearly evidenced by extremely large and sudden perturbations in the measurement time course) were rejected from the analysis. The impulse response functions were then averaged across all runs for each subject taken within the same scan session for each source-detector pair.
fMRI Acquisition and Processing
For study I, BOLD-fMRI measurements were preformed using a Siemens Allegra MR scanner (Siemens Medical Systems, Erlangen, Germany). fMRI data were collected using a gradient echo-planar imaging (EPI) sequence with five 6-mm oblique orientation slices (1-mm spacing) and 3.75-mm in-plane spatial resolution. Structural scans were performed using a T1-weighted magnetization-prepared rapid acquisition gradient echo (MPRAGE) sequence ( resolution, ]. To calculate the BOLD-based hemodynamic response functions, the functional images were first motion-corrected35 and spatially smoothed with a 6-mm Gaussian kernel. The response functions were then calculated by an OLS linear deconvolution. A third-order polynomial was included to remove drift effects. For each subject, the effect and standard deviation maps were input to a mixed effects analysis to generate a map of statistics.36 The mean of the response from 2 to was used to calculate the effect size.
For study II, ASL-fMRI was carried out at (same scanner as study I) using proximal inversion with control for off-resonance effects (PICORE) labeling geometry37 with Q2TIPS saturation38 to impose a controlled label duration. A postlabel delay of and a label duration of were used, with repetition and echo times of and , respectively . The PICORE labeling scheme allowed the collection of BOLD signals using the control phase of the acquisition. EPI was used to image five 6-mm slices (1-mm spacing) with 3.75-mm in-plane spatial resolution. Structural scans were also performed with the same scan prescriptions as study I. Although the original image acquisition was at , the stimulus had been jittered evenly on a 500-ms time step, which allowed for the response to be calculated at , albeit with a lower signal-to-noise ratio. To estimate the flow response, the ASL-functional scans were first separated into the control and negatively labeled tag images. These traces were then independently interpolated to upsample the data points to using a cubic spline model.39 Subtracting each negatively labeled tag scan from the immediately subsequent control scan generated the flow image series. This flow series was then deconvolved and the region of interest average calculated similar to its calculation for the BOLD response. The mean of the response from 1 to was used to calculate the effect size. Figure 2 depicts the processing stream.
Monte Carlo Simulations
The anatomical MPRAGE MRI scans taken at the start of each experimental session were resampled to 1-mm (isotropic) resolution and then segmented into white and gray matter using anatomical reconstruction tools available as part of the FREESURFER package from the Massachusetts General Hospital (MGH).40 This process first consisted of intensity normalization, skull stripping, and white-gray matter differentiation as described by Dale 41 and implemented in the RECON-ALL program as part of the MGH FREESURFER package. All parameters were set to their default values. To further differentiate the anatomical images into skin, skull, and cerebral spinal fluid (CSF) layers, a similar additional intensity-based initialization followed by sorting of ambiguous voxels based on the iterative expectation-maximization algorithm was applied to the unassigned regions of the intensity normalized image that were not accounted for after the skull-stripping algorithm.42, 43 For this, a custom written routine was used, which required manual supervision of this segmentation process. After full segmentation, tissue assignments were manually inspected. This process was repeated individually for all the subjects.
To simulate photon migration through the segmented brains, the probe positions were first determined from the location of the vitamin E fiducial markers that had been placed on the DOI probe in the MPRAGE scans. After locating the probe positioning, Monte Carlo–based simulations of light propagation were performed as described in Boas 16 The optical properties of each of the segmented tissues that were used in these simulations are listed in Table 1 .24 In these simulations, the optical properties of gray and white matter were taken to be the composite properties for brain tissue, which is optically dominated by gray matter. For each optode position, the trajectories of photons were simulated at both the 690-nm and 830-nm wavelengths. Because the amount of detected light decays exponentially with distance, only the three-point Green’s functions for the nearest-neighbor source-detectors were calculated. These simulations were also used to determine the differential path-length corrections for each wavelength and source-detector pair and applied to the calculation of concentration changes in the modified Beer-Lambert law used in Eqs. 5, 6.
Optical absorption and scattering properties of the tissue types used in the Monte Carlo simulations. Monte Carlo simulations of light propagation were performed on the tissue segmented structural images for each subject using the tissue optical properties at both wavelengths recorded (830 and 690nm ). Optical values for the various tissue types were taken from Ref. 24.
|Tissue||μa (mm−1)||μs′ (mm−1)||μa (mm−1)||μs′ (mm−1)|
|Brain (Gray and White)||0.0186||1.11||0.0178||1.25|
In order to test the spatial correspondence of fMRI activation with the DOI measurements, the fMRI hemodynamic responses were multiplied by the optical forward matrix to project the changes into source-detector space as described in Sec. 2.5. This was done by first registering and resampling the forward matrices from the 1-mm isotropic resolution of the anatomical images and Monte Carlo simulations to the lower spatial resolution of the functional scans. The forward matrix is then multiplied onto the data matrix holding the fMRI hemodynamic time courses from the BOLD or ASL scans. In order to reduce the noise contribution from the fMRI, the activation maps from the fMRI were used to first mask the fMRI data prior to projection. Although, only voxels under this criterion were used, this high choice of value should be considered quite lenient. The optical probe extended beyond the spatial coverage of the fMRI scan, which spanned roughly only axially [5 slices (6mm 1-mm gap . In order to avoid the introduction of partial volume errors, the projection and analysis were limited to only those source-detector pairs for which the forward operator was spanned by the functional MRI sample volume (i.e., the superior row of source and detector pairs).
In Fig. 3 , we illustrate the concept of the fMRI projection by showing the overlay of the three-point Green’s function (shown as contour lines in logarithmic scale) and the anatomical MRI volume. The functional activation region (shown only above half-maximum amplitude) is inset, showing peak activation approximately below this particular source-detector pair. The forward equation sums these changes over the volume using the weights of the sensitivity matrix (Green’s function) to predict the measurements between each source-detector pair.
Because the BOLD signal is proportional to deoxyhemoglobin (but this proportionality is a matter of debate),12, 44, 45 the resulting deoxyhemoglobin predictions need to be normalized. A single scaling factor, to normalize to unity the maximum change across the entire probe, was applied for each subject.
Although we report the results from using a spectrally defined forward operator [Eqs. 3, 4], we also examined the results using the Green’s functions at both 690 and . Due to the coarse spatial resolution of the fMRI relative to any wavelength-specific variations noted in the Green’s functions, the results were nearly identical with all three choices of forward operators. For consistency, we report only the spectral operator.
To compare the optical measurements and their predicted values from the fMRI, we used linear regression analysis to report correlation values for the spatial or temporal amplitude profiles. Pearson coefficients (i.e., values) describe the zeroth-order correlation between two vectors of values. The associated value for this test describes the probability that random chance alone could produce a correlation equal or better than this value.46
Because Pearson coefficients are nonlinear, to compare the correlations between hemodynamic parameters, we use a transformation of the value. The transformation is defined astest of their transformed values.46
Depth Sensitivity of Optical Imaging
The differential path-length factors (DPF) for each source-detector pairing were calculated from the Monte Carlo data from each subject. The DPF for a given source-detector pair was calculated as the actual mean path length traveled by all photons leaving the source and reaching the detector, divided by the linear source-detector separation.22, 23, 24
In addition to allowing us to forward model the fMRI signal into the optical imaging space, the determination of the three-point Green’s functions in this work allows us to comment on the depth sensitivity of the optical measurements. We calculated the mean activation depth of the BOLD and ASL responses for each subject by computing the first moment of distance to the nearest point on the head’s surface. The mean depth of the optical contrast was calculated in a similar way by weighting by the appropriate Green’s functions . That is,8, is the distance from volume position to the nearest point on the surface of the head . is the value of the Green’s function for the th source-detector pair and is the mean value of the contrast (from 2 to for BOLD or 1 to 6 for ASL) at volume position .
Study I: BOLD/DOI
In study I, we performed simultaneous optical and fMRI-BOLD recordings. Individual subject measurements were gathered and processed as described in the previous sections.
In Fig. 4 , we show the results from the projection of the BOLD signal through the optical forward equation, as described in Eq. 3, for each of the five subjects. The time course for each (nearest-neighbor) source-detector pair is shown for both the DOI HbR recording and the projected BOLD trace. Both of these time courses have been normalized to their respective maximal time points across all channels, such that the relative amplitudes between each pair are comparable. Both time courses are shown as positive changes for direct comparison.
In Fig. 4, we note that there is generally good spatial and temporal agreement between the BOLD and optical signals. In several of the subjects, we see the characteristics of poststimulus undershoots in both the BOLD and HbR data. Among the channels that did not significantly respond to the stimulus, we still see some degree of correlation in the time courses, suggesting that the noise in these channels has a physiological source common to both modalities.
In order to further investigate the spatial correlation between the optical and BOLD measurements, we calculated the maximum response amplitude for each source-detector pair. This is presented as a bar graph in Fig. 5 . We also calculated the mean response over the time period from 2 to for the HbR and BOLD responses, yielding results similar to those found using the maximum of the response (data not shown). The comparison of response amplitudes presented in Fig. 5 shows more clearly the spatial correspondence of the BOLD and optical signals. The distribution of response amplitudes across the probe is similar in both the measured optical data and predictions from the BOLD signal. In addition, we also performed a similar analysis to compare the BOLD and optically measured total hemoglobin (HbT) and responses. The results of the spatial correlation of BOLD with the optically measured hemoglobin parameters are given in Table 2 .
Spatial correlation (study I) results show the spatial correlation of BOLD to the optically measured HbR, HbO2 , and HbT parameters. The BOLD signal was projected through the HbR forward model (shown in Fig. 2) and then the mean amplitude for each channel was calculated (shown in Fig. 3). These amplitudes were then correlated between the two modalities. Displayed here are the correlation (Pearson’s) coefficients for each comparison for the five subjects used in study I. The probability ( p value) of each correlation is shown in brackets for each comparison.
The spatial correlation between the BOLD and optically measured HbR amplitudes was strongly significant ( , ) across all five subjects. Additionally, the BOLD/HbR spatial correlation was better than that observed between or BOLD/HbT . A paired test on these values allows us to test the probability that the correlation between BOLD/HbR is statistically better than that between and BOLD/HbT ( -transformed test).
Study II: ASL:DOI
In the second study, we repeated the same motor task experiment and simultaneously recorded ASL-based fMRI and optical signals. These measurements were analyzed similarly as in study I. In Fig. 6 , we present in optical measurement space, the time courses for the fMRI-measured blood flow, which has been projected from the ASL-fMRI data by the simulated optical forward matrix [Eq. 3]. The amplitudes have been scaled to match the single maximum amplitude of the fMRI and optical time courses, allowing for a comparison of the relative amplitudes across space and time. Here we see a good correspondence between those source-detector pairs predicted to have the largest response from the fMRI measurements and the optical data.
In Fig. 7 , we present the maximum response amplitudes across all source-detector pairs to clarify the spatial correlation between the fMRI and optical measurements. The calculated correlations are presented in Table 3 .
Spatial correlation (study II) results show the spatial correlation of ASL to the optically measured HbR, HbO2 , and HbT parameters. The ASL signal was projected through the HbO2 forward model (shown in Fig. 4) and then the mean amplitude for each channel was calculated (shown in Fig. 5). These amplitudes were then correlated between the two modalities. Displayed here are the correlation (Pearson’s) coefficients for each comparison for the five subjects used in study II. The probability ( p value) of each correlation is shown in brackets for each comparison.
The spatial correlation between the ASL and optically measured or HbT amplitudes is strongly significant ( , and , , respectively) across all five subjects. The ASL/HbT and spatial correlations are better than that observed between ASL/HbR ( , ). These correlations between and ASL/HbT were determined to be nearly significantly ( and ) better than the correlation between ASL/HbR ( -transform -test).
Since the pulsed ASL scan used in this study allows the extraction of BOLD time courses from the interweaved control images, we were able to repeat the analysis to examine the BOLD/HbR correlation from this set of subjects. Although not shown, these results are similar to those from study I, which used a dedicated BOLD scan. We found that the spatial correlation between BOLD and HbR was 0.43 . This result is slightly weaker than that reported in study I, due to the lower contrast-to-noise ratio of the BOLD signal inherent in the scan sequence . However, similar to the results presented in Table 2, the correlation between BOLD and HbR is higher than that between and BOLD/HbT ( -transform -test).
Comparison of BOLD versus ASL-based fMRI
Our results from Secs. 3.1, 3.2 suggest that the spatial profiles of HbR and are distinct, since these measurements correlated better to the BOLD and ASL activations, respectively. To further investigate this result, we examined the spatial comparison between BOLD and ASL directly. Because the pulsed ASL sequence used in study II interweaves unlabeled, control images among the inversion labeled flow images, this allowed the extraction of both BOLD and flow-weighted time courses. Statistical maps of the activations from both of these time-courses were calculated.
As reviewed by Detre and Wang,47 ASL and BOLD imaging have different spatiotemporal characteristics. These two fMRI methods, in particular BOLD imaging, depend on details of the vascular anatomy, which in turn affect the measurements and the resulting images. It was therefore not surprising that we observed a slight disparity between the two activation maps from these two fMRI modalities; this is consistent with similar observations by Silva, 48 Luh, 49 and Yang 50, 51 Shown in Fig. 8 are the ASL and BOLD activation maps for one such subject (F).
Spatial comparison of ASL and BOLD activation profiles. For each of these five subjects, we calculated the center of mass of the BOLD and ASL activations. This was done by calculating the positional moment weighted by the normalized p -value map as described in the text. The second radial moments were also calculated to describe the area of the activation.
|Subject||Activation Center (First Moment)||Activation Area (Second Moment)||Activation Center (First Moment)||Activation Area (Second Moment)||Displacement|
|F||[127, 56, 142]||[132, 54, 139]|
|G||[128, 52, 132]||[130, 47, 136]|
|H||[128, 67, 140]||[128, 65, 141]|
|I||[129, 52, 137]||[132, 48, 144]|
|J||[125, 35, 152]||[128, 35, 142]|
For each of these five subjects we calculated the “center of mass” of the BOLD and ASL activations (shown in Table 4 ). This was done by calculating the positional moment weighted by the normalized -value map (i.e., . The second radial moments were also calculated to describe the area of the activation . As shown in Table 4, we observed a displacement between the BOLD and ASL centers of activation (mean ).
Temporal Correlation between DOI and fMRI
In our previous work,12 we reported a strong temporal correlation between the region-of-interest averaged time courses of the optically measured HbR and BOLD signals and between the ASL and HbT signals. However, in that analysis, we did not take into account the effects of DOI sensitivity or partial volume effects across the probe. Previously, these region-of-interest averages had been calculated from the inclusion of all source-detector pairs or fMRI voxels that met a criterion. In this current approach, by using the DOI forward matrix to project a prediction of each source-detector measurement from the BOLD signal, we have the ability to better address these effects. The projection of the fMRI onto the optical probe using the optical forward model accounts for these partial volume effects across the probe. One aspect of this current analysis, however, is that the resulting time courses have a smaller contrast-to-noise ratio because the projection preferentially weights the fMRI response from a small number of superficial voxels.
The mean temporal response from all source-detector pairs from the top half of the DOI probe for all subjects was used to calculate the response for each modality. The resultant average time courses are shown in Fig. 9 . Consistent with the results from our previous work,12 the results presented in Table 5 show statistically strong temporal correlations between the BOLD and optically measured HbR responses ( , ) and again between the ASL and HbT responses ( , ). The BOLD/HbR correlation tends to be higher than that observed between the or BOLD/HbT. This trend, however, is not significant ( and , respectively) ( -transform -test). Similarly, the ASL and HbT correlation is stronger than that between ASL and HbR. This result was nearly significant in a test of the paired -transformed Pearson’s coefficients . The correlation between ASL and HbT is also higher than that between , although, this trend is again not significant ( -transform -test).
Temporal correlation results show the temporal correlation of BOLD and ASL to the optically measured HbR, HbO2 , and HbT parameters. The time courses from the projected fMRI and DOI were averaged over all source-detector channels and then correlated. Displayed are the correlation (Pearson’s) coefficients for each comparison for the five subjects used in study II. The probability ( p value) of each correlation is shown in brackets for each comparison.
Depth Sensitivity of Optical Imaging
We calculated the mean depth of contrast for the fMRI and DOI measurements as defined by Eq. 8. These values are reported in Table 6 . We found that the contrast depth measured by the fMRI was for the BOLD and for the ASL. The mean depth of the optical contrast was shallower at only for both the ASL- and BOLD-derived signals. This depth was approximately equal to the cortical depth (mean distance from the scalp to the cortex), which was . The thicknesses of the superficial layers (skin, skull, and CSF) are also shown in Table 6. The skull thickness and CSF layer varied considerably over the nine different anatomical volumes (B and F are the same subject). In some cases, the mean optical contrast was shallower than the cortical depth. This is due to the lower resolution of the functional MR images compared to the anatomical volume .
Contrast depth and sensitivity comparison between the depth of functional contrast from the surface of the head for the fMRI and optical measurements. These were calculated as described in the text in Eqs. 8a, 8b. This reflects the mean depth of the signal giving rise to each measurement. For comparison, we also show the cortical depth, defined as the mean distance of the cortex from the surface of the head. The thickness of the three superficial layers is given in millimeters. The mean DPF for each subject is shown with the standard deviation of the value over all source-detector pairs.
|Cortical Depth||Functional Depth||Photon Migration|
|Subject||Total [mm]||[Skin, Skull, CSF] thickness [mm]||fMRI Contrast [mm]||Optical Contrast [mm]||Differential Path-Length Factor|
|A||12.2||[2.6, 3.9, 5.7]||BOLD||24.9||13.1||5.46 (0.38)|
|B||14.4||[3.1, 4.6, 6.7]||BOLD||14.9||13.1||5.45 (0.38)|
|C||12.6||[2.7, 2.7, 7.2]||BOLD||22.4||21.0||5.63 (0.25)|
|D||11.7||[2.5, 3.6, 5.5]||BOLD||20.7||7.6||5.00 (0.30)|
|E||12.8||[2.7, 5.9, 4.2]||BOLD||26.8||14.6||4.79 (0.36)|
|F||14.4||[3.1, 4.6, 6.7]||BOLD||14.2||14.6||5.12 (0.33)|
|G||16.9||[3.1, 6.7, 7.0]||BOLD||19.5||19.3||5.53 (0.90)|
|H||13.6||[2.6, 3.1, 7.9]||BOLD||21.6||16.5||5.59 (0.23)|
|I||11.8||[1.8, 4.2, 5.9]||BOLD||12.2||9.9||5.93 (1.11)|
|J||10.8||[2.5, 3.2, 5.2]||BOLD||22.9||10.8||4.99 (0.47)|
|AII||13.0||[2.6, 4.2, 6.1]||BOLD||20.0||14.0||5.35 (0.35)|
Related to the mean sensitivity depth, in Table 6, we also show the mean differential path-length factor for each subject, averaged over all source-detector pairs. Consistent with previously published theoretical findings from anatomical head geometries,24, 30 we find that the DPF was approximately 5.4 and 5.3 for the 690- and 830-nm wavelengths. For subjects G and I, there was a fairly large variation in the DPF between optode positions, although the ratio of the DPFs for each wavelength showed less variation. From visual inspection, it appeared that the position of the probe on these two subjects was slightly more posterior than on the other subjects. A possible explanation of the higher variation in the DPF for these two subjects could be the variation in skull thickness, increasing toward the back of the head.
In this paper, we investigated the spatial relationships between optical and fMRI measures of the hemodynamic response to motor activity. By using Monte Carlo simulations of light propagation, we were able to calculate the sensitivity profiles for each optical source-detector pair on the DOI probe for each subject. Rather than performing image reconstruction of the optical data, we used these sensitivity profiles to forward project the fMRI data enabling us to predict the measurements between each optical channel. This method allows quantitative spatial comparison of optical and fMRI data without using the regularized inversion techniques standard to optical imaging.
Spatial Correlation between fMRI and DOI
In the first part of this paper, we compared diffuse optical measurements to simultaneously acquired BOLD-fMRI. The BOLD signal measures localized changes in deoxyhemoglobin (e.g., Buxton, 44 Obata 45). However, because MRI derives its signal from the hydrogen nuclei of water molecules, the BOLD signal reflects changes in the resonance lifetime of water rather than a direct measure of deoxyhemoglobin. The ferrous heme group in deoxyhemoglobin has a paramagnetic electron structure. This creates localized magnetic perturbations in and around blood vessels, which give rise to the shortening of the relaxation time of water molecules as they migrate through these regions. The change in relaxation times , in turn, gives rise to changes in signal intensity, which underlies the BOLD response. Because it derives from hydrogen nuclei, the BOLD signal contains contributions from both water molecules inside the blood vessels, which can interact directly with hemoglobin, and from water outside the blood vessels, which interact with the magnetic perturbations extending from the vessels into the parenchyma tissue.44, 45 The BOLD signal is thus an indirect measurement of HbR, depending on details of the fMRI acquisition sequence, baseline properties of the brain, and vascular sensitivities.52 Spatially, the BOLD signal should reflect the underlying vascular anatomy. It is generally believed that functional changes may be biased to larger venous structures.47, 53 This observation is the result of two factors. First, deoxyhemoglobin changes are believed to be larger in the veins than in the arterial and capillary compartments. Second, for gradient echo imaging, as used here, the extravascular signal is most sensitive to larger vessels, selectively biasing the BOLD signal to the draining pial veins.52 On the other hand, optical measures of HbR will similarly be biased to the venous areas because this is where these changes are the largest. Thus, one would expect that both the optically measured HbR and the BOLD signal would be biased distal to the center of neuronal activity.
The fMRI-ASL contrast is thought to arise mainly from blood flow changes in the arterial and capillary vessels, leading to a signal that is better localized to the cortex, as opposed to the larger vessel structures.47, 48 This result, along with the expectation of venous-weighted BOLD, helps explain the slight spatial discrepancies that have been noted between the BOLD and ASL modalities. 48, 49, 50, 51
In this experiment, as was expected based on the differential vascular compartments giving rise to each of the measurements, we found spatial agreement between HbR and BOLD and between the HbT and ASL profiles. In study I, we found that the amplitude of the BOLD signal was slightly better correlated to the optically measured HbR signal than to either the HbT or signals. In study II, when the complementary ASL and DOI comparisons were made, we found better correlation between the optically measured and ASL signals than between HbR and ASL.
Correlation of “Noise” between DOI and fMRI
Although the aim of this study was to compare the functional hemodynamic responses measured with DOI and fMRI, we observed that even some of the noise in the data was present in both the optical and fMRI time courses. This noise, which we consider to be the spurious early or anomalous undershoots in the time courses of the measurements, not typical of the traditional canonical responses of hemodynamic changes, can be seen in both data sets. While, conventionally, many of these features would be considered noise in the response and do not meet statistical significance in effects analysis, the clear presence in both modalities indicates that they most likely have true physiological underpinnings. While perhaps not a component of the functional response, these artifacts are indeed physiological and probably due to blood pressure or respiratory effects. Subject F from study II is an example of this feature showing both early dips in flow and volume and negative going response curves, which are present in both the measured optical and ASL data.
The correlation of these time courses between modalities and varying across regions of the probe indicates that this noise is of physiological origin. Such physiological fluctuations have been noted in previous literature in both fMRI and DOI studies (i.e., Toronov, 54 Bhattacharyya and Lowe,55 Wise 56). This suggests that accounting for physiological fluctuations may therefore improve future analysis of fMRI and/or optical data, because this is clear evidence of nonwhitened noise contributions to the signals. Because most popular statistical mapping methods used in fMRI or optical analysis often make assumptions of random, white, or signal noise, the presence of such strong contributions suggests more care should be taken to investigate these assumptions. Approaches to modeling these systemic signals have recently been described for analysis of optical data and fMRI and similar methods might be extended to combined multimodality analysis in the future. 57, 58, 59, 60
Furthermore, although optical imaging is very susceptible to superficial, systemic noise effects, such as cardiac or respiration signals, the presence of these correlated artifacts in the fMRI data suggests this noise is not entirely arising from the surface. Although, this superficial, systemic physiology clearly has an impact on the optical data, this alone cannot explain the correlated temporal noise between the fMRI and DOI data, because the fMRI signals are arising from mostly cortical tissues.
Generality of Results
Although this experimental design yielded evidence for the correlation of these two methodologies, we must be cautious in drawing conclusions as to the generality of these results to other experimental situations. For both DOI and fMRI, the observed hemodynamic parameters, such as blood volume, , or HbR, are spatially and temporally filtered versions of their true underlying changes in the brain.61, 62 In general, these filters are not expected to be identical between the DOI and fMRI. As a result, the measured results may differ because the filter functions may differ. Because these filters will depend on a number of factors, including the acquisition scheme employed, it is possible that other echo, scan parameters, or even brain regions will yield slightly differing results. At different magnetic field strengths, the extra- or intravascular weightings of the BOLD signal will vary, giving rise to images that are differentially biased to macrovascular structure.52, 63 Although this data clearly suggest that there is correspondence between the fMRI and DOI modalities, further investigation of the explicit details of these results should be pursued.
In this paper, we have shown that the ASL signal correlates with the spatial profile of blood volume changes in the brain. This observation is the result of the primarily arteriolar contributions to these two measurements. The correspondence of volume and flow changes will not be strictly obeyed and further investigation will be needed. In particular, the correspondence of ASL with optical measurements of blood flow, such as diffuse correlation spectroscopy,64, 65 should be explored in future work.
Implications toward Multimodality Integration
Understanding the spatial and temporal correspondence between fMRI and optical methods is an important prerequisite for the multimodality integration of these two techniques and a motivation of the current study. Incorporation of proper measurement models into the processing of multimodality data will enable the ability to exploit analysis schemes that combine the mutual information of all modalities appropriately. In the future, such approaches, designed to maximize the covariance between the two sets of measurements, might greatly improve our discrimination of functional contrast from measurement noise. Placed in a Bayesian formulation, calculating the hemodynamic response that maximizes the a posteriori estimate of the joint set of measurements, as opposed to treating the two modalities as independent measurements, could provide a more accurate estimate of the significance of the functional response.7, 66 Thus, analyzing these modalities together as measurements derived from a common physiological state variable, rather than as independent variables could allow us to better cancel instrument-associated noise factors. These combined analysis schemes are currently being pursued.
In addition, the spatial correspondence between optical and fMRI techniques may also have an impact on multimodality schemes for optical image reconstruction methods. Although DOI offers several advantages for studying the temporal dynamics of the hemodynamic response, as a tool for imaging or brain mapping, it has a number of limitations as an independent imaging tool. Because DOI is based on the measurement of back-reflected, diffuse light, optical imaging is restricted to superficial depths of 0.5 to into the brain.67 In addition, continuous-wave imaging systems, as used in this study, have virtually no ability to discern activation depth and are exponentially more sensitive to superficial activation. As a result, the varying activation depth across the region of interest creates partial-volume errors. This was the specific motivation for not using image reconstructions in this study.
With the addition of the high-spatial resolution anatomical and functional information that MRI can provide, the optical path-length and partial-volume corrections can be included to improve the quantitative accuracy of DOI (e.g., Intes, 5 Boas and Dale6). With functional MRI to provide prior information for the location of activation, image reconstruction schemes can make use of this to constrain solutions.7 Because this paper has directly shown that the fMRI data is indeed a possible solution to the inversion problem, we now have more confidence that these functionally constrained methods are valid. Although we did show that the spatial profiles of both modalities are consistent with one another, we found that the optical HbR more closely matched the BOLD data, whereas the more closely matched the ASL data. This suggests that prior spatial constraints used with these inversion methods should be used in a hemoglobin species-specific manner. While BOLD could provide a reliable spatial constraint on HbR, it does not correlate as well to the arterial-weighted or HbT changes. Likewise, ASL could provide a spatial prior to HbT, but might provide a less reliable constraint for HbR.
The synergistic information content provided by the DOI and fMRI methodologies opens up the possibility of future multimodality studies. It is clear that, although both methods measure aspects of the hemodynamic response, they do so using very different means. As a result, measurement sensitivity, data processing, and the very different physics of acquisition all can play important roles in interpretation of these results. Vital to the correct interpretation of such multimodality studies, proper control experiments must be first performed. By using the approach of projecting the fMRI signal through the calculated DOI sensitivity profile, we believe we can more quantitatively examine the spatial and temporal correlation between these measurements. In the first part of this study, the common measurement of HbR between the BOLD-fMRI and DOI recordings provided one such control. Our results here show a statistical temporal and spatial correlation between these measurements. In addition, in study II, we showed that optically measured and HbT were correlated with the ASL response both temporally ( , respectively) and spatially .
These results confirm the hypothesis that the optical and fMRI measurements arise from common hemodynamic physiological parameters. Furthermore, they suggest the measurements of , HbT, and ASL probe different vascular territories from HbR and the BOLD signal. These results have implications toward the future of the multimodality integration of DOI and fMRI.
We would like to thank Christiana Andre, Danny Joseph, Elizabeth Alf, and Monica Allen for helping in organizing and performing the experimental sessions. T. J. H. is funded by the Howard Hughes Medical Institute predoctorial fellowship program. This work was supported by the National Institutes of Health (Grants No. R01-EB002482, RO1-EB001954, and P41-RR14075).