Meta-analysis dataset

1 Main libraries and configuration

library(dplyr)
library(ExperimentHub)
library(GenomicRanges)
library(coMethDMR)

2 Meta-analysis of Genomic Regions

2.1 Paths

dir.result <- "meta_analysis_region_results/"
dir.result.meta.analysis <- file.path(dir.result, "step1_meta_analysis/")
dir.result.smoking_crosshyb <- file.path(dir.result, "step2_smoking_cross_anotation/")
dir.result.comp <- file.path(dir.result, "step3_comp/")
dir.result.cpg.vs.dmr <- file.path(dir.result, "step4_dmr_vs_cpgs/")
dir.result.lola <- file.path(dir.result, "step5_lola/")
dir.result.enrichment <- file.path(dir.result, "step6_enrichment/")
dir.result.pathway <- file.path(dir.result, "step7_pathway/")
data.final <- "meta_analysis_single_cpg_results/"
for(p in grep("dir",ls(),value = T)) dir.create(get(p),recursive = TRUE,showWarnings = FALSE)

2.2 Import datasets and pre-process for each cohort

library(dplyr)
library(tidyr)
preMeta <- function(cohort){
  
  ### Load data
  file <- dir(pattern = paste0(cohort,"_linear_df"),
              path =  "DATASETS/",
              recursive = TRUE,
              full.names = TRUE,
              ignore.case = TRUE)
  message("Reading file: ", file)
  linear.results.final <- readRDS(file)
  
  ### select the most sig cometh region for each input region
  pva.col <- grep("_pVal",colnames(linear.results.final),value = TRUE)
  colnames(linear.results.final)[grep("inputRegion",colnames(linear.results.final))] <- "inputRegion"
  
  # !! is used to unquote an input 
  # https://dplyr.tidyverse.org/articles/programming.html
  result_sig <- linear.results.final %>%
    group_by(inputRegion) %>%
    filter((!!as.symbol(pva.col)) == min(!!as.symbol(pva.col)))
  
  data.frame(result_sig, stringsAsFactors = FALSE)
}

Gasparoni <- preMeta(cohort = "Gasparoni") #dim: 39605  8
London_PFC <- preMeta(cohort = "London") #dim: 40000 8
MtSinai <- preMeta(cohort = "MtSinai") #dim:  38619 8
ROSMAP <- preMeta(cohort = "ROSMAP") #dim: 38562 8

2.3 Merge cohorts

### merge datasets
cohort_ls <- list(
  Gasparoni = Gasparoni,
  London_PFC = London_PFC,
  MtSinai = MtSinai,
  ROSMAP = ROSMAP
)

### outer join input region
multi_cohorts <- Reduce(
  function(x,y, ...) merge(x, y, by = "inputRegion", all = TRUE, ...),
  cohort_ls
)

2.4 Meta analysis

library(meta)

doParallel::registerDoParallel(cores = parallel::detectCores()/2)
meta_df <- plyr::adply(.data = multi_cohorts, 
                       .margins = 1, 
                       .fun =  function(rowOne_df){
                         
                         est <- rowOne_df[grep("estimate",colnames(rowOne_df))] %>% as.numeric
                         
                         direction <-  paste(ifelse(
                           is.na(est), ".",
                           ifelse(est > 0, "+", "-")
                         ),collapse = "")
                         
                         se <- rowOne_df[grep("se",colnames(rowOne_df))] %>% as.numeric
                         cohort <- gsub("_se","",grep("se",colnames(rowOne_df),value = T))
                         rowOne_reform_df <- data.frame(
                           cohort = cohort,
                           est = est,
                           se = se,
                           stringsAsFactors = FALSE
                         )
                         
                         f <- metagen(
                           TE = est,
                           seTE = se,
                           data = rowOne_reform_df
                         )
                         
                         tibble::tibble(
                           inputRegion = rowOne_df$inputRegion,
                           estimate = f$TE.fixed,
                           se = f$seTE.fixed,
                           pVal.fixed = f$pval.fixed,
                           pVal.random = f$pval.random,
                           pValQ = f$pval.Q,
                           direction = direction
                         )
                       }  , .progress = "time",
                       .parallel = TRUE,
                       .id = NULL
)

### create final pVal
meta_df$pVal.final <- ifelse(
  meta_df$pValQ > 0.05, meta_df$pVal.fixed, meta_df$pVal.random
)

### calculate fdr
meta_df$fdr <- p.adjust(meta_df$pVal.final, method = "fdr")

### order meta_df
meta_final_df <- meta_df[, c(grep("_",colnames(meta_df),invert = T),
                             grep("_",colnames(meta_df),invert = F))]
meta_final_ordered_df <- meta_final_df[order(meta_final_df$pVal.final),]

2.5 Add annotation to input regions

### find annotations
library(coMethDMR)

splited_input <- meta_final_ordered_df %>% 
  tidyr::separate(col = inputRegion,into =  c("chrom", "start", "end"),remove = FALSE)
splited_input_annot <- AnnotateResults(splited_input[,c("chrom", "start", "end")],nCores_int = 10) 

### merge annotation with meta analysis data
meta_ordered_withAnnot_df <- cbind(
  meta_final_ordered_df, splited_input_annot[, 4:7]
)

### order columns
meta_ordered_withAnnot_df <- meta_ordered_withAnnot_df %>% 
  dplyr::select(c(1,
                  (ncol(meta_ordered_withAnnot_df) - 3):ncol(meta_ordered_withAnnot_df),
                  2:(ncol(meta_ordered_withAnnot_df) - 4))
  )
### save dataset
write.csv(
  meta_ordered_withAnnot_df,
  paste0(dir.result.meta.analysis, "meta_analysis_ALL_df.csv"),
  row.names = FALSE
)

meta_all_sig <- meta_ordered_withAnnot_df[
  !is.na(meta_ordered_withAnnot_df$fdr) &
    (meta_ordered_withAnnot_df$fdr < 0.05), 
  ] #dim: 478 41
row.names(meta_all_sig) <- NULL

write.csv(
  meta_all_sig,
  paste0(dir.result.meta.analysis, "meta_analysis_ALL_sig_df.csv"),
  row.names = FALSE
)

2.6 Filtering coMethDMRs

2.6.1 Annotate coMethDMRs with crosshybrdizing probes chen et al. (2013) & smoking probes

2.6.1.1 Get crosshybrdizing probes

### call in all cross hybridizing probes
eh = ExperimentHub()

## snapshotDate(): 2020-04-03

query(eh, "DMRcate")

