## ----message=FALSE, warning=FALSE, include=FALSE------------------------------
library(sesame)
library(wheatmap)
library(dplyr)
options(rmarkdown.html_vignette.check_title = FALSE)

## -----------------------------------------------------------------------------
sesameDataCache("MM285")

## ----eval=FALSE---------------------------------------------------------------
#  res_grn = sesameDataDownload("204637490002_R05C01_Grn.idat")
#  res_red = sesameDataDownload("204637490002_R05C01_Red.idat")
#  pfx = sprintf("%s/204637490002_R05C01", res_red$dest_dir)

## ----eval=FALSE---------------------------------------------------------------
#  sset = readIDATpair(pfx)

## ----eval=FALSE---------------------------------------------------------------
#  openSesame(idat_dir)

## ----include=FALSE------------------------------------------------------------
sset = sesameDataGet('MM285.1.NOD.FrontalLobe')

## -----------------------------------------------------------------------------
sset_normalized = sset %>%
                  noob %>%
                  dyeBiasCorrTypeINorm

## -----------------------------------------------------------------------------
betas = sset_normalized %>%
        qualityMask %>%
        detectionMask %>%
        getBetas

## -----------------------------------------------------------------------------
sum(is.na(betas))
head(betas[grep('uk', names(betas))])

## -----------------------------------------------------------------------------
betas = sset_normalized %>%
        detectionMask %>%
        getBetas
sum(is.na(betas))
head(betas[grep('uk', names(betas))])

## -----------------------------------------------------------------------------
betas = sset_normalized %>% getBetas
sum(is.na(betas))

## -----------------------------------------------------------------------------
betas = sset_normalized %>%
        qualityMask %>%
        detectionMask %>%
        getBetas(mask = FALSE)
sum(is.na(betas))

## ----message=FALSE------------------------------------------------------------
betas = sesameDataGet("MM285.10.tissue")$betas
visualizeGene("Igf2", betas = betas, platform="MM285", refversion = "mm10")

## -----------------------------------------------------------------------------
sset <- sesameDataGet('MM285.1.NOD.FrontalLobe')

## -----------------------------------------------------------------------------
betas <- sset %>%
         noob %>%
         dyeBiasCorrTypeINorm %>%
         getBetas

## -----------------------------------------------------------------------------
vafs <- betaToAF(betas)

## -----------------------------------------------------------------------------
strain <- inferStrain(vafs)
strain$pval

## ----fig.width=6, fig.height=5------------------------------------------------
library(ggplot2)
df <- data.frame(strain=names(strain$probs), probs=strain$probs)
ggplot(data = df,  aes(x = strain, y = log(probs))) +
  geom_bar(stat = "identity", color="gray") +
  ggtitle("strain probabilities") +
  scale_x_discrete(position = "top") +
  theme(axis.text.x = element_text(angle = 90), legend.position = "none")

## -----------------------------------------------------------------------------
betas <- sesameDataGet("MM285.10.tissue")$betas[,1:2]

## ----fig.width=6, fig.height=5------------------------------------------------
compareMouseTissueReference(betas)

## -----------------------------------------------------------------------------
betas <- sesameDataGet('MM285.10.tissue')$betas

## -----------------------------------------------------------------------------
predictMouseAgeInMonth(betas[,1])

## ----message=FALSE------------------------------------------------------------
library(SummarizedExperiment)

## ----message=FALSE------------------------------------------------------------
se = sesameDataGet("MM285.10.tissues")[1:100,]
se_ok = (checkLevels(assay(se), colData(se)$sex) &
    checkLevels(assay(se), colData(se)$tissue))
se = se[se_ok,]

## -----------------------------------------------------------------------------
cf_list = summaryExtractCfList(DML(se, ~tissue + sex))

## -----------------------------------------------------------------------------
cf_list = DMR(se, cf_list$sexMale)
topSegments(cf_list) %>% dplyr::filter(Seg.Pval.adj < 0.05)