## ----echo = FALSE------------------------------------------------------------- knitr::opts_chunk$set(comment = "", message=FALSE, warning = FALSE) ## ----eval = FALSE------------------------------------------------------------- # if (!requireNamespace("BiocManager", quietly = TRUE)) # install.packages("BiocManager") # # BiocManager::install("methylGSA") ## ----------------------------------------------------------------------------- library(methylGSA) ## ----------------------------------------------------------------------------- library(IlluminaHumanMethylation450kanno.ilmn12.hg19) ## ----eval = FALSE------------------------------------------------------------- # library(IlluminaHumanMethylationEPICanno.ilm10b4.hg19) ## ----------------------------------------------------------------------------- data(cpgtoy) head(cpg.pval, 20) ## ----------------------------------------------------------------------------- res1 = methylglm(cpg.pval = cpg.pval, minsize = 200, maxsize = 500, GS.type = "KEGG") head(res1, 15) ## ----echo=FALSE--------------------------------------------------------------- glm_res = data.frame( "Column" = c("ID", "Description (N/A for user supplied gene set)", "Size", "pvalue", "padj"), "Explanation" = c("Gene set ID", "Gene set description", "Number of genes in gene set", "p-value in logistic regression", "Adjusted p-value") ) knitr::kable(glm_res) ## ----------------------------------------------------------------------------- library(org.Hs.eg.db) genes_04080 = select(org.Hs.eg.db, "04080", "SYMBOL", keytype = "PATH") head(genes_04080) ## ----eval=FALSE--------------------------------------------------------------- # # include all the IDs as the 2nd argument in select function # genes_all_pathway = select(org.Hs.eg.db, as.character(res1$ID), # "SYMBOL", keytype = "PATH") # head(genes_all_pathway) ## ----eval=FALSE--------------------------------------------------------------- # res2 = methylRRA(cpg.pval = cpg.pval, method = "ORA", # minsize = 200, maxsize = 210) # head(res2, 15) ## ----echo=FALSE--------------------------------------------------------------- ora_res = data.frame( "Column" = c("ID", "Description (N/A for user supplied gene set)", "Count", "overlap", "Size", "pvalue", "padj"), "Explanation" = c("Gene set ID", "Gene set description", "Number of significant genes in the gene set", "Names of significant genes in the gene set", "Number of genes in gene set", "p-value in ORA", "Adjusted p-value") ) knitr::kable(ora_res) ## ----eval=FALSE--------------------------------------------------------------- # res3 = methylRRA(cpg.pval = cpg.pval, method = "GSEA", # minsize = 200, maxsize = 210) # head(res3, 10) ## ----echo=FALSE--------------------------------------------------------------- gsea_res = data.frame( "Column" = c("ID", "Description (N/A for user supplied gene set)", "Size", "enrichmentScore", "NES", "pvalue", "padj"), "Explanation" = c("Gene set ID", "Gene set description", "Number of genes in gene set", "Enrichment score (see [3] for details)", "Normalized enrichment score (see [3] for details)", "p-value in GSEA", "Adjusted p-value") ) knitr::kable(gsea_res) ## ----eval=FALSE--------------------------------------------------------------- # res4 = methylgometh(cpg.pval = cpg.pval, sig.cut = 0.001, # minsize = 200, maxsize = 210) # head(res4, 15) ## ----------------------------------------------------------------------------- data(GSlisttoy) ## to make the display compact, only a proportion of each gene set is shown head(lapply(GS.list, function(x) x[1:30]), 3) ## ----eval=FALSE--------------------------------------------------------------- # library(BiocParallel) # res_p = methylglm(cpg.pval = cpg.pval, minsize = 200, # maxsize = 500, GS.type = "KEGG", parallel = TRUE) ## ----------------------------------------------------------------------------- data(CpG2Genetoy) head(CpG2Gene) ## ----------------------------------------------------------------------------- FullAnnot = prepareAnnot(CpG2Gene) ## ----------------------------------------------------------------------------- GS.list = GS.list[1:10] res5 = methylRRA(cpg.pval = cpg.pval, FullAnnot = FullAnnot, method = "ORA", GS.list = GS.list, GS.idtype = "SYMBOL", minsize = 100, maxsize = 300) head(res5, 10) ## ----eval=FALSE--------------------------------------------------------------- # res6 = methylglm(cpg.pval = cpg.pval, array.type = "450K", # GS.type = "Reactome", minsize = 100, maxsize = 110) # head(res6, 10) ## ----------------------------------------------------------------------------- barplot(res1, num = 8, colorby = "pvalue") ## ----sessionInfo-------------------------------------------------------------- sessionInfo()