## ExperimentHub with 9 records
## # snapshotDate(): 2020-04-03
## # $dataprovider: Ensembl, Springer, SickKids, Springer, Bioconductor
## # $species: Homo sapiens, Mus musculus
## # $rdataclass: GeneRegionTrack, GRanges, character, data.frame
## # additional mcols(): taxonomyid, genome, description,
## #   coordinate_1_based, maintainer, rdatadateadded, preparerclass, tags,
## #   rdatapath, sourceurl, sourcetype 
## # retrieve records with, e.g., 'object[["EH3129"]]' 
## 
##            title          
##   EH3129 | crosshyb       
##   EH3130 | snpsall        
##   EH3131 | XY.probes      
##   EH3132 | hg19.generanges
##   EH3133 | hg19.grt       
##   EH3134 | hg38.generanges
##   EH3135 | hg38.grt       
##   EH3136 | mm10.generanges
##   EH3137 | mm10.grt

crosshyb <- eh[["EH3129"]]

## see ?DMRcatedata and browseVignettes('DMRcatedata') for documentation

## loading from cache

2.6.1.2 Get significant smoking probes

smoking.file <- "https://www.ncbi.nlm.nih.gov/pmc/articles/PMC5267325/bin/NIHMS817273-supplement-001506_-_Supplemental_Tables.xlsx"

if(!file.exists(basename(smoking.file))) downloader::download(smoking.file,basename(smoking.file))

smoking <- readxl::read_xlsx(
  basename(smoking.file),
  sheet = "02",
  skip = 2
)
smoking.sig.probes <- smoking %>% dplyr::filter(`P-value` < 1*10^(-7)) %>% pull("Probe ID") 
length(smoking.sig.probes)

## [1] 2623

2.6.1.3 Load significant coMethDMRs identified in meta-analaysis

# Call in meta analysis final results
meta_all <- readr::read_csv(
  dir(dir.result.meta.analysis, pattern = "meta_analysis_ALL_sig_df.csv",full.names = TRUE),
  col_types = readr::cols()
)

# Find files with regions and probes
files <- dir(pattern = paste0(".*_residuals_cometh_input_ls.rds"),
             recursive = T,
             full.names = TRUE,
             ignore.case = T)
files <- grep("GASPARONI|LONDON_PFC|MtSinai|ROSMAP",files,value = TRUE,ignore.case = TRUE)

# Read files and Limit the cohort_ls to cohort_coMethRegion in meta_all 
cometh.probes.list <- lapply(files, function (f){
  print(f)
  all.clusters <- readRDS(f)$coMeth_ls # Read files
  cohort <- f %>% dirname %>% dirname %>% basename # get cohort from folder name
  
  col <- grep(paste0(cohort,"_coMethRegion"),
              colnames(meta_all),
              ignore.case = TRUE,
              value = TRUE) # get column with cohort sig regions
  
  # keep sig regions only
  all.clusters[names(all.clusters) %in% meta_all[[col]]]
})

## [1] "./DATASETS/GASPARONI/step10_residuals/Gasparoni_residuals_cometh_input_ls.RDS"
## [1] "./DATASETS/LONDON/step10_residuals/London_PFC_residuals_cometh_input_ls.RDS"
## [1] "./DATASETS/MtSinai/step10_residuals/MtSinai_residuals_cometh_input_ls.RDS"
## [1] "./DATASETS/ROSMAP/step10_residuals/ROSMAP_residuals_cometh_input_ls.RDS"

names(cometh.probes.list) <-  files %>% dirname %>% dirname %>% basename

lapply(cometh.probes.list,length)

## $GASPARONI
## [1] 476
## 
## $LONDON
## [1] 478
## 
## $MtSinai
## [1] 469
## 
## $ROSMAP
## [1] 474

2.6.1.4 Map probes in each region to smoking and crosshybrdizing

extractCrosHybridization <- function(probes.list){
  crosshyb.extracted <- plyr::laply(probes.list,function(probes){
    paste(intersect(probes, crosshyb), 
          collapse = ", "
    )
  })
  smoking.extracted <- plyr::laply(probes.list,function(probes){
    paste(intersect(probes, smoking.sig.probes), 
          collapse = ", "
    )
  })
  tibble::tibble(
    "coMethRegion" = names(probes.list),
    "crossHyb" = crosshyb.extracted, 
    "crossHyb_bi" = ifelse(crosshyb.extracted == "",0,1),
    "smoke" = smoking.extracted,
    "smoke_bi" = ifelse(smoking.extracted == "",0,1)
  )
}

Gasparoni_crossHyb_df <- extractCrosHybridization(cometh.probes.list$GASPARONI)
colnames(Gasparoni_crossHyb_df) <- paste0("GASPARONI_",colnames(Gasparoni_crossHyb_df))
plyr::count(Gasparoni_crossHyb_df, vars = grep("_bi",colnames(Gasparoni_crossHyb_df),value = TRUE))

London_PFC_crossHyb_df <- extractCrosHybridization(cometh.probes.list$LONDON)
colnames(London_PFC_crossHyb_df) <- paste0("London_",colnames(London_PFC_crossHyb_df))
plyr::count(London_PFC_crossHyb_df, vars = grep("_bi",colnames(London_PFC_crossHyb_df),value = TRUE))

MtSinai_crossHyb_df <- extractCrosHybridization(cometh.probes.list$MtSinai)
colnames(MtSinai_crossHyb_df) <- paste0("MtSinai_",colnames(MtSinai_crossHyb_df))
plyr::count(MtSinai_crossHyb_df, vars = grep("_bi",colnames(MtSinai_crossHyb_df),value = TRUE))

ROSMAP_crossHyb_df <- extractCrosHybridization(cometh.probes.list$ROSMAP)
colnames(ROSMAP_crossHyb_df) <- paste0("ROSMAP_",colnames(ROSMAP_crossHyb_df))
plyr::count(ROSMAP_crossHyb_df, vars = grep("_bi",colnames(ROSMAP_crossHyb_df),value = TRUE))

2.6.1.5 Merge smoking and crossHyb probes information with meta analysis results

meta_all_final <- meta_all %>% left_join(Gasparoni_crossHyb_df)  %>% 
  left_join(London_PFC_crossHyb_df) %>% 
  left_join(MtSinai_crossHyb_df) %>% 
  left_join(ROSMAP_crossHyb_df)

## Joining, by = "GASPARONI_coMethRegion"

## Joining, by = "London_coMethRegion"

## Joining, by = "MtSinai_coMethRegion"

## Joining, by = "ROSMAP_coMethRegion"

### Add information with input regions with any cross hybridization in cohorts 
meta_all_final$crossHyb_bi <- rowSums(meta_all_final[,grep("crossHyb_bi",colnames(meta_all_final))],na.rm = TRUE) > 0
meta_all_final$smoke_bi <- rowSums(meta_all_final[,grep("smoke_bi",colnames(meta_all_final))],na.rm = TRUE) > 0

