## ----setup, include=FALSE--------------------------------------------------
knitr::opts_chunk$set(message = FALSE, warning = FALSE)

## ----Load_Packages, message = FALSE----------------------------------------
library(MEAL)
library(brgedata)
library(MultiDataSet)
library(missMethyl)
library(minfi)
library(GenomicRanges)
library(ggplot2)

## ----Methylation_Data------------------------------------------------------
data(brge_methy)
brge_methy
colData(brge_methy)

## ----Expression_Data-------------------------------------------------------
data(brge_gexp)
brge_gexp
lapply(pData(brge_gexp), table)

## --------------------------------------------------------------------------
targetRange <- GRanges("chr11:102600000-103300000")

## ----Meth_Analysis---------------------------------------------------------
cellCounts <- colnames(colData(brge_methy))[3:9]
methRes <- runPipeline(set = brge_methy, variable_names = "sex", covariable_names = cellCounts, range = targetRange)
methRes
names(methRes)

## ----Plot QQ 1-------------------------------------------------------------
plot(methRes, rid = "DiffMean", type = "qq")

## ----Meth_Analysis SVA, eval = FALSE---------------------------------------
#  methRes <- runPipeline(set = brge_methy, variable_names = "sex", range = targetRange, sva = TRUE)
#  methRes

## ----Manhattans------------------------------------------------------------
targetRangeNum <- GRanges("11:102600000-103300000")
plot(methRes, rid = "DiffMean", type = "manhattan", main = "Differences in Means", highlight = targetRangeNum)
plot(methRes, rid = "DiffVar", type = "manhattan", main = "Differences in Variances", highlight = targetRangeNum)

## ----plotRDA---------------------------------------------------------------
plotRDA(object = methRes, 
        pheno = colData(brge_methy)[, "sex", drop = FALSE])

## ----getRDAresults---------------------------------------------------------
getRDAresults(methRes)

## ----Regional plot 1-------------------------------------------------------
plotRegion(rset = methRes, range = targetRange)

## ----Exp show--------------------------------------------------------------
targetRange <- GRanges("chr11:102600000-103300000")
gexpRes <- runPipeline(set = brge_gexp, variable_names = "sex", range = targetRange)
names(gexpRes)

## ----Plot QQ exp 1---------------------------------------------------------
plot(gexpRes, rid = "DiffMean", type = "qq")

## ----Volcano gexp----------------------------------------------------------
plot(gexpRes, rid = "DiffMean", type = "volcano") + ggtitle("Differences in Means")
plot(gexpRes, rid = "DiffVar", type = "volcano") + ggtitle("Differences in Variances")

## ----Manhattans gexp-------------------------------------------------------
targetRange <- GRanges("chr11:102600000-103300000")
plot(gexpRes, rid = "DiffMean", type = "manhattan", main = "Differences in Means", highlight = targetRangeNum)
plot(gexpRes, rid = "DiffVar", type = "manhattan", main = "Differences in Variances", highlight = targetRangeNum)

## ----plotRDA gexp----------------------------------------------------------
brge_gexp$sex <- as.factor(brge_gexp$sex)
plotRDA(object = gexpRes, pheno = pData(brge_gexp)[, "sex", drop = FALSE])

## ----getRDAresults gexp----------------------------------------------------
getRDAresults(gexpRes)

## ----Regional plot exp-----------------------------------------------------
plotRegion(rset = gexpRes, range = targetRange)

## ----Regional plot 2 exp---------------------------------------------------
plotRegion(rset = methRes, rset2 = gexpRes, range = targetRange, results = c("bumphunter", "DiffMean"))

## ----New Multi Meth Exp----------------------------------------------------
multi <- createMultiDataSet()
multi <- add_genexp(multi, brge_gexp)
multi <- add_methy(multi, brge_methy)

## --------------------------------------------------------------------------
multi.filt <- multi[, , targetRange]

## ----Corr Meth Exp---------------------------------------------------------
methExprs <- correlationMethExprs(multi.filt)
head(methExprs)

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