## ----setup, include = FALSE---------------------------------------------------
knitr::opts_chunk$set(
    collapse = TRUE,
    comment = "##"
)

## ----eval = FALSE-------------------------------------------------------------
# if (!requireNamespace("BiocManager", quietly = TRUE)) {
#     install.packages("BiocManager")
# }
# BiocManager::install("batchCorr")

## ----message=FALSE, warning=FALSE---------------------------------------------
library(batchCorr)
data("ThreeBatchData")

## -----------------------------------------------------------------------------
peaks <- SimpleList(t(PTnofill), t(PTfill))
sampleData <- meta
featureData <- peakInfo(PT = PTnofill, sep = "@", start = 3)
rownames(featureData) <- rownames(peaks[[1]])

se <- SummarizedExperiment(assays = peaks, colData = sampleData,
    rowData = featureData)
names(assays(se)) <- c("nofill", "fill")

## -----------------------------------------------------------------------------
# Extract peakinfo (i.e. m/z and rt of features),
# These column names have 2 leading characters describing LC-MS mode
# -> start at 3
peakIn <- peakInfo(PT = PTnofill, sep = "@", start = 3)
# Perform multi-batch alignment
alignBat <- alignBatches(
    peakInfo = peakIn, PeakTabNoFill = PTnofill,
    PeakTabFilled = PTfill, batches = meta$batch,
    sampleGroups = meta$grp, selectGroup = "QC",
    report = FALSE
)
# Extract new peak table
PT <- alignBat$PTalign

## -----------------------------------------------------------------------------
se <- alignBatches(PeakTabNoFill = se, PeakTabFilled = se, batches = "batch",
    sampleGroups = "grp", report = FALSE, assay.type1 = "nofill",
    assay.type2 = "fill", name = "aligned", rt_col = "rt", mz_col = "mz")

## ----message = FALSE----------------------------------------------------------
# Batch B
batchB <- getBatch(
    peakTable = PT, meta = meta,
    batch = meta$batch, select = "B"
)
BCorr <- correctDrift(
    peakTable = batchB$peakTable,
    injections = batchB$meta$inj,
    sampleGroups = batchB$meta$grp, QCID = "QC",
    G = seq(5, 35, by = 3), modelNames = c("VVE", "VEE"),
    report = FALSE
)
# Batch F
batchF <- getBatch(
    peakTable = PT, meta = meta,
    batch = meta$batch, select = "F"
)
FCorr <- correctDrift(
    peakTable = batchF$peakTable,
    injections = batchF$meta$inj,
    sampleGroups = batchF$meta$grp,
    QCID = "QC", G = seq(5, 35, by = 3),
    modelNames = c("VVE", "VEE"),
    report = FALSE
)
# Batch H
batchH <- getBatch(
    peakTable = PT, meta = meta,
    batch = meta$batch, select = "H"
)
HCorr <- correctDrift(
    peakTable = batchH$peakTable,
    injections = batchH$meta$inj,
    sampleGroups = batchH$meta$grp,
    QCID = "QC", G = seq(5, 35, by = 3),
    modelNames = c("VVE", "VEE"),
    report = FALSE
)

HCorr$BIC

## ----message = FALSE----------------------------------------------------------
batch_labels <- unique(colData(se)$batch)
batches <- lapply(batch_labels, function(batch_label) {
    se[, colData(se)$batch == batch_label]
})

batches <- lapply(batches, correctDrift, injections = "inj", 
    sampleGroups = "grp", RefID = "Ref", G = seq(5, 35, by = 3), 
    modelNames = c("VVE", "VEE"), report = FALSE, 
    assay.type = "aligned", name = "corrected")

## ----warning = FALSE----------------------------------------------------------
mergedData <- mergeBatches(list(BCorr, FCorr, HCorr), qualRatio = 0.5)
normData <- normalizeBatches(
    peakTableCorr = mergedData$peakTableCorr,
    batches = meta$batch, sampleGroup = meta$grp,
    refGroup = "Ref", population = "all"
)
PTnorm <- normData$peakTable

## ----results = FALSE----------------------------------------------------------
se <- do.call(cbind, batches)
se <- se[which(apply(assay(se), 1, function(x) any(is.na(x)))), ]

se <- normalizeBatches(peakTableCorr = se, 
    batches = "batch", sampleGroup = "grp", refGroup = "Ref",
    population = "all", assay.type = "corrected", name = "normalized")

## ----echo = FALSE-------------------------------------------------------------
sessionInfo()