Synopsis
Functional Magnetic Resonance Imaging (fMRI) provides a popular approach to imaging cerebral activation in response to stimuli. Reliably detecting activation is, however, not an easy proposition because only a very small proportions of voxels show true activation. These truly activated voxels are known to be spatially localized, yet incorporating this information is challenging to implement practically. We provide a model-based approach that incorporates spatial context in a practical and methodologically sound manner while postulating our a priori expectation that a certain proportion of voxels is truly active. Results on simulation experiments for different noise levels are uniformly encouraging. The methodology is also illustrated on a sports imagination experiment and shows its potential in making possible the adoption of fMRI as a clinical tool to detect awareness and improve treatment in individual patients in persistent vegetative state, such as traumatic brain injury survivors. Purpose
To improve activation detection accuracy in fMRI by a novel approach that explicitly incorporates the well-known view that no more than only 1-3% of voxels are expected to be a priori truly active in a typical fMRI experiment, and to do so in a methodologically sound but practical manner while also preservingspatial context in the context of detecting activated voxels. (Note that requiring exactly 1-3% activated voxels in the activation maps is not an accurate representation of our prior knowledge that 1-3% of voxels are activated on average and would increase the chance of missing pathologies, thus mis-diagnosing anomalies in a clinical setting.)
Background
One of the most important aspects of fMRI is the construction of voxel-wise maps to identify regions of neural activation. A common strategy fits a general linear model1 at each voxel to the time series of the BOLD response and uses Statistical Parametric Mapping2 to obtain a test statistic at each voxel and its p-value summarizing the association between the BOLD response and the stimulus, and its statistical significance, respectively. These p-values are thresholded to determine activation. Because of scanner and between-subject variability, inherent signal unreliability and,subject motion3-5, the several-seconds delay in the onset of the BOLD response and above all,by the fact that activation is expected to occur in no more than 1-3% of the voxels in the brain. Consequently,identified activation regions vary greatly from one dataset to another even when the same subject is scanned under the same paradigm.6-7
Methods
9
Imaging: The publicly-available8 dataset was from a study of a healthy female volunteer imaged while alternately resting and imagining playing tennis,blocks of 30 seconds each. AFNI9 was used to preprocess the data, after which analysis was done using the authors' R/C code.
Statistical:
Let $$$p=(p_1,p_2,\ldots,p_n)$$$ be the p-values obtained at each voxel. A general approach
10-11 models these
p-values as a mixture of beta distributions, with the first component as the uniform distribution and representing inactivated voxels. The remaining components follow the beta distribution with the first parameter substantially smaller than the second. Thus, each
pi is drawn from the mixture model $$f(p_i;\pi,\alpha,\beta)=\pi_0+\sum_{k=1}^K\pi_kb(p_i;\alpha_k,\beta_k)$$ where $$$b(\cdot;\alpha_k,\beta_k)$$$ is the beta density with shapeparameters $$$\alpha_k$$$ and $$$\beta_k$$$ with constraints that $$$(\pi_0,\pi_1,\ldots,\pi_K)$$$ sum to 1, and that $$$\pi_0\geq\delta$$$ where $$$1-\delta$$$ denotes the maximum proportion of voxels a priori expected to be activated. Thus, for given
K (denoting the different regions and strengths of activation), the loglikelihood equation is given by $$\sum_{i=1}^n\log[\pi_0+\sum_{k=1}^K\pi_kb(p_i;\alpha_k,\beta_k)]$$. To incorporate spatial context in the results, we add a penalty term involving the voxel coordinate to the loglikelihood as follows:$$$\sum_{i\neq j =1}^n\sum_{k=0}^K W_{i,k}W_{j,k}\left[-\frac{(v_{x,i} - v_{x,j})^2}{h_{x,k}^2}-\frac{(v_{y,i}-v_{y,j})^2}{h_{y,k}^2}-\frac{(v_{z,i} - v_{z,j})^2}{h_{z,k}^2}-\log h_{x,k}h_{y,k}h_{z,k}\right]$$$ where $$$(v_{x,i},v_{y,i},v_{z,i})^\prime$$$ are
ith voxel coordinates in the
x-, y- and
z-axes. The Expectation-Maximization algorithm
13 estimates the parameters, while ICL-BIC
14 estimates
K.
Results
Figure 1 displays the results of activation regons detected using our method -- here we have used a maximum true expected activation rate of 1%. For comparison purposes, we have also included results obtained using FDR thresholding (q<0.05) Figure 214 and cluster thresholding15 in Figure 3. The cluster thresholding was obtained using first-order nearest neighbors and minimum cluster sizes of 10 voxels (estimated via simulation by AFNI after accounting for the intrinsic smoothness level for this dataset and at an overall significance level of 0.05).Figure 1 shows detected activation in the pre-SMA, the SMA and the parietal cortex of the brain. A larger study with the same experimental paradigm over 14 normal subjects16 also reported similar findings (and, because of the larger-scale nature of this study can be considered to be a gold standard). Therefore, it is encouraging to note that we were able to get similar findings from our single-subject study. In comparison, however, FDR thresholding (Figure 2) was unable to detect activation in the SMA and very weak activation in the other regions. Clusterwise-thresholding (Figure 3) was able to detect faint activation in the SMA but stronger activation in the other regions of the brain.
Discussion
We have provided a novel model-based approach to image segmentation in fMRI Activation studies. Our approach sets the maximum proportion of activated voxels expected a priori to be activated (but importantly does not fix it to be so). Because of its encouraging performance in a single-subject experiment, it has the potential to be applied to detect pathologies and be adapted for use in a clinical setting, such as in the case of traumatic brain injury (TBI) survivors.
Acknowledgements
The research was supported in part by the National Institute of Biomedical Imaging and Bioengineering (NIBIB) of the National Institutes of Health (NIH) under its award no. R21EB016212. The content of this paper however is solely the responsibility of the authors and does not represent the official views of the NIH.References
1. Friston, K. J., Holmes, A. P., Worsley, K. J., Poline, J.-B., Frith, C. D., and Frackowiak, R. S. J. (1995), “Statistical parametricmaps in functional imaging: A general linear approach,” Human Brain Mapping, 2, 189–210.
2. Friston, K. J., Frith, C. D., Liddle, P. F., Dolan, R. J., Lammertsma, A. A., and Frackowiak, R. S. J. (1990), “The relationship between global and local changes in PET scans,” Journal of Cerebral Blood Flow and Metabolism, 10, 458–466.
3. Biswal, B., DeYoe, A. E., and Hyde, J. S. (1996), “Reduction of physiological fluctuations in fMRI using digital filters.” Magnetic Resonance in Medicine, 35, 107–113.
4. Genovese, C. R. (2000), “A Bayesian Time-Course Model for Functional Magnetic Resonance Imaging Data (with dis-cussion),” Journal of the American Statistical Association, 95, 691–719.
5. Hajnal, J. V., Myers, R., Oatridge, A., Schweiso, J. E., Young, J. R., and Bydder, G. M. (1994), “Artifacts due to stimulus-correlated motion in functional imaging of the brain,” Magnetic Resonance in Medicine, 31, 283–291.
6. Maitra, R., Roys, S. R., and Gullapalli, R. P. (2002), “Test-retest reliability estimation of functional MRI data,” Magnetic Resonance in Medicine, 48, 62–70.
7. Maitra, R. (2010), “A Re-defined and Generalized Percent-Overlap-of-Activation Measure for Studies of fMRI Reproducibilityand its Use in Identifying Outlier Activation Maps,” Neuroimage, 50, 124–135.
8. Tabelow, K. and Polzehl, J. (2011), “Statistical Parametric Maps for Functional MRI Experiments in R: The Package fmri,” Journal of Statistical Software, 44, 1–21.
9. Cox, R. W. (1996), “AFNI: software for analysis and visualization of functional magnetic resonance neuroimages,” Cox, R. W. (1996), “AFNI: software for analysis and visualization of functional magnetic resonance neuroimages,” Computers and biomedical research, an international journal, 29, 162–173.
10. Parker, R. A. and Rothenberg, R. B. (1988), “Identifying important results from multiple statistical tests,” Statistics in Medicine, 7, 1031–1043.
11. Allison, D. B., Gadbury, G. L., Heo, M., Ferna´ndez, J. R., Les, C.-K., Prolla, J. A., and Weindruch, R. (2002), “A mixture model approach for the analysis of microarray gene expression data,” Computational Statistics and Data Analysis, 39,1–20.
12. McLachlan, G. J. and Krishnan, T. (2008), The EM Algorithm and Extensions, Wiley.
13. Biernacki, C., Celeux, G., and Gold, E. M. (2000), “Assessing a mixture model for clustering with the integrated completed likelihood,” IEEE Transactions on Pattern Analysis and Machine Intelligence, 22, 719–725.
14. Genovese, C. R., Lazar, N. A., and Nichols, T. (2002), “Thresholding of statistical maps in functional neuroimagingusing the false discovery rate,” Neuroimage, 15, 870–878.
15. Heller, R., Stanley, D., Yekutieli, D., Rubin, N., and Benjamini, Y. (2006), “Cluster-based analysis of fMRI data,” Neuroimage, 33, 599–608.
16. Bardin, J. C., Fins, J. J., Katz, D. I., Hersh, J., Heier, L. A., Tabelow, K., Dyke, J. P., Ballon, D. J., Schiff, N. D., and Voss,H. U. (2011), “Dissociations between behavioural and functional magnetic resonance imaging-based evaluationsof cognitive function after brain injury,”
Brain,
134, 769–782.