R/methylationGLMM_T1T2.R
methylationGLMM_T1T2.RdmethylationGLMM_T1T2() is the high-level coordinator for the longitudinal
linear mixed-effects stage of the dnaEPICO workflow. It prepares the merged
phenotype-plus-beta input, fits one mixed-effects model per CpG for each
requested phenotype, extracts phenotype-specific coefficient summaries,
optionally collects significant interaction tables, generates diagnostic plots,
annotates the combined summary table, and optionally writes legacy-style
outputs to disk. The default behavior is now in-memory and quiet, which makes
the function easier to compose with other package functions and more aligned
with typical Bioconductor usage.
methylationGLMM_T1T2(
inputPheno = "rData/preprocessingPheno/mergeData/phenoBetaT1T2.RData",
outputLogs = "logs",
outputRData = "rData/methylationGLMM_T1T2/models",
outputPlots = "figures/methylationGLMM_T1T2",
personVar = "person",
timeVar = "Timepoint",
phenotypes = c("DASS_Depression", "DASS_Anxiety", "DASS_Stress", "PCL5_TotalScore",
"MHCSF_TotalScore", "BRS_TotalScore"),
covariates = "Sex,Age,Ethnicity,TraumaDefinition,Leukocytes,Epithelial.cells",
factorVars = "Sex,Ethnicity,TraumaDefinition,Timepoint",
lmeLibs = "lme4,lmerTest",
prsMap = NULL,
libPath = NULL,
cpgPrefix = "cg",
cpgLimit = NA,
nCores = 32,
summaryPval = NA,
plotWidth = 2000,
plotHeight = 1000,
plotDPI = 150,
interactionTerm = NULL,
saveSignificantInteractions = TRUE,
significantInteractionDir = "preliminaryResults/cpgs/methylationGLMM_T1T2",
significantInteractionPval = 0.05,
saveTxtSummaries = TRUE,
chunkSize = NULL,
summaryTxtDir = "preliminaryResults/summary/methylationGLMM_T1T2/lmer",
fdrThreshold = 0.05,
padjmethod = "fdr",
annotationPackage = "IlluminaHumanMethylationEPICv2anno.20a1.hg38",
annotationCols = c("Name", "chr", "pos", "UCSC_RefGene_Group", "UCSC_RefGene_Name",
"Relation_to_Island", "GencodeV41_Group"),
annotatedLMEOut = "data/methylationGLMM_T1T2",
display = FALSE,
verbose = FALSE,
logs = FALSE,
saveOutputs = FALSE
)Character. Path to the merged longitudinal phenotype-plus-beta
.RData or .rds object created by preprocessingPheno(). The default
points to the combined timepoint object produced by the package workflow.
Character. Directory used for optional log files.
Character. Directory used for optional serialized mixed-model and summary outputs.
Character. Directory used for optional TIFF diagnostic plots.
Character. Subject identifier variable used for the random
intercept. When this column is missing, it is derived from SID using the
package's existing sample naming convention.
Character. Name of the longitudinal time variable included as a fixed effect in every model.
Character vector or comma-separated phenotype variables to model.
Character. Comma-separated fixed-effect covariates included in every mixed model.
Character. Comma-separated variables that should be coerced
to factors before modeling. This usually includes categorical covariates and
timeVar.
Character. Comma-separated package names to validate on worker
processes. The default is "lme4,lmerTest".
Character or NULL. Optional phenotype-to-PRS mapping in the
form "Phenotype1:PRS_1,Phenotype2:PRS_2".
Character vector or NULL. Optional library paths forwarded to
worker processes. By default, the current .libPaths() are used.
Character. Prefix used to identify methylation columns in the
merged phenotype-plus-beta input object. The default is "cg".
Integer or NA. Maximum number of CpGs to analyse. Use NA
to keep all CpGs matching cpgPrefix.
Integer. Number of worker processes to use while fitting models and extracting summaries.
Numeric or NA. Optional p-value threshold applied to the
returned longitudinal CpG summary tables. Use NA to keep all summary rows.
Integer. TIFF width in pixels when plots are written to disk.
Integer. TIFF height in pixels when plots are written to disk.
Integer. TIFF resolution in DPI when plots are written to disk.
Character or NULL. Optional interaction term. When
supplied and present in the input data, the phenotype is modeled together
with its interaction against this variable.
Logical. If TRUE, collect coefficient
tables for CpGs passing significantInteractionPval in the returned object
and optionally write them to disk when saveOutputs = TRUE.
Character. Directory used for optional significant-interaction coefficient tables.
Numeric. P-value threshold used to collect or write significant interaction coefficient tables.
Logical. If TRUE and saveOutputs = TRUE, write
tab-delimited summary tables to summaryTxtDir.
Integer or NULL. Number of CpGs processed per summary
extraction chunk. NULL chooses a value automatically.
Character. Directory used for optional tab-delimited LME summary tables.
Numeric. False-discovery-rate threshold used to highlight CpGs in the residual-significance diagnostic plots.
Character. P-value adjustment method passed to
stats::p.adjust(). The default is "fdr".
Character. Annotation package or object name passed to
minfi::getAnnotation(), for example
"IlluminaHumanMethylationEPICv2anno.20a1.hg38".
Character vector or comma-separated annotation columns to append to the combined LME summary table. Available columns depend on the selected annotation package.
Character. Directory used for the optional annotated LME summary CSV file.
Logical. If TRUE, draw diagnostic plots on the active
graphics device.
Logical. If TRUE, emit progress messages with message(). The
default is FALSE, so the function is quiet unless requested.
Logical. If TRUE, write the same progress messages to
file.path(outputLogs, "log_methylationGLMM_T1T2.txt").
Logical. If TRUE, write optional serialized model files,
summary tables, significant-interaction tables, annotated results, and TIFF
plots to the requested output directories. The default is FALSE, so the
function returns in-memory results without writing files.
A list with class "dnaEPICO_methylationGLMM_T1T2".
Object returned by prepareMethylationGLMM_T1T2Data()
containing the merged longitudinal phenotype-plus-beta analysis table and
modeling metadata.
Object returned by fitMethylationGLMM_T1T2Models()
containing the per-phenotype CpG mixed-effects model fits.
Object returned by
summarizeMethylationGLMM_T1T2Models() containing the combined CpG summary
tables used for reporting and annotation.
Object returned by
collectSignificantInteractionsMethylationGLMM_T1T2() containing optional
phenotype-specific significant-interaction tables.
Object returned by
plotMethylationGLMM_T1T2Diagnostics() describing the diagnostic plot
objects and any written TIFF files.
Object returned by
annotateMethylationGLMM_T1T2Summaries() containing the annotated combined
summary table.
Object returned by
writeMethylationGLMM_T1T2Outputs() when saveOutputs = TRUE, otherwise
NULL.
See dnaEPICO_methylationGLMM_T1T2 for a class-level overview.
if (
requireNamespace("IlluminaHumanMethylation450kanno.ilmn12.hg19", quietly = TRUE) &&
requireNamespace("lmerTest", quietly = TRUE)
) {
tmp <- tempdir()
toy_path <- file.path(tmp, "phenoBetaT1T2.RData")
phenoBT1T2 <- data.frame(
SID = c("P1A", "P1B", "P2A", "P2B", "P3A", "P3B", "P4A", "P4B"),
person = c(1, 1, 2, 2, 3, 3, 4, 4),
Timepoint = factor(c("1", "2", "1", "2", "1", "2", "1", "2")),
score = c(10, 12, 9, 11, 13, 14, 8, 9),
sex = factor(c("F", "F", "M", "M", "F", "F", "M", "M")),
cg00000029 = c(0.25, 0.27, 0.20, 0.22, 0.30, 0.31, 0.18, 0.20),
cg00000108 = c(0.50, 0.53, 0.55, 0.57, 0.48, 0.49, 0.60, 0.61),
check.names = FALSE
)
save(phenoBT1T2, file = toy_path)
result <- methylationGLMM_T1T2(
inputPheno = toy_path,
phenotypes = "score",
covariates = "sex",
factorVars = "sex,Timepoint",
cpgLimit = 2,
nCores = 1,
summaryPval = 1,
annotationPackage = "IlluminaHumanMethylation450kanno.ilmn12.hg19",
annotationCols = "Name,chr,pos",
display = FALSE,
verbose = FALSE,
logs = FALSE,
saveOutputs = FALSE
)
class(result)
}
#> [1] "dnaEPICO_methylationGLMM_T1T2"