## ----setup_data, message=FALSE, warning=FALSE------------------------------
library(dplyr)
library(biotmle)
library(biotmleData)
data(illuminaData)
suppressMessages(library(SummarizedExperiment))
"%ni%" = Negate("%in%")

## ----clean_data------------------------------------------------------------
# discretize "age" in the phenotype-level data
colData(illuminaData) <- colData(illuminaData) %>%
  data.frame %>%
  dplyr::mutate(age = as.numeric(age > median(age))) %>%
  DataFrame

# specify column index of treatment/exposure variable of interest
varInt_index <- which(names(colData(illuminaData)) %in% "benzene")

## ----biomarkerTMLE_eval, eval=FALSE----------------------------------------
#  biomarkerTMLEout <- biomarkertmle(se = illuminaData,
#                                    varInt = varInt_index,
#                                    g_lib = c("SL.glmnet", "SL.randomForest",
#                                              "SL.polymars", "SL.mean"),
#                                    Q_lib = c("SL.glmnet", "SL.randomForest",
#                                              "SL.nnet", "SL.mean")
#                                   )

## ----load_biomarkerTMLE_result, echo=FALSE---------------------------------
data(biomarkertmleOut)

## ----limmaTMLE_eval--------------------------------------------------------
limmaTMLEout <- modtest_ic(biotmle = biomarkerTMLEout)

## ----pval_hist_limma_adjp--------------------------------------------------
plot(x = limmaTMLEout, type = "pvals_adj")

## ----pval_hist_limma_rawp--------------------------------------------------
plot(x = limmaTMLEout, type = "pvals_raw")

## ----heatmap_limma_results-------------------------------------------------
varInt_index <- which(names(colData(illuminaData)) %in% "benzene")
designVar <- as.data.frame(colData(illuminaData))[, varInt_index]
designVar <- as.numeric(designVar == max(designVar))

# build heatmap
heatmap_ic(x = limmaTMLEout, left.label = "none", scale = TRUE,
           clustering.method = "hierarchical", row.dendrogram = TRUE,
           design = designVar, FDRcutoff = 0.05, top = 25)

## ----volcano_plot_limma_results--------------------------------------------
volcano_ic(biotmle = limmaTMLEout)

## ----sessionInfo, echo=FALSE-----------------------------------------------
sessionInfo()