# Sort by region meta analysis FDR
# Cluster columns of the projects together
meta_all_final <- meta_all_final[
  order(meta_all_final$fdr),
  c(
    grep("Gasparoni|MtSinai|London|ROSMAP",colnames(meta_all_final),ignore.case = TRUE,invert = TRUE),
    grep("Gasparoni",colnames(meta_all_final),ignore.case = TRUE),
    grep("MtSinai",colnames(meta_all_final),ignore.case = TRUE),
    grep("London",colnames(meta_all_final),ignore.case = TRUE),
    grep("ROSMAP",colnames(meta_all_final),ignore.case = TRUE)
  )
  ]
str(meta_all_final)

## tibble [478 × 59] (S3: tbl_df/tbl/data.frame)
##  $ inputRegion            : chr [1:478] "chr19:49220102-49220485" "chr7:27153580-27153944" "chr7:27146237-27146445" "chr7:27154262-27155548" ...
##  $ UCSC_RefGene_Group     : chr [1:478] "1stExon;5'UTR;Body;TSS1500;TSS200" "1stExon;5'UTR;TSS1500;TSS200" "3'UTR" "5'UTR;TSS1500" ...
##  $ UCSC_RefGene_Accession : chr [1:478] "NM_001130915;NM_182574" "NM_030661;NM_153631;NM_153632" "NM_030661;NM_153631;NM_153632" "NM_030661;NM_153631;NM_153632" ...
##  $ UCSC_RefGene_Name      : chr [1:478] "MAMSTR" "HOXA3" "HOXA3" "HOXA3" ...
##  $ Relation_to_Island     : chr [1:478] "N_Shelf" "Island;N_Shore" "Island" "Island;N_Shore;S_Shore" ...
##  $ estimate               : num [1:478] 0.0607 0.0838 0.0659 0.0472 0.0449 ...
##  $ se                     : num [1:478] 0.00751 0.01085 0.00862 0.00619 0.00634 ...
##  $ pVal.fixed             : num [1:478] 6.02e-16 1.16e-14 2.23e-14 2.37e-14 1.35e-12 ...
##  $ pVal.random            : num [1:478] 4.24e-14 1.45e-11 3.03e-11 4.09e-10 1.35e-12 ...
##  $ pValQ                  : num [1:478] 0.341 0.286 0.291 0.237 0.741 ...
##  $ direction              : chr [1:478] "++++" "++++" "++++" "++++" ...
##  $ pVal.final             : num [1:478] 6.02e-16 1.16e-14 2.23e-14 2.37e-14 1.35e-12 ...
##  $ fdr                    : num [1:478] 2.41e-11 2.32e-10 2.37e-10 2.37e-10 1.08e-08 ...
##  $ crossHyb_bi            : logi [1:478] FALSE FALSE FALSE FALSE FALSE FALSE ...
##  $ smoke_bi               : logi [1:478] FALSE FALSE FALSE FALSE FALSE FALSE ...
##  $ GASPARONI_nCoMethRegion: num [1:478] 1 1 1 1 1 0 1 1 1 1 ...
##  $ GASPARONI_coMethRegion : chr [1:478] "chr19:49220102-49220235" "chr7:27153580-27153847" "chr7:27146237-27146445" "chr7:27154387-27155548" ...
##  $ GASPARONI_nCpGs        : num [1:478] 4 6 4 14 4 3 4 5 3 8 ...
##  $ Gasparoni_estimate     : num [1:478] 0.0409 0.0664 0.0559 0.0261 0.0487 ...
##  $ Gasparoni_se           : num [1:478] 0.0264 0.0336 0.0345 0.0189 0.0218 ...
##  $ Gasparoni_pVal         : num [1:478] 0.1289 0.054 0.1121 0.1723 0.0306 ...
##  $ Gasparoni_fdr          : num [1:478] 0.996 0.996 0.996 0.996 0.996 ...
##  $ GASPARONI_crossHyb     : chr [1:478] "" "" "" "" ...
##  $ GASPARONI_crossHyb_bi  : num [1:478] 0 0 0 0 0 0 0 0 0 0 ...
##  $ GASPARONI_smoke        : chr [1:478] "" "" "" "" ...
##  $ GASPARONI_smoke_bi     : num [1:478] 0 0 0 0 0 0 0 0 1 0 ...
##  $ MtSinai_nCoMethRegion  : num [1:478] 1 1 1 1 1 0 1 1 1 1 ...
##  $ MtSinai_coMethRegion   : chr [1:478] "chr19:49220102-49220235" "chr7:27153580-27153847" "chr7:27146237-27146445" "chr7:27154720-27155548" ...
##  $ MtSinai_nCpGs          : num [1:478] 4 6 4 10 4 3 4 5 3 7 ...
##  $ MtSinai_estimate       : num [1:478] 0.0846 0.118 0.1051 0.0662 0.0604 ...
##  $ MtSinai_se             : num [1:478] 0.0167 0.027 0.0225 0.0131 0.0158 ...
##  $ MtSinai_pVal           : num [1:478] 1.45e-06 2.57e-05 7.49e-06 1.45e-06 2.08e-04 ...
##  $ MtSinai_fdr            : num [1:478] 0.0188 0.0499 0.0322 0.0188 0.1232 ...
##  $ MtSinai_crossHyb       : chr [1:478] "" "" "" "" ...
##  $ MtSinai_crossHyb_bi    : num [1:478] 0 0 0 0 0 0 0 0 0 0 ...
##  $ MtSinai_smoke          : chr [1:478] "" "" "" "" ...
##  $ MtSinai_smoke_bi       : num [1:478] 0 0 0 0 0 0 0 0 1 0 ...
##  $ London_nCoMethRegion   : num [1:478] 1 1 1 2 1 1 1 1 1 1 ...
##  $ London_coMethRegion    : chr [1:478] "chr19:49220102-49220235" "chr7:27153580-27153847" "chr7:27146237-27146445" "chr7:27154262-27154387" ...
##  $ London_nCpGs           : num [1:478] 4 6 4 4 4 3 3 5 3 7 ...
##  $ London_estimate        : num [1:478] 0.0497 0.099 0.0632 0.0526 0.0415 ...
##  $ London_se              : num [1:478] 0.01287 0.01991 0.0137 0.01241 0.00945 ...
##  $ London_pVal            : num [1:478] 2.07e-04 3.10e-06 1.31e-05 5.41e-05 3.05e-05 ...
##  $ London_fdr             : num [1:478] 0.0958 0.0297 0.0439 0.0631 0.0474 ...
##  $ London_crossHyb        : chr [1:478] "" "" "" "" ...
##  $ London_crossHyb_bi     : num [1:478] 0 0 0 0 0 0 0 0 0 0 ...
##  $ London_smoke           : chr [1:478] "" "" "" "" ...
##  $ London_smoke_bi        : num [1:478] 0 0 0 0 0 0 0 0 1 0 ...
##  $ ROSMAP_nCoMethRegion   : num [1:478] 1 1 1 1 1 0 1 1 1 1 ...
##  $ ROSMAP_coMethRegion    : chr [1:478] "chr19:49220102-49220235" "chr7:27153580-27153847" "chr7:27146237-27146445" "chr7:27154720-27155548" ...
##  $ ROSMAP_nCpGs           : num [1:478] 4 6 4 10 3 3 3 5 3 7 ...
##  $ ROSMAP_estimate        : num [1:478] 0.0622 0.0649 0.0555 0.0393 0.0408 ...
##  $ ROSMAP_se              : num [1:478] 0.01224 0.01641 0.01374 0.00955 0.01147 ...
##  $ ROSMAP_pVal            : num [1:478] 4.94e-07 8.51e-05 6.00e-05 4.39e-05 4.03e-04 ...
##  $ ROSMAP_fdr             : num [1:478] 0.00955 0.12481 0.10655 0.10456 0.23625 ...
##  $ ROSMAP_crossHyb        : chr [1:478] "" "" "" "" ...
##  $ ROSMAP_crossHyb_bi     : num [1:478] 0 0 0 0 0 0 0 0 0 0 ...
##  $ ROSMAP_smoke           : chr [1:478] "" "" "" "" ...
##  $ ROSMAP_smoke_bi        : num [1:478] 0 0 0 0 0 0 0 0 1 0 ...

