## ----setup, echo = FALSE, warning = FALSE-------------------------------------
suppressMessages(library(knitr))
suppressMessages(library(GenomicRanges))
suppressMessages(library(limma))
suppressMessages(library(minfi))
suppressMessages(library(ExperimentHub))
suppressMessages(library(recountmethylation))
knitr::opts_chunk$set(echo = TRUE, eval = TRUE, warning = FALSE, 
  message = FALSE)

## ---- echo = TRUE, message = TRUE---------------------------------------------
sm <- as.data.frame(smfilt(get_servermatrix()))
if(is(sm, "data.frame")){knitr::kable(sm, align = "c")}

## -----------------------------------------------------------------------------
hub = ExperimentHub::ExperimentHub() # connect to the hubs
rmdat <- AnnotationHub::query(hub, "recountmethylation") # query the hubs
rmdat

## -----------------------------------------------------------------------------
eid <- "EH3778" # h5 test file id
fpath <- rmdat[[eid]] # download with default caching
rhdf5::h5ls(fpath)

## -----------------------------------------------------------------------------
dn <- "remethdb-h5se_gr-test_0-0-1_1590090412"
path <- system.file("extdata", dn, package = "recountmethylation")
h5se.test <- HDF5Array::loadHDF5SummarizedExperiment(path)

## -----------------------------------------------------------------------------
class(h5se.test)

## -----------------------------------------------------------------------------
dim(h5se.test)

## -----------------------------------------------------------------------------
summary(h5se.test)

## -----------------------------------------------------------------------------
h5se.md <- minfi::pData(h5se.test)
dim(h5se.md)

## -----------------------------------------------------------------------------
colnames(h5se.md)

## -----------------------------------------------------------------------------
h5se.bm <- minfi::getBeta(h5se.test)
dim(h5se.bm)

## -----------------------------------------------------------------------------
colnames(h5se.bm) <- h5se.test$gsm
knitr::kable(head(h5se.bm), align = "c")

## -----------------------------------------------------------------------------
an <- minfi::getAnnotation(h5se.test)
dim(an)

## -----------------------------------------------------------------------------
colnames(an)

## -----------------------------------------------------------------------------
ant <- as.matrix(t(an[c(1:4), c(1:3, 5:6, 9, 19, 24, 26)]))
knitr::kable(ant, align = "c")

## -----------------------------------------------------------------------------
dn <- "remethdb-h5_rg-test_0-0-1_1590090412.h5"
h5.test <- system.file("extdata", "h5test", dn, 
                    package = "recountmethylation")

## -----------------------------------------------------------------------------
h5.rg <- getrg(dbn = h5.test, all.gsm = TRUE)

## -----------------------------------------------------------------------------
h5.red <- minfi::getRed(h5.rg)
h5.green <- minfi::getGreen(h5.rg)
dim(h5.red)

## -----------------------------------------------------------------------------
knitr::kable(head(h5.red), align = "c")

## -----------------------------------------------------------------------------
knitr::kable(head(h5.green), align = "c")

## -----------------------------------------------------------------------------
identical(rownames(h5.red), rownames(h5.green))

## ---- eval = FALSE------------------------------------------------------------
#  # download from GEO
#  dlpath <- tempdir()
#  gsmv <- c("GSM1038308", "GSM1038309")
#  geo.rg <- gds_idat2rg(gsmv, dfp = dlpath)
#  colnames(geo.rg) <- gsub("\\_.*", "", colnames(geo.rg))

## ---- eval = FALSE------------------------------------------------------------
#  geo.red <- minfi::getRed(geo.rg)
#  geo.green <- minfi::getGreen(geo.rg)

## ---- eval = FALSE------------------------------------------------------------
#  int.addr <- intersect(rownames(geo.red), rownames(h5.red))
#  geo.red <- geo.red[int.addr,]
#  geo.green <- geo.green[int.addr,]
#  geo.red <- geo.red[order(match(rownames(geo.red), rownames(h5.red))),]
#  geo.green <- geo.green[order(match(rownames(geo.green), rownames(h5.green))),]
#  identical(rownames(geo.red), rownames(h5.red))
#  identical(rownames(geo.green), rownames(h5.green))
#  class(h5.red) <- "integer"
#  class(h5.green) <- "integer"

## ---- eval = FALSE------------------------------------------------------------
#  identical(geo.red, h5.red)

## ---- eval = FALSE------------------------------------------------------------
#  identical(geo.green, h5.green)

## ---- eval = FALSE------------------------------------------------------------
#  geo.gr <- minfi::preprocessNoob(geo.rg)

## ---- eval = FALSE------------------------------------------------------------
#  geo.bm <- as.matrix(minfi::getBeta(geo.gr))

## ---- eval = FALSE------------------------------------------------------------
#  h5se.bm <- as.matrix(h5se.bm)
#  int.cg <- intersect(rownames(geo.bm), rownames(h5se.bm))
#  geo.bm <- geo.bm[int.cg,]
#  geo.bm <- geo.bm[order(match(rownames(geo.bm), rownames(h5se.bm))),]

## ---- eval = FALSE------------------------------------------------------------
#  identical(summary(geo.bm), summary(h5se.bm))

## ---- eval = FALSE------------------------------------------------------------
#  identical(rownames(geo.bm), rownames(h5se.bm))

## ---- eval = FALSE------------------------------------------------------------
#  detectionP(rg[,1:50])
#  "Error in .local(Red, Green, locusNames, controlIdx, TypeI.Red, TypeI.Green, dim(Red_grid) == dim(detP_sink_grid) are not all TRUE"

## ---- eval = FALSE------------------------------------------------------------
#  preprocessFunnorm(rg)
#  "Error: 'preprocessFunnorm()' only supports matrix-backed minfi objects.""

## ---- eval = FALSE------------------------------------------------------------
#  rg.h5se <- loadHDF5SummarizedExperiment(rg.path) # full h5se RGChannelSet
#  rg.sub <- rg.h5se[,c(1:20)] # subset samples of interest
#  # re-make as non-DA object
#  rg.new <- RGChannelSet(Red = getRed(rg.sub),
#                         Green = getGreen(rg.sub),
#                         annotation = annotation(rg.sub))
#  gr <- preprocessFunnorm(rg.new) # repeat preprocessing

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