### R code from vignette source 'vignettes/minfi/inst/doc/minfi.Rnw'

###################################################
### code chunk number 1: options
###################################################
options(width=70)


###################################################
### code chunk number 2: load
###################################################
require(minfi)
require(minfiData)


###################################################
### code chunk number 3: baseDir
###################################################
baseDir <- system.file("extdata", package = "minfiData")
list.files(baseDir)


###################################################
### code chunk number 4: baseDir
###################################################
list.files(file.path(baseDir, "5723646052"))


###################################################
### code chunk number 5: sampleSheet
###################################################
targets <- read.450k.sheet(baseDir)
targets


###################################################
### code chunk number 6: BasenameColumn
###################################################
sub(baseDir, "", targets$Basename)


###################################################
### code chunk number 7: paths
###################################################
RGset <- read.450k.exp(base = baseDir, targets = targets)


###################################################
### code chunk number 8: pData
###################################################
RGset
pd <- pData(RGset)
pd[,1:4]


###################################################
### code chunk number 9: read2
###################################################
RGset2 = read.450k.exp(file.path(baseDir, "5723646052"))
RGset3 = read.450k.exp(baseDir, recursive = TRUE)


###################################################
### code chunk number 10: sampleSheet2
###################################################
targets2 <- read.csv(file.path(baseDir, "SampleSheet.csv"), 
                     stringsAsFactors = FALSE, skip = 7)
targets2


###################################################
### code chunk number 11: Basename
###################################################
targets2$Basename <- file.path(baseDir, targets2$Sentrix_ID, 
                               paste0(targets2$Sentrix_ID, 
                                      targets2$Sentrix_Position))


###################################################
### code chunk number 12: qcReport-quick (eval = FALSE)
###################################################
## qcReport(RGset, sampNames = pd$Sample_Name, 
##          sampGroups = pd$Sample_Group, pdf = "qcReport.pdf")


###################################################
### code chunk number 13: qcReport-density
###################################################
densityPlot(RGset, sampGroups = pd$Sample_Group, 
            main = "Beta", xlab = "Beta")


###################################################
### code chunk number 14: qcReport-bean
###################################################
par(oma=c(2,10,1,1))
densityBeanPlot(RGset, sampGroups = pd$Sample_Group, 
                sampNames = pd$Sample_Name)


###################################################
### code chunk number 15: qcReport-stripplot
###################################################
controlStripPlot(RGset, controls="BISULFITE CONVERSION II", 
                 sampNames = pd$Sample_Name)


###################################################
### code chunk number 16: Msetraw
###################################################
MSet.raw <- preprocessRaw(RGset)


###################################################
### code chunk number 17: allMsets
###################################################
MSet.norm <- preprocessIllumina(RGset, bg.correct = TRUE,
                                 normalize = "controls", reference = 2)


###################################################
### code chunk number 18: MSet
###################################################
getMeth(MSet.raw)[1:4,1:3]
getUnmeth(MSet.raw)[1:4,1:3]
getBeta(MSet.raw, type = "Illumina")[1:4,1:3]
getM(MSet.raw)[1:4,1:3]


###################################################
### code chunk number 19: qcReport-mdsplot2
###################################################
mdsPlot(MSet.norm, numPositions = 1000, sampGroups = pd$Sample_Group, 
	sampNames = pd$Sample_Name)


###################################################
### code chunk number 20: preprocessSwan
###################################################
Mset.swan <- preprocessSWAN(RGsetEx, MsetEx)


###################################################
### code chunk number 21: plotBetaType
###################################################
par(mfrow=c(1,2))
plotBetasByType(MsetEx[,1], main = "Raw")
plotBetasByType(Mset.swan[,1], main = "SWAN")


###################################################
### code chunk number 22: subset-mset
###################################################
mset <- MSet.norm[1:20000,]


###################################################
### code chunk number 23: dmpFinder-categorical
###################################################
table(pd$Sample_Group)
M <- getM(mset, type = "beta", betaThreshold = 0.001)
dmp <- dmpFinder(M, pheno=pd$Sample_Group, type="categorical")
head(dmp)


###################################################
### code chunk number 24: plot-dmps-categorical
###################################################
cpgs <- rownames(dmp)[1:4]
par(mfrow=c(2,2))
plotCpg(mset, cpg=cpgs, pheno=pd$Sample_Group)


###################################################
### code chunk number 25: set-seed
###################################################
set.seed(123)


###################################################
### code chunk number 26: sim-pheno
###################################################
continuousPheno <- rnorm(nrow(pd))


###################################################
### code chunk number 27: dmpFinder-continuous
###################################################
dmp <- dmpFinder(mset, pheno=continuousPheno, type="continuous")
dmp[1:3,]


###################################################
### code chunk number 28: filter-dmp
###################################################
dmp <- subset(dmp, abs(beta)>1)


###################################################
### code chunk number 29: plot-dmps-continuous
###################################################
cpgs <- rownames(dmp)[1:4]
par(mfrow=c(2,2))
plotCpg(mset, cpg=cpgs, type="continuous",
        pheno=continuousPheno, xlab="Phenotype 1")


###################################################
### code chunk number 30: manifest
###################################################
IlluminaHumanMethylation450kmanifest
head(getProbeInfo(IlluminaHumanMethylation450kmanifest, type = "I"), n = 3)
head(getProbeInfo(IlluminaHumanMethylation450kmanifest, type = "II"), n = 3)
head(getProbeInfo(IlluminaHumanMethylation450kmanifest, type = "Control"), n = 3)


###################################################
### code chunk number 31: sessionInfo
###################################################
toLatex(sessionInfo())