2.6.1.6 Save

write.csv(
  meta_all_final,
  paste0(dir.result.smoking_crosshyb, "meta_analysis_ALL_sig_add_crossHyb_df.csv"),
  row.names = FALSE
)

meta_all_final %>% 
  dplyr::filter(smoke_bi == 0 & crossHyb_bi == 0)  %>% 
  DT::datatable(filter = 'top',
                style = "bootstrap",
                extensions = 'Buttons',
                options = list(scrollX = TRUE, 
                               dom = 'Bfrtip',
                               buttons = I('colvis'),
                               keys = TRUE, 
                               pageLength = 10), 
                rownames = FALSE,
                caption = "meta-analysis results")

2.6.2 Overlap with comb-p DMRs

2.6.2.1 Import datasets

library(GenomicRanges)

linear_sig <- readr::read_csv(
  dir(dir.result.smoking_crosshyb, pattern = "meta_analysis_ALL_sig_add_crossHyb_df.csv",full.names = TRUE),
  col_types = readr::cols()
)
linear_sig_input_gr <- linear_sig %>% 
  tidyr::separate("inputRegion",c("chrom","start","end")) %>%
  makeGRangesFromDataFrame()
linear_sig_input_gr

## GRanges object with 478 ranges and 0 metadata columns:
##         seqnames              ranges strand
##            <Rle>           <IRanges>  <Rle>
##     [1]    chr19   49220102-49220485      *
##     [2]     chr7   27153580-27153944      *
##     [3]     chr7   27146237-27146445      *
##     [4]     chr7   27154262-27155548      *
##     [5]     chr7   27179161-27179432      *
##     ...      ...                 ...    ...
##   [474]    chr16   87886871-87886933      *
##   [475]    chr13 113422497-113422671      *
##   [476]    chr12 121476549-121476792      *
##   [477]    chr19       836606-836898      *
##   [478]     chr1   22978556-22978660      *
##   -------
##   seqinfo: 22 sequences from an unspecified genome; no seqlengths

combp <- readxl::read_xlsx("../../Michael/1_7_2020_001-200/1_7_2020_001-200/combp_results.xlsx") # get comb-p results
colnames(combp)[1] <- "chrom"
combp$chrom <- paste0("chr",combp$chrom)
combp_sig <- combp %>% filter(z_sidak_p < 0.05 & n_probes > 2)

combp_sig %>% 
  DT::datatable(filter = 'top',
                style = "bootstrap",
                extensions = 'Buttons',
                options = list(scrollX = TRUE, 
                               dom = 'Bfrtip',
                               buttons = I('colvis'),
                               keys = TRUE, 
                               pageLength = 10), 
                rownames = FALSE,
                caption = "comb-p results")

combp_sig_gr <- makeGRangesFromDataFrame(combp_sig)
combp_sig_gr

## GRanges object with 187 ranges and 0 metadata columns:
##         seqnames              ranges strand
##            <Rle>           <IRanges>  <Rle>
##     [1]     chr7   27153580-27153848      *
##     [2]    chr19   10736006-10736356      *
##     [3]     chr7   27146237-27146446      *
##     [4]     chr7   27170313-27171402      *
##     [5]    chr19   49220102-49220236      *
##     ...      ...                 ...    ...
##   [183]     chr6 149806077-149806132      *
##   [184]     chr6 152702330-152702388      *
##   [185]    chr19   46999224-46999290      *
##   [186]     chr4   24797140-24797177      *
##   [187]     chr6   28984234-28984270      *
##   -------
##   seqinfo: 22 sequences from an unspecified genome; no seqlengths

2.6.2.2 Save overlapping results

overlapping.results <- linear_sig[
  queryHits(findOverlaps(linear_sig_input_gr,combp_sig_gr)),
  ] # nrow: 146

overlapping.combp <- combp_sig[
  subjectHits(findOverlaps(linear_sig_input_gr,combp_sig_gr)),
  ]
colnames(overlapping.combp) <- paste0("combp_", colnames(overlapping.combp))

overlapping.results <- cbind(
  overlapping.results, overlapping.combp
)

overlapping.results.unique <- overlapping.results %>%
  group_by(inputRegion) %>%
  filter(combp_z_sidak_p == min(combp_z_sidak_p))


write.csv(
  overlapping.results.unique,
  paste0(dir.result.comp, "meta_analysis_ov_comb_p.csv"),
  row.names = FALSE
)


write.csv(
  overlapping.results.unique %>% dplyr::filter(smoke_bi == 0 & crossHyb_bi == 0),
  paste0(dir.result.comp, "meta_analysis_no_crossHyb_smoking_ov_comb_p.csv"),
  row.names = FALSE
)


linear_sig_input_no_crossHyb_smoking_gr <- linear_sig %>% 
  dplyr::filter(smoke_bi == 0 & crossHyb_bi == 0) %>% 
  tidyr::separate("inputRegion",c("chrom","start","end")) %>%
  makeGRangesFromDataFrame()
linear_sig_input_no_crossHyb_smoking_gr

## GRanges object with 421 ranges and 0 metadata columns:
##         seqnames              ranges strand
##            <Rle>           <IRanges>  <Rle>
##     [1]    chr19   49220102-49220485      *
##     [2]     chr7   27153580-27153944      *
##     [3]     chr7   27146237-27146445      *
##     [4]     chr7   27154262-27155548      *
##     [5]     chr7   27179161-27179432      *
##     ...      ...                 ...    ...
##   [417]    chr16   87886871-87886933      *
##   [418]    chr13 113422497-113422671      *
##   [419]    chr12 121476549-121476792      *
##   [420]    chr19       836606-836898      *
##   [421]     chr1   22978556-22978660      *
##   -------
##   seqinfo: 22 sequences from an unspecified genome; no seqlengths

