combat.enigma (version 1.1)

lmm_fit: Fit LMM harmonization and obtain model coefficients

Description

Fit linear mixed-effects models (LMM) to test magnetic resonance imaging (MRI) data from different sites. Briefly, this function calls lme (from the nlme package) massively to fit LMM for each ROI/voxel/vertex with the variables in cov as fixed-effects factors and site as a random intercept. Importantly, it also addresses specific brain imaging details (e.g., prescaling fMRI measures if needed or accounting for varying collinearity due to differing brain coverage; see "Neural correlates of human fear conditioning and sources of variability: A mega-analysis and normative modeling study of fMRI data from 2,199 individuals").

Usage

lmm_fit(dat, site, cov = NULL, n.min = 10, impute_missing_cov = FALSE,
        prescaling = FALSE, lin.hyp = NULL, verbose = TRUE)

Value

A list of parameters plus the results of the LMM.

Arguments

dat

matrix or data.frame with the MRI data (e.g., ROIs, voxels, vertexes, etc.).

site

factor specifying the site of each individual.

cov

matrix or data.frame with the fixed-effects covariates.

n.min

(optional) number specifying the minimum size of a site to be analyzed.

impute_missing_cov

(optional) logical, whether to impute missing covariates.

prescaling

(optional) logical, whether to prescale the sites' data before conducting ComBat. See 'Details'.

lin.hyp

A hypothesis vector or matrix giving linear combinations of coefficients by rows as described in linearHypothesis (from the car package).

verbose

(optional) logical, whether to print some messages during execution.

Author

Joaquim Radua

Details

In several situations, the use of ComBat is problematic, such as when all participants of some sites have missing data for some ROIs or voxels. In these cases, a linear mixed effects model (LMM) may be preferable. On another note, setting prescaling equal to TRUE may be especially beneficial when the sites use different scales, as may be easily the case in fMRI mega-analyses.

References

Neural correlates of human fear conditioning and sources of variability: A mega-analysis and normative modeling study of fMRI data from 2,199 individuals, submitted.

See Also

combat_fit, prescale_fit, and linearHypothesis

Examples

raw_mri = combat_example[,6:19]
site = factor(combat_example$site)
mod = as.matrix(combat_example[,c("disorder", "age", "sex")])

# Estimate the effects of disorder with simple linear models
# WITHOUT harmonizing MRI data across sites
Out_raw = t(apply(raw_mri, 2, function (y) {
  m = summary(lm(y ~ mod, data = combat_example))
  coef(m)[2,c(1,3,4)]
}))

# Estimate the effects of disorder with simple linear models
# AFTER ComBat harmonizing MRI data across sites
combat = combat_fit(raw_mri, site, mod)
harmonized_mri = combat_apply(combat, raw_mri, site, mod)$dat.combat
Out_combat = t(apply(harmonized_mri, 2, function (y) {
  m = summary(lm(y ~ mod, data = combat_example))
  coef(m)[2,c(1,3,4)]
}))

# Estimate the effects of disorder with LMM harmonizing MRI data across sites
lmm = lmm_fit(raw_mri, site, mod)
Out_lmm = data.frame(
  b = lmm$b_map[[2]],
  t = lmm$t_map[[2]],
  p = 2 * pnorm(-abs(lmm$z_map[[2]]))
)
rownames(Out_lmm) = colnames(raw_mri)

# Results without harmonizing, with combat_fit and with lmm_fit:
cat("\nRaw results:\n")
print(round(Out_raw, 3))
cat("\nComBat results:\n")
print(round(Out_combat, 3))
cat("\nLMM results:\n")
print(round(Out_lmm, 3))

# Correlations between the three methods:
cat("\nCorrelation in coefficients:\n")
print(round(cor(cbind(Out_raw[,1], Out_combat[,1], Out_lmm[,1])), 3))
cat("\nCorrelation in p-values:\n")
print(round(cor(cbind(Out_raw[,3], Out_combat[,3], Out_lmm[,3])), 3))