Advanced Medical Image AnalysisView this Special Issue
Enhancing the Detection of BOLD Signal in fMRI by Reducing the Partial Volume Effect
Purpose. To investigate the advantages of reducing the partial volume effect (PVE) to enhance the detection of the BOLD signal in fMRI. Methods. A linear phase term was added in k-space to obtain half-voxel shifting of 64 × 64 -weighted echo-planar images. Three sets of image data shifted in the x, y, and diagonal direction, respectively, are combined with the original 64 × 64 data to form the 128 × 128 voxel-shifted interpolated data. Results. A simulation of a synthetic fMRI dataset shows that the voxel-shifted interpolation (VSI) can increase the t-score up to 50% in single-voxel activations. An fMRI study () demonstrates that 20.4% of the interpolated voxels have higher t-scores than their nearest neighboring voxels in the original maps. The average increase of the t-score in these interpolated voxels is 13.3%. Conclusion. VSI yields increased sensitivity in detecting voxel-size BOLD activations, improved spatial accuracy of activated regions, and improved detection of the peak BOLD signal of an activated region. VSI can potentially be used as an alternative to the high-resolution fMRI studies in which reduction in SNR and increase in imaging time become prohibitive.
Functional MRI (fMRI) studies with the blood oxygenation level-dependent (BOLD) contrast [1, 2] typically use a spatial resolution of a few millimeters (e.g., with an in-plane resolution of 2–4 mm and through-plane resolution of 3–6 mm). It is desirable to increase the spatial resolution in fMRI for at least two reasons. First, the spatial localization of neuronal activation can be improved with a higher spatial resolution [3, 4]. The gray matter (GM) typically has a thickness of 2-3 mm in human brains. Even within a cortical layer of GM, the distribution of venous vessels is heterogeneous and, therefore, the BOLD contrast can be increased by using a voxel size as small as 1 mm3 [5–7]. Secondly, the BOLD response originating from sub-voxel activation can be increased by reducing the voxel size. For example, it was demonstrated that BOLD activation in the hippocampal formation is often undetectable at a resolution of 4 × 4 × 5 mm3 but can be reliably detected at a higher resolution of 2 × 2 × 3 mm3 . In fMRI studies designed to examine midbrain nuclei , single-voxel BOLD activation is not unusual because of the small size of these nuclei.
The ability to detect small region (e.g., single-voxel size) BOLD activation is highly dependent on the location of the activation. The detected BOLD response is highest when the activation site is located at the center of a voxel and is lowest when the activation is located at the corners of 4 neighboring voxels in two-dimensional data (see Figure 1(a)) or 8 neighboring voxels in three-dimensional data. Similarly, the detected BOLD response at the edge of an activated region is dependent on the relative location of the edge to the voxels. This dependency of detected BOLD response on the relative location of the BOLD activation to the grid of voxels is determined by the voxel sensitivity function  and is referred to as the partial volume effect (PVE). PVE occurs in MRI when a voxel is only partially filled with the source of imaging contrast.
PVE can be reduced by acquiring images with high spatial resolution. High spatial resolution is commonly achieved by increasing the acquisition matrix of fMRI data, because any reduction of field-of-view (FOV) is limited by the dimension of the imaging objects or region of interest. Several studies have demonstrated the advantages of higher spatial resolution fMRI, including the increased spatial specificity of BOLD activations [4–6] and suppression of the BOLD signal originating from large vessels . The acquisition of high-resolution fMRI data, however, inevitably reduces spatial SNR and prolongs acquisition time. With the current echo-planar imaging (EPI) or spiral imaging technology, a low resolution (e.g., with a 64 × 64 matrix) whole-brain dataset can be acquired with a repetition time (TR) on the order of 2 seconds . Acquiring whole-brain fMRI with a higher spatial resolution (e.g., with a matrix size larger than 128 × 128) will considerably reduce the temporal resolution and spatial SNR, resulting in a greatly reduced statistical power in fMRI data analysis. This reduced statistical power can make the high-resolution approach unfavorable or even prohibitive to most fMRI studies. Furthermore, event-related fMRI studies require a TR of 2 seconds or shorter in order to adequately sample the hemodynamic response, making high-resolution acquisition even less useful for these studies. Additionally, image distortion induced by field inhomogeneity becomes more problematic in single-shot high-resolution echo-planar images .
In fMRI, PVE can be the result of three conditions: (1) when the size of the activated region is smaller than the voxel size, (2) when the peak of the true BOLD activation is mismatched from the location of any voxel, and (3) when the voxels at the boundary of activated regions are only partially occupied by the activation. When the dimension of an activation site is similar to the size of a voxel and the location of the activation is not centered on any one voxel (Figure 1(a)), the second and third conditions above can become the primary sources of PVE. In this case, the averaging of the BOLD signal with nonactivated partial volume greatly dilutes the BOLD contrast. When the BOLD contrast is diluted below certain threshold, the activated region becomes undetectable.
In this study, we hypothesize that spatial resampling with a half-voxel shift in the image space can greatly reduce PVE (Figure 1(b)) and enhance the ability to detect the activation through a better spatial match. In this paper, both simulation and fMRI experiments demonstrate that reducing PVE with voxel-shifted interpolation can substantially improve the detection of BOLD signal in fMRI.
2. Materials and Methods
2.1. Voxel-Shifted Interpolation
In the voxel-shift interpolation (VSI), three image sets with a 64 × 64 matrix were reconstructed with a half-voxel shift along the , , and diagonal direction, respectively, in addition to the original 64 × 64 image. According to the Fourier shift theorem, these shifts were implemented by applying a linear phase term in the complex -space data in the corresponding directions: where half-integer and is the matrix size of the original image. Finally, the voxels in the three shifted image matrices were interspersed with those in the original image, based on the relative shift that was applied to each of the shifted image, to form an interpolated image with a 128 × 128 matrix size (see Figure 2).
VSI can also be implemented by filling zeroes in the k-space and is thus also referred to as zero-filled interpolation  or by sinc interpolation with complex images. Although the names of VSI, zero-filled interpolation, and sinc interpolation suggest different approaches of implementation, they are nevertheless mathematically equivalent. We use the term VSI in this paper, as it provides a pictorial description of the interpolation that can better illustrate the reduction of PVE.
A synthetic dataset mimicking fMRI activation was generated for simulation. The data set had a matrix size of 128 × 128, with 240 time points, alternating “ON” states for 30 time points and “OFF” states for 30 time points. The background intensity was set at 10,000 in the real component and zero in the imaginary component. The data set had an SNR of 500 with Gaussian noise introduced in both real and imaginary components. The SNR of the synthetic data was intentionally set higher than that typically used in fMRI scans in order to better demonstrate the improvement in detecting “activation” through reducing PVE. In the “ON” state, a single-voxel “activation,” with a signal of 1%, 2.5%, and 5%, respectively, was placed at 4 locations. The 2D coordinates of these 4 locations are in the combination of odd, odd; even, odd; odd, even; and even, even, respectively, as shown in the upper part of Figure 3. The same “activated” voxels were duplicated at the lower part of the figure, with surrounding “activated” voxels being added to these single-voxel “activations.” The “activation” signal of each of the surrounding voxels was half of that in the central voxel. The simulation of the single-voxel “activation” was to mimic the activations in mid-brain nuclei, while the simulation of the multivoxel “activation” was to demonstrate the effect of PVE in cortical activations.
A two-dimensional Fourier transform was applied to the dataset to generate the -space data. The central portion of the -space data with a size of 64 × 64 was used to reconstruct the 64 × 64 original images. The same data set was also used to reconstruct the shifted images using VSI to form the 128 × 128 interpolated images. -scores of each voxel in the original and the interpolated images were calculated using the general linear model in SPM2b software package (Wellcome Department of Cognitive Neurology, University College London, UK). A -score of 3.18 was used to threshold the “activation” maps.
2.3. fMRI Experiments
Functional MRI data were acquired using -weighted gradient-echo EPI on a 3T MRI scanner (General Electric, Waukesha, WI), while subjects (, male/female = 3/4, age = years old) were performing a smooth pursuit eye movement (SPEM) task [14–16], under a protocol approved by a local institutional review board. The imaging parameters were TE/TR/flip angle = 32 ms/2500 ms/77°, 24 slices, field-of-view = 22 cm, slice thickness = 4 mm, and matrix size = 64 × 64. The task was run twice; each consisted of task (30 s)/rest (30 s) for 4 cycles.
Functional images were reconstructed offline from the -space data using a software package developed by our lab. The original image set with a 64 × 64 matrix was reconstructed directly from the acquired data. VSI was applied to the -space data to obtain the 128 × 128 interpolated images.
VSI was also applied to the original magnitude images for comparison. The magnitude VSI is an approximation of the -space VSI because the magnitude fMRI data are more accessible in scanners. The magnitude VSI can also be used to retrospectively interpolate the fMRI data acquired without saving the complex -space data. In this approach, an original 64 × 64 magnitude image was used as the real component of the complex image and the imaginary component of the complex image was set to zero. Fourier transform was applied to the complex images to generate the 64 × 64 -space data prior to VSI to a 128 × 128 matrix.
2.4. fMRI Data Analysis
fMRI data analysis was performed using the SPM2b software package. A -score was calculated for each individual voxel. The same procedure was used in the data analysis to ensure the comparability of results between the original images with a 64 × 64 matrix and the interpolated images. The procedure is not involved in any spatial interpolation or smoothing. In order to focus on the comparison between the interpolated and original fMRI data and avoid unnecessary variables in the comparison, no other preprocessing procedures, such as motion correction or normalization to standard space, were performed prior to the fMRI data analysis. These steps involve additional spatial interpolation and can further complicate the comparisons of interest in this study. The data from one of the subjects was discarded due to excessive motion during data acquisition. A -score of 3.18 () was used as the threshold of activation and the same scale of -score (i.e., from 3.18 to 22.83) was used in the display of the statistical parametric mapping results.
Simulation results show that the activated regions have different degrees of blurring in the original 64 × 64 images (see Figure 4(a)), depending on whether their - and -coordinates are even or odd numbers. The activations located at odd, odd coordinates (left column) have the least blurring and more likely have the highest -scores. The activations located at even, even coordinates (right column) have the most severe blurring and more likely have the lowest -scores. The activations located at even, odd coordinates are blurred horizontally and the activations located at odd, even coordinates are blurred vertically (see the middle two columns). These activations likely have intermediate -scores. In the 64 × 64 images shifted along -, -, and -directions, similar spatial dependency of the -scores on the relative location of the activations to the voxels is observed, except that the peak -scores are at the second column from the left in the -shifted image (Figure 4(b)), the third column from the left in the y-shifted image (Figure 4(c)), and the right column in the -shifted image (Figure 4(d)), respectively. In the interpolated 128 × 128 image (Figure 4(e)), this spatial dependency of -scores is largely removed. The residual variation of the -scores in different columns is caused by the added noise. As expected, all the activations were spatially blurred due to the 64 × 64 truncation of the -space data prior to the interpolation.
Figure 5 is the scatter plot demonstrating the increase of peak -scores in the interpolated 128 × 128 images (-axis) compared to the peak -score in the original and shifted 64 × 64 images (-axis) for each of the 24 activations as in Figure 3. Up to a 50% increase in the maximum t-score of the 1% single-voxel activation was observed (see the lower-left corner of the plot). Substantial increase of the peak -score can also occur with VSI in the single-voxel activations with higher BOLD signal and multivoxel activation with lower BOLD signal, as shown in the middle part of the plot. This increase of peak -score with VSI becomes modest in multivoxel activations with high BOLD signal, as shown in the upper-right corner of the plot. These simulation results indicate that VSI can substantially improve the detection of small activated regions with a low BOLD signal. This feature is quite relevant in real fMRI studies, because small activated regions usually have a lower BOLD signal and are more likely undetectable in fMRI analysis.
Human fMRI experiments also show substantial improvement using VSI to detect BOLD signal response. An example of the original (a) and interpolated (b) activation map at the visual cortex is shown in Figure 6. In this example, a single voxel activation in the original map can be a single-voxel (circle), a group of 6 voxels (short arrow), or a strip of 7 voxels (long arrow) in the interpolated map. Two apparently separated activated areas in the original map are connected in the interpolated map (ellipse). On the other hand, two similar 6-voxel activated regions in the interpolated map (long and medium arrows) appear as single-voxel and 3-voxel activated regions in the original map. In Figure 6(b), the voxels in the rows and columns indicated by the black arrows are obtained through VSI. The other voxels, one out of four in any 2 × 2 cluster, are from the original image.
The activation maps at the frontal eye fields (FEF) in 6 of the 7 subjects are shown in Figure 7. In this figure, the original maps are at the top, the maps interpolated with -space VSI are in the second row, and the maps interpolated with magnitude VSI are at the third row. In subject 1, two apparently connected activated voxels in the original map are separated in the interpolated map (large blue circle). A single voxel activation in the original maps can be a single voxel activation (small blue circle in subject 1) or a group of 6 activated voxels (small green circle in subject 2) in the interpolated maps. In subject 3, the activated voxels apparently form an enclosed loop in the original map (large green circle) but not in the interpolated map. In general, the activation with VSI confirms better to the known gyral anatomical location of the FEF acquired with a high-resolution SPEM study .
Figure 8 shows the evaluation of the -scores in the -space interpolated and original activation maps in the FEF region in six subjects. There were 304 voxels that had -scores above the threshold and thus considered as activated voxels. An interpolated voxel with a -score higher than any of its two nearest horizontal or vertical original neighbors (labeled as Os in the inset drawing) is referred to as an A-type voxel (labeled as A), while an interpolated voxel with a higher -score than any of its four nearest diagonal original neighbors (labeled as B) is referred to as a B-type voxel. Each circle in this plot shows an interpolated voxel that has a -score (-axis) higher than the highest t-score of its nearest original neighbors (-axis). In the 6 interpolated maps, 20.4% of the voxels (i.e., 62 out of 304 voxels) are either A- (47) or B-type (15) voxels, with an average of 13.3% increase in -score. In each of the 6 A-type interpolated voxels shown as dark solid circles, none of its nearest original neighbors has a t-score above the threshold of 3.18. These 6 interpolated voxels have -scores above the threshold and become part of the multivoxel activated regions. Single-voxel activations in the interpolated maps are excluded in the statistics. Similar evaluation of the -scores was applied to the magnitude interpolated activation maps (see Figure 9), with 311 activated voxels. Note that the signal in each voxel is slightly different between the k-space interpolated data and magnitude interpolated data. This evaluation demonstrates that 15.8% of the voxels (i.e., 49 out of 311 voxels) are either A- (38) or B-type (11) voxels, with an average of 14.4% increase of the -scores. The location of 76% A-type (31) and B-type voxels (6) in the magnitude interpolated maps are overlapped with that in the -space interpolated maps.
4.1. VSI versus Other Spatial Interpolation Approaches
The shifted and original images are independently reconstructed from the same -space data using the same algorithm. The only difference in the reconstruction of the shifted images is the added linear phase in -space prior to the inverse Fourier transform. The shifted images can be considered to carry the same amount of information and as “true” as the original images. Therefore, the interpolated fMRI data is not expected to increase false positive detection when voxel-by-voxel fMRI analysis is performed. Also, the shifted images have the same SNR and physical voxel dimension as the original images. In comparison, other spatial interpolation methods use the intensity of the neighboring voxels to generate the point data at locations not originally sampled. As such, the intensities of these spatially interpolated voxels are different than the true intensities if they were actually acquired. It is noteworthy that -space VSI does not introduce additional blurring, because it does not increase the physical voxel dimension, although the displayed interpolated images appear much less pixelated. In comparison, other interpolation methods usually introduce additional spatial blurring.
Complex -space VSI offers several advantages compared to the truncated sinc interpolation in the image space. First, the sinc interpolation applied to the magnitude image is only an approximation of the complex -space VSI, due to the lack of phase information presented in the image data. Secondly, the kernel size of the sinc interpolation is usually severely truncated in order to perform the interpolation with reasonable computation time. Because the equivalent kernel size of VSI is equal to the matrix size of the image, both complex -space VSI and magnitude VSI are superior to severely truncated sinc interpolation.
4.2. PVE and Its Reduction in BOLD fMRI
VSI provides a unique alternative for the detection of highly localized activation with a weak BOLD response. A voxel-size BOLD activation usually has a weak BOLD signal. For example, a previous study of the point spread function of the BOLD activation in human V1 at 4 Tesla demonstrated that both the spatial extent and the amplitude of BOLD signal decrease when the stimuli with a constant amplitude decrease in size . PVE can severely hamper the detection of a voxel-size BOLD activation, due to the small size and low BOLD signal change of the activation. Figure 4 shows that PVE can thus severely reduce the spatial accuracy in depicting the activated regions in the original map, especially at the boundaries of activated regions or when the size of the activated region is in the order of voxel size. Although PVE can be reduced by increasing the spatial resolution, the reduction of SNR in high-resolution fMRI can outweigh the benefits of reduced PVE or even become prohibitive in detecting the weak BOLD signal response in some studies.
The shape of the activated region is better depicted and the apparent location of activated region, which can be defined by the center of mass, can be more precisely determined with VSI. The increased precision in locating the activated regions can potentially improve the fMRI mapping of visual cortices , tonotopy [20–22], and somatotopy . VSI may be particularly useful in presurgical mapping studies, where the precision of anatomical location is of paramount importance.
In the simulation, the spatial dependency of the shape and maximum -score on the relative location of an activation site on the grid of voxels in the low resolution images (i.e., 64 × 64) clearly demonstrate the effect of PVE on the ability to detect BOLD responses in fMRI. In comparing the spatially dependent variation of the maximum -score for each of 24 activated regions in simulation, only the second condition of PVE discussed in the introduction (i.e., the peak of the BOLD activation is mismatched from the location of any voxel) was involved. The results of this simulation show that VSI is very effective in reducing PVE arising from this condition, especially when the activation is weak and highly localized.
The increased -score in the A- and B-type voxels in the human fMRI experiment is the result of reducing PVE through better matching the voxel sensitivity function with the spatial profile of the BOLD signal. In six A-type voxels shown as dark solid circles in Figure 8, the -scores of their nearest original neighbors are below the threshold. These A-type voxels, as indicated by the blue arrows in Figure 7, demonstrates the reduction of PVE arising from the third condition (i.e., the voxels at the boundary of activated regions are partially occupied by the activation) and improvement in delineating the boundary of activated regions.
In the magnitude interpolated maps, the A- and B-type voxels consist of 15.8% of the total activated voxels, comparable to 20.4% in -space interpolated maps. Furthermore, 76% of these A- and B-type voxels are at the same location as the A- or B-type voxels in the k-space interpolated maps. These comparisons suggest that the magnitude interpolation is a good approximation of the k-space interpolation and thus puts forward this VSI method as a viable postprocessing option for conventional fMRI image data.
4.3. Other Considerations
VSI also improves the delineation of GM and white matter (WM) in the functional images, as shown in Figure 10. The thickness of GM in human brains is on the order of 2 mm, smaller than the typical voxel size used in fMRI. The shape of the GM regions is less continuous in the original EPI data due to PVE or the lack of sampling density in the image space. After VSI, the GM regions in the EPI data became clearer and more continuous. This improved delineation between GM and WM can lead to (1) more accurate coregistration between the EPI volumes at different time points during motion correction, (2) more accurate coregistration between the EPI data and the high-resolution anatomical data, and (3) improved warping of the EPI data onto standard brain template, such as the Montreal Neurological Institute (MNI) template.
Further increasing the matrix size with an interpolation factor of four (e.g., from 64 × 64 to 256 × 256) or higher may slightly reduce PVE and improve the detection of fMRI signal. However, the additional improvement is expected to be much less substantial than the improvement obtained with an interpolation factor of two (e.g., from 64 × 64 to 128 × 128) . The penalty in computation time and the increased size of a dataset at a higher interpolation factors can easily outweigh the benefit of a more reduced PVE.
A direct comparison between the VSI activation maps and the ones obtained from high-resolution acquisition would be helpful for the validation of the proposed technique. Since high-resolution fMRI scans were not performed in the current study, such a comparison should be considered in future studies.
VSI improves the spatial accuracy of activation in fMRI through reducing PVE. Additionally, VSI can improve the sensitivity in detecting BOLD signal, especially when the size of activated regions is on the order of voxel size. VSI achieves these advantages of higher resolution fMRI without penalties in SNR, temporal resolution, and image artifacts.
Conflict of Interests
The authors declare that there is no conflict of interests regarding the publication of this paper.
This work was partially supported by the National Key Basic Research Program of China (2013CB329501), the National Natural Science Foundation of China (81371518, 81261120411), and the Ministry of Science and Technology (2011BAI12B01). The Brain and Behavior Research Foundation, and the Blowitz-Ridgeway Foundation.
P. A. Bandettini, E. C. Wong, R. S. Hinks, R. S. Tikofsky, and J. S. Hyde, “Time course EPI of human brain function during task activation,” Magnetic Resonance in Medicine, vol. 25, no. 2, pp. 390–397, 1992.View at: Google Scholar
R. S. Menon and B. G. Goodyear, “Submillimeter functional localization in human striate cortex using BOLD contrast at 4 Tesla: implications for the vascular point-spread function,” Magnetic Resonance in Medicine, vol. 41, no. 2, pp. 230–235, 1999.View at: Google Scholar
J. P. Szaflarski, S. K. Holland, V. J. Schmithorst, R. S. Dunn, and M. D. Privitera, “High-resolution functional MRI at 3 T in healthy and epilepsy subjects: Hippocampal activation with picture encoding task,” Epilepsy and Behavior, vol. 5, no. 2, pp. 244–252, 2004.View at: Publisher Site | Google Scholar
Y. P. Du, D. L. Parker, W. L. Davis, and G. Cao, “Reduction of partial-volume artifacts with zero-filled interpolation in three-dimensional MR angiography,” Journal of Magnetic Resonance Imaging, vol. 4, no. 5, pp. 733–741, 1994.View at: Google Scholar
B. Luna and J. A. Sweeney, “Cognitive functional magnetic resonance imaging at very-high-field: eye movement control,” Topics in Magnetic Resonance Imaging, vol. 10, no. 1, pp. 3–15, 1999.View at: Google Scholar
C. Rosano, C. M. Krisky, J. S. Welling et al., “Pursuit and saccadic eye movement subregions in human frontal eye field: a high-resolution fMRI investigation,” Cerebral Cortex, vol. 12, no. 2, pp. 107–115, 2002.View at: Google Scholar
K. Ueno, R. A. Waggoner, K. Tanaka, and K. Cheng, “Spatial precision of BOLD-fMRI in human V1: point spread function measured at 4T with spatially localized and size-varied stimuli,” in Proceedings of 11th Annual Meeting of the Organization for Human Brain Mapping (OHBM '05), p. 124, Toronto, Canada, 2005.View at: Google Scholar
M. I. Sereno, A. M. Dale, J. B. Reppas et al., “Borders of multiple visual areas in humans revealed by functional magnetic resonance imaging,” Science, vol. 268, no. 5212, pp. 889–893, 1995.View at: Google Scholar
P. Hluštík, A. Solodkin, R. P. Gullapalli, D. C. Noll, and S. L. Small, “Somatotopy in human primary motor and somatosensory hand representations revisited,” Cerebral Cortex, vol. 11, no. 4, pp. 312–321, 2001.View at: Google Scholar