2.6.2.3 Create Venn Diagram

library(ChIPpeakAnno)
library(ggplot2)

methods      <- c("linear", "combp")
methodsLabel <- c("coMethDMR linear", "combp")


n <- length(methods)

gg_color_hue <- function(n){
  hues = seq(15, 375, length = n + 1)
  hcl(h = hues, l = 65, c = 100)[1:n]
}

cols = gg_color_hue(n)

ranges.list <- list(linear_sig_input_no_crossHyb_smoking_gr, combp_sig_gr) 

# overlap for 1 repetition
p <- makeVennDiagram(
  ranges.list,
  NameOfPeaks = methodsLabel[1:n],
  totalTest = 37159, 
  by = "region",
  main = paste0(" "),
  col = c("#440154ff", '#21908dff'),
  fill = c(alpha("#440154ff",0.3), alpha('#21908dff',0.3)),
  cex = 1,
  fontfamily = "sans",
  cat.cex = 1,
  cat.default.pos = "outer",
  cat.pos = c(180, 180),
  cat.dist = c(-0.055, -0.055),
  cat.fontfamily = "sans",
  cat.col = c("#440154ff", '#21908dff'))

print(p)

## $p.value
##      coMethDMR.linear combp          pval
## [1,]                1     1 2.491331e-189
## 
## $vennCounts
##      coMethDMR.linear combp Counts
## [1,]                0     0  36672
## [2,]                0     1     66
## [3,]                1     0    302
## [4,]                1     1    119
## attr(,"class")
## [1] "VennCounts"

pdf(
  file = paste0(dir.result.comp,"venn_coMethDMR_linear_vs_combp.pdf"),
  width = 5, height = 5
) 
makeVennDiagram(
  ranges.list,
  NameOfPeaks = methodsLabel[1:n],
  totalTest = 37159, 
  by = "region",
  main = paste0(" "),
  col = c("#440154ff", '#21908dff'),
  fill = c(alpha("#440154ff",0.3), alpha('#21908dff',0.3)),
  cex = 1,
  fontfamily = "sans",
  cat.cex = 1,
  cat.default.pos = "outer",
  cat.pos = c(180, 180),
  cat.dist = c(-0.055, -0.055),
  cat.fontfamily = "sans",
  cat.col = c("#440154ff", '#21908dff'))
dev.off()

3 Meta-analysis Single-cpgs

3.1 Import datasets and pre-process for each cohort

library(dplyr)
dir.create(data.final,recursive = TRUE,showWarnings = FALSE)
results.files <- dir("DATASETS/",
                     pattern = "single_cpg_pVal_df.csv",
                     recursive = TRUE,
                     full.names = TRUE,
                     ignore.case = TRUE)

for(i in results.files){
  data <- readr::read_csv(i)
  dataset <- unlist(stringr::str_split(i,"//|/"))[2]  %>% as.character()
  aux <- paste0(dataset,c("_estimate", "_se", "_pValue", "_fdr"))
  colnames(data) <- c("cpg", aux)
  assign(dataset,data)
}

3.2 Create a merged final dataset

cohort_ls <- list(
  Gasparoni = GASPARONI,
  London_PFC = LONDON,
  MtSinai = MtSinai,
  ROSMAP = ROSMAP
)

### outer join input region
multi_cohorts <- Reduce(
  function(x,y, ...) merge(x, y, by = "cpg", all = TRUE, ...),
  cohort_ls
) #dim: 437713 17
dim(multi_cohorts)

3.3 Meta analysis

### calculate meta analysis z scores and p values
library(meta)

doParallel::registerDoParallel(cores = parallel::detectCores()/2)
meta_df <- plyr::adply(.data = multi_cohorts, 
                       .margins = 1, 
                       .fun =  function(rowOne_df){
                         
                         est <- rowOne_df[grep("estimate",colnames(rowOne_df))] %>% as.numeric
                         
                         direction <-  paste(ifelse(
                           is.na(est), ".",
                           ifelse(est > 0, "+", "-")
                         ),collapse = "")
                         
                         se <- rowOne_df[grep("se",colnames(rowOne_df))] %>% as.numeric
                         cohort <- gsub("_se","",grep("se",colnames(rowOne_df),value = T))
                         rowOne_reform_df <- data.frame(
                           cohort = cohort,
                           est = est,
                           se = se,
                           stringsAsFactors = FALSE
                         )
                         
                         f <- metagen(
                           TE = est,
                           seTE = se,
                           data = rowOne_reform_df
                         )
                         
                         tibble::tibble(
                           cpg = rowOne_df$cpg,
                           estimate = f$TE.fixed,
                           se = f$seTE.fixed,
                           pVal.fixed = f$pval.fixed,
                           pVal.random = f$pval.random,
                           pValQ = f$pval.Q,
                           direction = direction
                         )
                       }  , .progress = "time",
                       .parallel = TRUE,
                       .id = NULL
)

### create final pVal
meta_df$pVal.final <- ifelse(
  meta_df$pValQ > 0.05, meta_df$pVal.fixed, meta_df$pVal.random
)

### calculate fdr
meta_df$fdr <- p.adjust(meta_df$pVal.final, method = "fdr")

### order meta_df
meta_final_df <- meta_df[, c(grep("_",colnames(meta_df),invert = T),
                             grep("_",colnames(meta_df),invert = F))]
meta_final_ordered_df <- meta_final_df[order(meta_final_df$pVal.final),]

3.4 Add annotation to input cpgs

library(sesame)
probes.info <- sesameDataGet("HM450.hg19.manifest")
probes.info <- probes.info[meta_final_ordered_df$cpg %>% as.character()] %>% as.data.frame %>% dplyr::select(c("seqnames","start","end"))

result_annot_df <- merge(
  y = meta_final_ordered_df,
  x = probes.info,
  by.y = "cpg",
  by.x = "row.names",
  all.y = TRUE,
  sort = FALSE
)

colnames(result_annot_df)[1] <- "cpg"
### final raw data
write.csv(
  result_annot_df  %>% as.data.frame(),
  paste0(data.final, "meta_analysis_single_cpg_df.csv"),
  row.names = FALSE
)

3.5 Delete cross-hybridizing and smoking probes from sig probes

### Exclude non-significant probes
result_annot_sig_df <- result_annot_df %>% filter(fdr < 0.05) #dim:3979 24

write.csv(
    result_annot_sig_df %>% as.data.frame(),
    paste0(data.final, "meta_analysis_single_cpg_sig_df.csv"),
    row.names = FALSE
)

### Get crosshybrdizing probes
library(ExperimentHub)

eh = ExperimentHub()
query(eh, "DMRcate")
crosshyb <- eh[["EH3129"]]


### Get significant smoking probes
smoking.file <- "https://www.ncbi.nlm.nih.gov/pmc/articles/PMC5267325/bin/NIHMS817273-supplement-001506_-_Supplemental_Tables.xlsx"

