## -----------------------------------------------------------------------------
#| label: setup
#| include: false
knitr::opts_chunk$set(tidy = FALSE, cache = TRUE, collapse = TRUE,
                      comment = "#>")
options(width = 80)


## -----------------------------------------------------------------------------
#| label: loadLibraries
suppressPackageStartupMessages({
  library(methylumi)
  library(TCGAMethylation450k)
  library(FDb.InfiniumMethylation.hg19)
})


## -----------------------------------------------------------------------------
#| label: loadData
## read in 10 BRCA IDATs
idatPath <- system.file("extdata/idat", package = "TCGAMethylation450k")
mset450k <- methylumIDAT(getBarcodes(path = idatPath), idatPath = idatPath)
sampleNames(mset450k) <- paste0("TCGA", seq_along(sampleNames(mset450k)))
show(mset450k)


## -----------------------------------------------------------------------------
#| label: fig-controls
#| fig-width: 5
#| fig-height: 7
#| cache: false
#| fig-cap: "Some of the controls on the 450k chip."
library(ggplot2)
## for larger datasets, the by.type argument can be set to FALSE;
## positional effects will manifest as a wave-like pattern
p <- qc.probe.plot(mset450k, by.type = TRUE)
print(p)


## -----------------------------------------------------------------------------
#| label: preprocess
mset450k.proc <- stripOOB(normalizeMethyLumiSet(methylumi.bgcorr(mset450k)))


## -----------------------------------------------------------------------------
#| label: fig-controls2
#| fig-width: 5
#| fig-height: 7
#| cache: false
#| fig-cap: "Controls after preprocessing."
p2 <- qc.probe.plot(mset450k.proc, by.type = TRUE)
print(p2)


## -----------------------------------------------------------------------------
#| label: coerceLumi
suppressPackageStartupMessages(library(lumi))
mset450k.lumi <- as(mset450k.proc, "MethyLumiM")
show(mset450k.lumi)


## -----------------------------------------------------------------------------
#| label: coerceBack
mset450k.andBack <- as(mset450k.lumi, "MethyLumiSet")
show(mset450k.andBack)


## -----------------------------------------------------------------------------
#| label: coerceMinfi
suppressPackageStartupMessages(library(FDb.InfiniumMethylation.hg19))
rgSet450k <- as(mset450k, "RGChannelSet")
show(rgSet450k)


## -----------------------------------------------------------------------------
#| label: coerceMinfi2
suppressPackageStartupMessages({
  library(minfi)
  library(IlluminaHumanMethylation450kanno.ilmn12.hg19)
})
grSet450k <- mapToGenome(mset450k.andBack)

sexChroms <- GRanges(seqnames = c("chrX", "chrY"),
                     IRanges(start = c(1, 1),
                             end = c(155270560, 59373566)),
                     strand = c("*", "*"))
summary(subsetByOverlaps(grSet450k, sexChroms))
dim(subsetByOverlaps(grSet450k, sexChroms))


## -----------------------------------------------------------------------------
#| label: subsetMinfi
## perhaps more topical:
suppressPackageStartupMessages({
  library(TxDb.Hsapiens.UCSC.hg19.knownGene)
  library(Homo.sapiens)
})
txdb <- TxDb.Hsapiens.UCSC.hg19.knownGene

KDM6AEntrezID <- org.Hs.egSYMBOL2EG[["KDM6A"]]
txs.KDM6A <- transcriptsBy(txdb, "gene")[[KDM6AEntrezID]]
tss.KDM6A <- unique(resize(txs.KDM6A, 1, fix = "start")) ## two start sites
promoters.KDM6A <- flank(tss.KDM6A, 100) ## an arbitrary distance upstream
show(subsetByOverlaps(grSet450k, promoters.KDM6A)) ## probes in this window


## -----------------------------------------------------------------------------
sessionInfo()