Merge MSBB Data

This script uses ComBat to merge MSBB data.

Dependencies

Load requisite packages.

library(data.table)
library(sva)

Read Data

Define brain regions.

brain_regions = c("IFG", "STG", "PHG", "FP")

Read expression data.

# read mData
MSBB_mData = lapply(brain_regions, function(brain_region){
  # process expression file
  cat("processing", brain_region, "...\n")
  load(paste0("../Data/MSBB-", brain_region,"-24.Rdata"))
  return(expSet$mData)
})
all(unlist(lapply(MSBB_mData, function(x) all(rownames(x) == rownames(MSBB_mData[[1]])))))
MSBB_mData = do.call(cbind, MSBB_mData)

Read covariate data.

# load cov data
MSBB_cov = lapply(brain_regions, function(brain_region){
  # process expression file
  cat("processing", brain_region, "...\n")
  load(paste0("../Data/MSBB-", brain_region,"-24.Rdata"))
  expSet$cov$brain_region = brain_region
  return(expSet$cov)
})
MSBB_cov = do.call(rbind, MSBB_cov)

Confirm sample order.

all(colnames(MSBB_mData) == MSBB_cov$individualIdentifier)

Remove Batch Effects

First, create the model matrix for the adjustment variables, including the variable of interest. Note that you do not include batch when creating this model matrix; rather, batch will be included in the ComBat function call in the following chunk. In this case, there are no other adjustment variables so we simply fit an intercept term.

modcombat = model.matrix(~1, data=MSBB_cov)

Now, remove batch effects using the ComBat approach from the sva package. The input data is assumed to be cleaned and normalized before batch effect removal.

MSBB_combat = ComBat(dat=MSBB_mData, batch=MSBB_cov$brain_region, mod=modcombat)

Save Results

expSet = list(mData = MSBB_combat,
              cov = MSBB_cov)
save(expSet, file = "../Data/MSBB-combated.Rdata")