if(!file.exists(basename(smoking.file))) downloader::download(smoking.file,basename(smoking.file))

smoking <- readxl::read_xlsx(
    basename(smoking.file),
    sheet = "02",
    skip = 2
)
smoking.sig.probes <- smoking %>% dplyr::filter(`P-value` < 1*10^(-7)) %>% pull("Probe ID")


### Exclude cross-hybridizing and smoking probes
result_annot_sig_df$cpg <- as.character(result_annot_sig_df$cpg)

result_annot_sig_no_crossHyb_smoking_df <- result_annot_sig_df[
    !(result_annot_sig_df$cpg %in% c(crosshyb, smoking.sig.probes)),
] #dim: 3751 28

### Add annotation
library(coMethDMR)

result_annot_sig_no_crossHyb_smoking_df$chrom <- as.character(
    result_annot_sig_no_crossHyb_smoking_df$seqnames
)

result_final <- AnnotateResults(result_annot_sig_no_crossHyb_smoking_df)

result_final_ordered <- result_final[
    ,c(
        1:4, 30:33, 5:12,
        grep("GASPARONI", colnames(result_final), ignore.case = TRUE, invert = FALSE),
        grep("LONDON", colnames(result_final), ignore.case = TRUE, invert = FALSE),
        grep("MTSINAI", colnames(result_final), ignore.case = TRUE, invert = FALSE),
        grep("ROSMAP", colnames(result_final), ignore.case = TRUE, invert = FALSE)
    )
]

write.csv(
    result_final_ordered %>% as.data.frame(),
    paste0(data.final, "meta_analysis_single_cpg_sig_no_crossHyb_smoking_df.csv"),
    row.names = FALSE
)

readr::read_csv(
  paste0(data.final, "meta_analysis_single_cpg_sig_no_crossHyb_smoking_df.csv"),
  col_types = readr::cols()
)

4 Results from single CPG analysis and DMR analysis

readr::read_csv(
  paste0(dir.result.cpg.vs.dmr, "meta_analysis_sig_no_crossHyb_smoking_ov_comb_p_with_sig_single_cpgs.csv"),
  col_types = readr::cols()
)

4.1 Venn plot of single cpg sig vs. probes in DMRs

overlapping.results <-  readr::read_csv(paste0(dir.result.comp, "meta_analysis_no_crossHyb_smoking_ov_comb_p.csv"),
                                        col_types = readr::cols())
probes.in.all.regions <- coMethDMR::GetCpGsInAllRegion(overlapping.results$inputRegion)
probes.in.all.regions <- probes.in.all.regions %>% unique %>% unlist 

single.cpg.results <- readr::read_csv(
  "meta_analysis_single_cpg_results/meta_analysis_single_cpg_sig_no_crossHyb_smoking_df.csv",
  col_types = readr::cols()
) 

single.cpg.results.probes <- single.cpg.results %>% pull(cpg) %>% as.character 

# library
library(VennDiagram)

## Loading required package: grid

## Loading required package: futile.logger

library(ggplot2)

# Make the plot
venn <- venn.diagram(
  x = list(
    probes.in.all.regions %>% unique ,    
    single.cpg.results.probes %>% unique  
  ), 
  category.names = c("DMR analysis probes" , "Single cpg analysis probes" ),
  filename = file.path(dir.result.cpg.vs.dmr, '/venn_DMR_cpg.png'),
  output = TRUE ,
  imagetype = "png" ,
  height = 700 ,
  width = 700 ,
  resolution = 300,
  compression = "lzw",
  lwd = 1,
  col = c("#440154ff", '#21908dff'),
  fill = c(alpha("#440154ff",0.3), alpha('#21908dff',0.3)),
  cex = 0.5,
  cat.cex = 0.3,
  cat.default.pos = "outer",
  cat.pos = c(0, 180),
  cat.dist = c(0.01, 0.01)
)

 venn.diagram(
  x = list(
    probes.in.all.regions %>% unique ,    
    single.cpg.results.probes %>% unique  
  ), 
  category.names = c("DMR analysis probes" , "Single cpg analysis probes" ),
  filename = NULL,
  output = TRUE ,
  imagetype = "png" ,
  height = 700 ,
  width = 700 ,
  resolution = 300,
  compression = "lzw",
  lwd = 1,
  col = c("#440154ff", '#21908dff'),
  fill = c(alpha("#440154ff",0.3), alpha('#21908dff',0.3)),
  cex = 1,
  cat.cex = 1,
  cat.default.pos = "outer",
  cat.pos = c(0, 180),
  cat.dist = c(0.01, 0.01)
) %>% grid.draw()

5 Session information

devtools::session_info()

