Introduction:

Task-based fMRI has been widely used to determine brain activations during cognitive tasks of different populations. Popular group-level analysis is based on the univariate general linear model (GLM)^1,2. However, univariate methods are known to suffer from low sensitivity for a given specificity because the spatial covariance structure at each voxel is not taken entirely into account³. At the subject-level analysis, multivariate canonical correlation analysis applied to a small neighborhood has been shown to improve performance over univariate methods^{4, 5}. In this study, a spatially constrained multivariate model is introduced for group-level analysis to improve sensitivity at a given specificity for activation detection.

Methods:

Model: Subject-level GLM model can be written as y=Xβ+ε, where vector y represents a time course from a single voxel, X represents functions used to model the BOLD response, and β is the subject effect. Univariate group analysis is modeled as β _u=X_Gβ_Gu+η_Gu where β_u represents effects from N subjects at a single voxel and X_G(Nx1) is the design matrix modeling within-group activation detection. Multivariate group analysis is then modeled as a multivariate multiple regression problem according to Bα_G=X_Gβ_G+η_G , where B=(β⁽¹⁾_u,...,β^(m)_u) is a matrix that represents 1^st level effects of m voxels in a local neighborhood from N subjects, α_G is a spatial weight vector and β_G is the combined group effect. A previous study⁶ suggests that β^(k)_u of the k^th voxel from all N subjects follows a multivariate Gaussian distribution with mean X_Gβ_Gu and covariance V_Gu. As stated above, Bα_G in the multivariate model is a linear combination of β^(k)_u in a local neighborhood; therefore, Bα_G also follows a multivariate Gaussian distribution with mean X_Gβ_G and covariance V_G. We estimate the group inference β_G, spatial weight vector α_G, and the between group variance term in V_G by maximizing the log-likelihood of the probability density function of the multivariate Gaussian distribution. Statistics of β_G is finally computed and assigned to the center voxel. A constraint is further added to the model to guarantee the dominance of the center voxel in a local neighborhood so that false positives are limited. Imaging: fMRI data (TR/TE/resolution=2s/30ms/1.7x1.7x5mm³, 25 slices, coronal oblique, 288 time frames) from 16 subjects each consisting of a resting-state data set and an episodic memory task data set were analyzed. The memory task involved viewing faces paired with occupations and contained instruction, encoding, recognition and control periods. Validation with simulated data: 500 5x5 neighborhoods with active center voxels and 500 neighborhoods with inactive center voxels for 16 subjects were simulated. The distribution of active neighbors in each local neighborhood followed the empirical distribution of real fMRI data analyzed with the univariate method. Time-series for simulated neighborhoods were obtained from neighborhoods in real data with the same activation patterns (both resting-state and task fMRI). Specifically, wavelet-resampled resting-state time-series were added to the task time-series with different noise fractions to simulate time courses at different noise levels. Both subject-level and group-level univariate and multivariate techniques were applied to the simulated data and the receiver operating characteristic (ROC) method was used to evaluate the performance of each method. Validation with real data: The standard mass-univariate method was applied in fMRI 1^st level analysis and a voxel-wise effect map of contrast encoding v/s control was computed for each subject. Both univariate and multivariate techniques were applied to the 2^nd level analysis. Statistical significance was determined from null distributions, which were generated by repeating the exact analysis on wavelet-resampled resting-state time series at least 200 times.

Results:

Table 1 summarizes optimization models used in the univariate⁶ and the proposed multivariate group-analysis methods. Table 2 lists the methods applied and compared using simulated data, with the details of subject-level smoothness and both subject-level and group-level analysis methods. Simulation results are shown in Fig.1. Areas under the ROC curves (AUC) at multiple noise levels demonstrate optimum performance of the proposed multivariate methods, especially at high noise levels. Applying both levels of multivariate analysis methods will further increase the AUC by 5% at noise levels of 0.85. Whole brain within-group activation maps for real fMRI data with contrast encoding v/s control are shown in Fig.2 (p<0.05, corrected) and improved activations are observed in fusiform-gyrus and para-hippocampal areas.

Discussion and Conclusion:

We have introduced a constrained multivariate method to incorporate local neighboring voxels in an fMRI group-level analysis. Using simulation, we have demonstrated better performance in activation detection of the proposed method over univariate techniques. Applying the proposed method to real fMRI episodic memory data, larger within-group activations in fusiform gyrus and para-hippocampal areas were observed.

Acknowledgements

The study is supported by the National Institutes of Health (grant number 1R01EB014284 and P20GM109025).