## ----include = FALSE----------------------------------------------------------
startTime <- Sys.time()
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>",
  crop = NULL
)

library(dnaEPICO)

can_run_preprocessing <- requireNamespace("minfiData", quietly = TRUE) &&
  requireNamespace("IlluminaHumanMethylation450kmanifest", quietly = TRUE) &&
  requireNamespace("IlluminaHumanMethylation450kanno.ilmn12.hg19", quietly = TRUE)

can_run_sva <- requireNamespace("minfiData", quietly = TRUE)

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

## ----eval = FALSE, message = FALSE--------------------------------------------
# library("dnaEPICO")

## ----preprocessing-local-inputs, eval = can_run_preprocessing, message = FALSE, warning = FALSE----
preprocessing_inputs <- dnaEPICO:::exampleMinfiIdatInputsDnaEpico(n = 6L)

names(preprocessing_inputs)
preprocessing_inputs$tempDir
basename(preprocessing_inputs$phenoFile)
head(basename(list.files(preprocessing_inputs$idatFolder, full.names = TRUE)), 4)
basename(preprocessing_inputs$probeExclusionPath)

## ----preprocessing-local, eval = can_run_preprocessing, message = FALSE, warning = FALSE----
preprocessResult <- preprocessingMinfiEwasWater(
  phenoFile = preprocessing_inputs$phenoFile,
  idatFolder = preprocessing_inputs$idatFolder,
  outputLogs = file.path(preprocessing_inputs$tempDir, "logs"),
  nSamples = 6,
  SampleID = "Sample_Name",
  arrayType = preprocessing_inputs$arrayType,
  annotationVersion = preprocessing_inputs$annotationVersion,
  scriptLabel = "preprocessingMinfiEwasWater",
  baseDataFolder = file.path(preprocessing_inputs$tempDir, "rData"),
  figureBaseDir = file.path(preprocessing_inputs$tempDir, "figures"),
  detPThreshold = 1,
  normMethods = "quantile",
  sexColumn = "Sex",
  pvalThreshold = 1,
  chrToRemove = "",
  snpsToRemove = "SBE",
  mafThreshold = 1,
  probeExclusionPath = preprocessing_inputs$probeExclusionPath,
  plotGroupVar = "Sex",
  lcRef = "saliva",
  phenoOrder = "Sample_Name;Sex;Basename;Sentrix_ID;Sentrix_Position",
  lcPhenoDir = preprocessing_inputs$tempDir,
  display = FALSE,
  verbose = FALSE,
  logs = FALSE,
  saveOutputs = FALSE
)

class(preprocessResult)
names(preprocessResult)
head(preprocessResult$targets[, c(
  "Sample_Name",
  "Sex",
  "Sentrix_ID",
  "Sentrix_Position"
)])

## ----preprocessing-local-inspect, eval = can_run_preprocessing, message = FALSE, warning = FALSE----
class(preprocessResult$RGSet)
names(preprocessResult$metricsData)
head(preprocessResult$lcData$phenoLC[, 1:6])

## ----sva-local-inputs, eval = can_run_sva, message = FALSE, warning = FALSE----
sva_inputs <- dnaEPICO:::exampleMinfiBaseDataDnaEpico()
sva_temp_dir <- tempdir()
sva_targets_file <- file.path(sva_temp_dir, "sva_targets.csv")
sva_rgset_file <- file.path(sva_temp_dir, "sva_RGSet.RData")

names(sva_inputs)
class(sva_inputs$RGSet)
utils::write.csv(sva_inputs$targets, sva_targets_file, row.names = FALSE)
RGSet <- sva_inputs$RGSet
save(RGSet, file = sva_rgset_file)

basename(c(sva_targets_file, sva_rgset_file))
file.exists(c(sva_targets_file, sva_rgset_file))

## ----sva-local, eval = can_run_sva, message = FALSE, warning = FALSE----------
svaResult <- svaEnmix(
  phenoFile = sva_targets_file,
  rgsetData = sva_rgset_file,
  outputLogs = file.path(tempdir(), "logs"),
  SampleID = "Sample_Name",
  arrayType = "IlluminaHumanMethylation450k",
  annotationVersion = "ilmn12.hg19",
  SentrixIDColumn = "Sentrix_ID",
  SentrixPositionColumn = "Sentrix_Position",
  ctrlSvaPercVar = 0.90,
  ctrlSvaFlag = 1,
  scriptLabel = "svaEnmix",
  display = FALSE,
  verbose = FALSE,
  logs = FALSE,
  saveOutputs = FALSE
)

class(svaResult)
names(svaResult)
dim(svaResult$svaData$sva)
head(svaResult$mergedPheno[, 1:min(6, ncol(svaResult$mergedPheno))])

## ----sva-local-inspect, eval = can_run_sva, message = FALSE, warning = FALSE----
class(svaResult$analysisData)
names(svaResult$analysisData)

## ----preprocessingPheno-local-inputs, message = FALSE, warning = FALSE--------
pheno_inputs <- dnaEPICO:::examplePreprocessingPhenoStateDnaEpico()

names(pheno_inputs)
pheno_inputs$tempDir
basename(c(
  pheno_inputs$phenoPath,
  pheno_inputs$betaPath,
  pheno_inputs$mPath,
  pheno_inputs$cnPath
))

## ----preprocessingPheno-local, message = FALSE, warning = FALSE---------------
phenoResult <- preprocessingPheno(
  phenoFile = pheno_inputs$phenoPath,
  betaPath = pheno_inputs$betaPath,
  mPath = pheno_inputs$mPath,
  cnPath = pheno_inputs$cnPath,
  SampleID = "Sample_Name",
  timeVar = "Timepoint",
  timepoints = "1,2",
  combineTimepoints = "1,2",
  outputPheno = file.path(pheno_inputs$tempDir, "data", "preprocessingPheno"),
  outputRData = file.path(
    pheno_inputs$tempDir,
    "rData",
    "preprocessingPheno",
    "metrics"
  ),
  outputRDataMerge = file.path(
    pheno_inputs$tempDir,
    "rData",
    "preprocessingPheno",
    "mergeData"
  ),
  sexColumn = "Sex",
  outputLogs = file.path(pheno_inputs$tempDir, "logs"),
  outputDir = file.path(pheno_inputs$tempDir, "clockFoundation"),
  verbose = FALSE,
  logs = FALSE,
  saveOutputs = FALSE
)

class(phenoResult)
names(phenoResult)
names(phenoResult$combinedData)
head(phenoResult$combinedData$phenoMethylation[, 1:6])

## ----preprocessingPheno-local-inspect, message = FALSE, warning = FALSE-------
names(phenoResult$clockFoundation)
head(phenoResult$clockFoundation$phenoCF)

## ----echo = FALSE-------------------------------------------------------------
Sys.time()

## ----echo = FALSE-------------------------------------------------------------
totalTime <- diff(c(startTime, Sys.time()))
round(totalTime, digits = 3)

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