## ─ Session info ───────────────────────────────────────────────────────────────
##  setting  value                                             
##  version  R Under development (unstable) (2020-02-25 r77857)
##  os       macOS Catalina 10.15.3                            
##  system   x86_64, darwin15.6.0                              
##  ui       X11                                               
##  language (EN)                                              
##  collate  en_US.UTF-8                                       
##  ctype    en_US.UTF-8                                       
##  tz       America/New_York                                  
##  date     2020-04-26                                        
## 
## ─ Packages ───────────────────────────────────────────────────────────────────
##  package                                       * version     date       lib
##  acepack                                         1.4.1       2016-10-29 [1]
##  ade4                                            1.7-15      2020-02-13 [1]
##  annotate                                        1.65.1      2020-01-27 [1]
##  AnnotationDbi                                   1.49.1      2020-01-25 [1]
##  AnnotationFilter                                1.11.0      2019-11-06 [1]
##  AnnotationHub                                 * 2.19.11     2020-04-16 [1]
##  askpass                                         1.1         2019-01-13 [1]
##  assertthat                                      0.2.1       2019-03-21 [1]
##  backports                                       1.1.6       2020-04-05 [1]
##  base64                                          2.0         2016-05-10 [1]
##  base64enc                                       0.1-3       2015-07-28 [1]
##  beanplot                                        1.2         2014-09-19 [1]
##  Biobase                                         2.47.3      2020-03-16 [1]
##  BiocFileCache                                 * 1.11.6      2020-04-16 [1]
##  BiocGenerics                                  * 0.33.3      2020-03-23 [1]
##  BiocManager                                     1.30.10     2019-11-16 [1]
##  BiocParallel                                    1.21.2      2019-12-21 [1]
##  BiocVersion                                     3.11.1      2019-11-13 [1]
##  biomaRt                                         2.43.5      2020-04-02 [1]
##  Biostrings                                    * 2.55.7      2020-03-24 [1]
##  biovizBase                                      1.35.1      2019-12-03 [1]
##  bit                                             1.1-15.2    2020-02-10 [1]
##  bit64                                           0.9-7       2017-05-08 [1]
##  bitops                                          1.0-6       2013-08-17 [1]
##  blob                                            1.2.1       2020-01-20 [1]
##  boot                                            1.3-24      2019-12-20 [1]
##  BSgenome                                        1.55.4      2020-03-19 [1]
##  bumphunter                                      1.29.0      2019-11-07 [1]
##  callr                                           3.4.3       2020-03-28 [1]
##  cellranger                                      1.1.0       2016-07-27 [1]
##  checkmate                                       2.0.0       2020-02-06 [1]
##  ChIPpeakAnno                                  * 3.21.6      2020-03-13 [1]
##  cli                                             2.0.2       2020-02-28 [1]
##  cluster                                         2.1.0       2019-06-19 [1]
##  codetools                                       0.2-16      2018-12-24 [1]
##  colorspace                                      1.4-1       2019-03-18 [1]
##  coMethDMR                                     * 0.0.0.9001  2020-03-24 [1]
##  crayon                                          1.3.4       2017-09-16 [1]
##  crosstalk                                       1.1.0.1     2020-03-13 [1]
##  curl                                            4.3         2019-12-02 [1]
##  data.table                                      1.12.9      2020-02-26 [1]
##  DBI                                             1.1.0       2019-12-15 [1]
##  dbplyr                                        * 1.4.3       2020-04-19 [1]
##  DelayedArray                                    0.13.12     2020-04-10 [1]
##  DelayedMatrixStats                              1.9.1       2020-03-30 [1]
##  desc                                            1.2.0       2018-05-01 [1]
##  devtools                                        2.3.0       2020-04-10 [1]
##  dichromat                                       2.0-0       2013-01-24 [1]
##  digest                                          0.6.25      2020-02-23 [1]
##  DMRcatedata                                   * 2.5.0       2019-10-31 [1]
##  doRNG                                           1.8.2       2020-01-27 [1]
##  dplyr                                         * 0.8.99.9002 2020-04-02 [1]
##  DT                                              0.13        2020-03-23 [1]
##  ellipsis                                        0.3.0       2019-09-20 [1]
##  ensembldb                                       2.11.4      2020-04-17 [1]
##  evaluate                                        0.14        2019-05-28 [1]
##  ExperimentHub                                 * 1.13.7      2020-04-16 [1]
##  fansi                                           0.4.1       2020-01-08 [1]
##  farver                                          2.0.3       2020-01-16 [1]
##  fastmap                                         1.0.1       2019-10-08 [1]
##  foreach                                         1.5.0       2020-03-30 [1]
##  foreign                                         0.8-78      2020-04-13 [1]
##  formatR                                         1.7         2019-06-11 [1]
##  Formula                                         1.2-3       2018-05-03 [1]
##  fs                                              1.4.1       2020-04-04 [1]
##  futile.logger                                 * 1.4.3       2016-07-10 [1]
##  futile.options                                  1.0.1       2018-04-20 [1]
##  genefilter                                      1.69.0      2019-11-06 [1]
##  generics                                        0.0.2       2018-11-29 [1]
##  GenomeInfoDb                                  * 1.23.17     2020-04-13 [1]
##  GenomeInfoDbData                                1.2.3       2020-04-20 [1]
##  GenomicAlignments                               1.23.2      2020-03-24 [1]
##  GenomicFeatures                                 1.39.7      2020-03-19 [1]
##  GenomicRanges                                 * 1.39.3      2020-04-08 [1]
##  GEOquery                                        2.55.1      2019-11-18 [1]
##  ggplot2                                       * 3.3.0       2020-03-05 [1]
##  ggpubr                                          0.2.5       2020-02-13 [1]
##  ggsignif                                        0.6.0       2019-08-08 [1]
##  glue                                            1.4.0       2020-04-03 [1]
##  GO.db                                           3.10.0      2020-03-30 [1]
##  graph                                           1.65.3      2020-04-14 [1]
##  gridExtra                                       2.3         2017-09-09 [1]
##  gtable                                          0.3.0       2019-03-25 [1]
##  Gviz                                            1.31.12     2020-03-05 [1]
##  HDF5Array                                       1.15.18     2020-04-10 [1]
##  Hmisc                                           4.4-0       2020-03-23 [1]
##  hms                                             0.5.3       2020-01-08 [1]
##  htmlTable                                       1.13.3      2019-12-04 [1]
##  htmltools                                       0.4.0       2019-10-04 [1]
##  htmlwidgets                                     1.5.1       2019-10-08 [1]
##  httpuv                                          1.5.2       2019-09-11 [1]
##  httr                                            1.4.1       2019-08-05 [1]
##  idr                                             1.2         2014-09-04 [1]
##  IlluminaHumanMethylation450kanno.ilmn12.hg19    0.6.0       2020-03-24 [1]
##  IlluminaHumanMethylationEPICanno.ilm10b2.hg19   0.6.0       2020-03-24 [1]
##  IlluminaHumanMethylationEPICanno.ilm10b4.hg19   0.6.0       2020-04-09 [1]
##  illuminaio                                      0.29.0      2019-11-06 [1]
##  interactiveDisplayBase                          1.25.0      2019-11-06 [1]
##  IRanges                                       * 2.21.8      2020-03-25 [1]
##  iterators                                       1.0.12      2019-07-26 [1]
##  jpeg                                            0.1-8.1     2019-10-24 [1]
##  jsonlite                                        1.6.1       2020-02-02 [1]
##  knitr                                           1.28        2020-02-06 [1]
##  lambda.r                                        1.2.4       2019-09-18 [1]
##  later                                           1.0.0       2019-10-04 [1]
##  lattice                                         0.20-41     2020-04-02 [1]
##  latticeExtra                                    0.6-29      2019-12-19 [1]
##  lazyeval                                        0.2.2       2019-03-15 [1]
##  lifecycle                                       0.2.0       2020-03-06 [1]
##  limma                                           3.43.8      2020-04-14 [1]
##  lme4                                            1.1-23      2020-04-07 [1]
##  lmerTest                                        3.1-2       2020-04-08 [1]
##  locfit                                          1.5-9.4     2020-03-25 [1]
##  magrittr                                        1.5         2014-11-22 [1]
##  MASS                                            7.3-51.5    2019-12-20 [1]
##  Matrix                                          1.2-18      2019-11-27 [1]
##  matrixStats                                     0.56.0      2020-03-13 [1]
##  mclust                                          5.4.6       2020-04-11 [1]
##  memoise                                         1.1.0       2017-04-21 [1]
##  mime                                            0.9         2020-02-04 [1]
##  minfi                                           1.33.1      2020-03-05 [1]
##  minqa                                           1.2.4       2014-10-09 [1]
##  multtest                                        2.43.1      2020-03-12 [1]
##  munsell                                         0.5.0       2018-06-12 [1]
##  nlme                                            3.1-147     2020-04-13 [1]
##  nloptr                                          1.2.2.1     2020-03-11 [1]
##  nnet                                            7.3-13      2020-02-25 [1]
##  nor1mix                                         1.3-0       2019-06-13 [1]
##  numDeriv                                        2016.8-1.1  2019-06-06 [1]
##  openssl                                         1.4.1       2019-07-18 [1]
##  pillar                                          1.4.3       2019-12-20 [1]
##  pkgbuild                                        1.0.6       2019-10-09 [1]
##  pkgconfig                                       2.0.3       2019-09-22 [1]
##  pkgload                                         1.0.2       2018-10-29 [1]
##  plyr                                            1.8.6       2020-03-03 [1]
##  png                                             0.1-7       2013-12-03 [1]
##  preprocessCore                                  1.49.2      2020-02-01 [1]
##  prettyunits                                     1.1.1       2020-01-24 [1]
##  processx                                        3.4.2       2020-02-09 [1]
##  progress                                        1.2.2       2019-05-16 [1]
##  promises                                        1.1.0       2019-10-04 [1]
##  ProtGenerics                                    1.19.3      2019-12-25 [1]
##  ps                                              1.3.2       2020-02-13 [1]
##  purrr                                           0.3.4       2020-04-17 [1]
##  quadprog                                        1.5-8       2019-11-20 [1]
##  R6                                              2.4.1       2019-11-12 [1]
##  rappdirs                                        0.3.1       2016-03-28 [1]
##  RBGL                                            1.63.1      2019-11-06 [1]
##  RColorBrewer                                    1.1-2       2014-12-07 [1]
##  Rcpp                                            1.0.4.6     2020-04-09 [1]
##  RcppZiggurat                                    0.1.5       2018-06-10 [1]
##  RCurl                                           1.98-1.2    2020-04-18 [1]
##  readr                                           1.3.1       2018-12-21 [1]
##  readxl                                          1.3.1       2019-03-13 [1]
##  regioneR                                        1.19.2      2020-03-28 [1]
##  remotes                                         2.1.1       2020-02-15 [1]
##  reshape                                         0.8.8       2018-10-23 [1]
##  reshape2                                        1.4.4       2020-04-09 [1]
##  Rfast                                           1.9.9       2020-03-08 [1]
##  rhdf5                                           2.31.10     2020-04-02 [1]
##  Rhdf5lib                                        1.9.3       2020-04-15 [1]
##  rlang                                           0.4.5.9000  2020-03-20 [1]
##  rmarkdown                                       2.1         2020-01-20 [1]
##  rngtools                                        1.5         2020-01-23 [1]
##  rpart                                           4.1-15      2019-04-12 [1]
##  rprojroot                                       1.3-2       2018-01-03 [1]
##  Rsamtools                                       2.3.7       2020-03-18 [1]
##  RSQLite                                         2.2.0       2020-01-07 [1]
##  rstudioapi                                      0.11        2020-02-07 [1]
##  rtracklayer                                     1.47.0      2019-11-06 [1]
##  S4Vectors                                     * 0.25.15     2020-04-04 [1]
##  scales                                          1.1.0       2019-11-18 [1]
##  scrime                                          1.3.5       2018-12-01 [1]
##  seqinr                                          3.6-1       2019-09-07 [1]
##  sessioninfo                                     1.1.1       2018-11-05 [1]
##  shiny                                           1.4.0.2     2020-03-13 [1]
##  siggenes                                        1.61.0      2019-11-06 [1]
##  statmod                                         1.4.34      2020-02-17 [1]
##  stringi                                         1.4.6       2020-02-17 [1]
##  stringr                                         1.4.0       2019-02-10 [1]
##  SummarizedExperiment                            1.17.5      2020-03-27 [1]
##  survival                                        3.1-12      2020-04-10 [1]
##  testthat                                        2.3.2       2020-03-02 [1]
##  tibble                                          3.0.1       2020-04-20 [1]
##  tidyr                                           1.0.2       2020-01-24 [1]
##  tidyselect                                      1.0.0       2020-01-27 [1]
##  usethis                                         1.6.0       2020-04-09 [1]
##  VariantAnnotation                               1.33.4      2020-04-09 [1]
##  vctrs                                           0.2.99.9010 2020-04-02 [1]
##  VennDiagram                                   * 1.6.20      2018-03-28 [1]
##  withr                                           2.1.2       2018-03-15 [1]
##  xfun                                            0.13        2020-04-13 [1]
##  XML                                             3.99-0.3    2020-01-20 [1]
##  xml2                                            1.3.1       2020-04-09 [1]
##  xtable                                          1.8-4       2019-04-21 [1]
##  XVector                                       * 0.27.2      2020-03-24 [1]
##  yaml                                            2.2.1       2020-02-01 [1]
##  zlibbioc                                        1.33.1      2020-01-24 [1]
##  source                                     
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  Bioconductor                               
##  Bioconductor                               
##  Bioconductor                               
##  Bioconductor                               
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  Bioconductor                               
##  Bioconductor                               
##  Bioconductor                               
##  CRAN (R 4.0.0)                             
##  Bioconductor                               
##  Bioconductor                               
##  Bioconductor                               
##  Bioconductor                               
##  Bioconductor                               
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  Bioconductor                               
##  Bioconductor                               
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  Bioconductor                               
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  Bioconductor                               
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  local                                      
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  Bioconductor                               
##  Bioconductor                               
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  Bioconductor                               
##  CRAN (R 4.0.0)                             
##  Github (tidyverse/dplyr@affb977)           
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  Bioconductor                               
##  CRAN (R 4.0.0)                             
##  Bioconductor                               
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  Bioconductor                               
##  CRAN (R 4.0.0)                             
##  Bioconductor                               
##  Bioconductor                               
##  Bioconductor                               
##  Bioconductor                               
##  Github (Bioconductor/GenomicRanges@70e6e69)
##  Bioconductor                               
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  Bioconductor                               
##  Bioconductor                               
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  Bioconductor                               
##  Bioconductor                               
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  Bioconductor                               
##  Bioconductor                               
##  Bioconductor                               
##  Bioconductor                               
##  Bioconductor                               
##  Bioconductor                               
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  Bioconductor                               
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  Bioconductor                               
##  CRAN (R 4.0.0)                             
##  Bioconductor                               
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  Bioconductor                               
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  Bioconductor                               
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  Bioconductor                               
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  Bioconductor                               
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  Bioconductor                               
##  Bioconductor                               
##  Github (r-lib/rlang@a90b04b)               
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  Bioconductor                               
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  Bioconductor                               
##  Bioconductor                               
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  Bioconductor                               
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  Bioconductor                               
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  Bioconductor                               
##  Github (r-lib/vctrs@fd24927)               
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  CRAN (R 4.0.0)                             
##  Bioconductor                               
##  CRAN (R 4.0.0)                             
##  Bioconductor                               
## 
## [1] /Library/Frameworks/R.framework/Versions/4.0/Resources/library