Detecting Differentially Methylated Regions with coMethDMR
2022-07-23
1 Overview
We present the workflow to analyse DNA methylation data to detect co-methylated genomics regions with differential methylation. There are two main components to the workflow:
- Detecting genomic regions of co-methylation, and
- Testing if methylation in these co-methylated regions are related to the phenotype of interest.
We will use the Tidyverse package suite for data wrangling.
You will need the following packages to complete the exercises in this example:
# Data packages; install but don't load
# BiocManager::install("IlluminaHumanMethylationEPICanno.ilm10b4.hg19")
# BiocManager::install("sesameData")
# Packages to load for this workshop
library(corrplot)
library(coMethDMR)
library(pathwayPCA)
library(kableExtra)
library(tidyverse)2 Data Inputs and Setup
2.1 Phenotype Data Structure
An emerging field of research suggests that many disease conditions and environmental factors may accelerate cellular age (AKA epigenetic age) relative to chronological age. Epigenetic aging is thought to play a direct role in age-related illnesses and its biomarkers alone can significantly predict mortality. One of the most valid and versatile epigenetic biomarkers compares predicted epigenetic age (based on quantification of DNA methylation) to actual chronological age in order to observe evidence of accelerated epigenetic aging. The authors sought to investigate the association between long-term illicit opioid use and accelerated DNA methylation among individuals of African ancestry who use heroin.
pheno_df <-
read_csv("data/pheno_clean_20220225.csv") %>%
# the coMethDMR package expects the sample ID column to be named "Sample"
rename(Sample = SampleID) %>%
# the linear mixed model (as currently coded) cannot take non-numeric
# phenotypes; this is a limitation we are working on resolving. See Issue 17
# https://github.com/TransBioInfoLab/coMethDMR/issues/17
mutate(is_case = is_case + 0) %>%
# We have reason to believe that the age effect of may be non-linear for
# certain regions of the genome: https://doi.org/10.1038/s41598-021-88504-0
mutate(Age_sq = Age ^ 2) %>%
select(Sample, Tissue, starts_with("is"), everything())
pheno_df[1:5, ]2.2 DNA Methylation Data Structure
DNA isolated from peripheral leukocytes was obtained from male and female participants of African ancestry who confirmed heroin as their drug of choice. Venous blood was collected in 8.5 ml ACD vacutainer tubes. Within 48 hrs of their collection, blood samples were transferred to Columbia University’s Human Genetics Research Core where DNA was isolated and stored at -20 C. Once extracted, DNA samples were checked for purity at an OD280 ratio of between 1.8 - 2.0. Methylation analysis was performed with the Illumina Infinium Human Methylation EPIC BeadChip that measures bisulfite-conversion-based, single-CpG resolution DNAm levels at 866,836 CpG sites in the human genome. The standard protocol of Illumina methylation assays quantifies methylation levels by the beta value using the ratio of intensities between methylated (signal A) and un-methylated (signal B) alleles. Specifically, the beta value is calculated from the intensity of the methylated (M corresponding to signal A) and un-methylated (U corresponding to signal B) alleles, as the ratio of fluorescent signals:
\[ \text{beta} := \frac{Max(M,0)}{Max(M,0)+Max(U,0)+100}. \]
Thus, beta values range from 0 (completely un-methylated) to 1 (completely methylated) (Weinhold et al., 2016)[https://www.doi.org/10.1186/s12859-016-1347-4].
Due to computational and time constraints, we will only show examples from chromosome 16. Also, because the real data is currently embargoed by NIDA, we show the real data plus small synthetic random perturbations (while preserving correlational structure between the methylation probes).
DNAm_df <-
read_csv("data/dnaM_chr16_synthetic_20220716.csv") %>%
# to preserve row names; note that many Tidyverse functions do not support row
# names for methylation data (and data frames in general)
as.data.frame() %>%
column_to_rownames(var = "Sample")
DNAm_df[1:5, 1:5]2.3 Pre-Calculated List of Regions with Contiguous CpGs
The coMethDMR package include pre-computed genomic regions that include CpGs located close to each other. These regions are in different files. The “Gene_3_200” file include genic regions with at least 3 CpG probes, with the probes located within 200 bp to each other.
Again, for computational concerns, we will only use data from the genic regions of chromosome 16 in this workshop’s example, dropping the other chromosomes from our list of close-by probes.
closeByGeneAll_ls <- readRDS("data/EPIC_10b4_Gene_3_200.rds")
chr16_lgl <- str_detect(names(closeByGeneAll_ls), "chr16")
# table(chr16_lgl)
closeByGeneChr16_ls <- closeByGeneAll_ls[chr16_lgl]
rm(closeByGeneAll_ls, chr16_lgl)
closeByGeneChr16_ls[1:5]## $`chr16:6068831-6069198`
## [1] "cg03586879" "cg19378133" "cg02230017" "cg08545287" "cg03934713"
##
## $`chr16:6069731-6070292`
## [1] "cg02119348" "cg03072621" "cg08582701" "cg08015447" "cg05107500"
##
## $`chr16:6533549-6533789`
## [1] "cg03169041" "cg27256528" "cg12310850" "cg00853103"
##
## $`chr16:6623995-6624293`
## [1] "cg16427239" "cg01914555" "cg06560026" "cg07891288"
##
## $`chr16:7354254-7354432`
## [1] "cg00590139" "cg00982312" "cg00470955"2.4 Example: Visualizing Methylation in a Single Region
2.4.1 Methylation Data for the Region
We pick the first region in this list and inspect the region IDs and the corresponding data.
# Region Name
testRegionName_char <- "chr16:6068831-6069198"
# Region
testRegion_char <- closeByGeneChr16_ls[[testRegionName_char]]
testRegion_char## [1] "cg03586879" "cg19378133" "cg02230017" "cg08545287" "cg03934713"# DNA Methylation Data for that Region (5 Samples)
DNAm_df[testRegion_char, 1:5]2.4.2 A Boxplot of Methylation by Probe
What does the methylation look like for these probes? We first pivot the data so that each row is one methylation value.
testMethPivot_df <-
DNAm_df[testRegion_char, ] %>%
pathwayPCA::TransposeAssay(omeNames = "rowNames") %>%
rownames_to_column("Sample") %>%
pivot_longer(!Sample, names_to = "CpG") %>%
arrange(Sample)Now we plot the methylation beta values in this example region for these participants.
ggplot(data = testMethPivot_df) +
aes(x = CpG, y = value) +
labs(
title = "Unadjusted DNA Methylation Beta Values for 96 Participants",
subtitle = testRegionName_char,
y = "Methylation Beta Value"
) +
scale_y_continuous(limits = c(0, 1)) +
geom_boxplot(position = "dodge")
3 Detect Regions of Co-Methylation with CoMethAllRegions()
We now have an example region, but we would like to know if the probes in this region are co-methylated (i.e., that their methylation patterns are correlated within this region). The steps to do this are as follows:
- Regress out covariate effects from the methylation, then
- Find correlated clusters of CpGs within a region
Because this step does not use the phenotype information, it is unsupervised.
3.1 Step 1: Adjust Methylation for Covariates with GetResiduals()
From Lily: People might ask - why not include the covariates in one model?
We first find the methylation M-value residuals, then fit the following model separately for each probe \(j\):
\[ M_{i,j} = \beta_0 + \beta_1 * \text{age}_i + \beta_2 * \text{age}^2_i + \beta_3 * \text{sex}_i + \varepsilon_{i,j}. \]
We then use the residuals, \(\boldsymbol{\varepsilon}\), from these independent linear models as the inputs for the region-based clustering algorithm in Step 2.
Using data from chromosome 16 only, this code runs in 1.6 minutes on a 2016 MacBook Pro laptop with 16Gb RAM and 1.7 minutes on a 2020 Dell Latitude 5501 with 16 RAM.
system.time(
residByAgeSex_mat <- GetResiduals(
dnam = DNAm_df,
pheno_df = pheno_df,
betaToM = TRUE,
covariates_char = c("Age", "Age_sq", "is_female"),
# on my Mac, 1 core worked best; on Windows, 3 cores worked best
nCores_int = 1L
)
)
residByAgeSex_df <- as_tibble(residByAgeSex_mat, rownames = "probeID")
# saveRDS(
# residByAgeSex_df,
# file = "results/chr16_residuals_age_sex_20220716.rds"
# )If you don’t want to wait, we have these results calculated, so you can load them.
residByAgeSex_df <-
readRDS(file = "results/chr16_residuals_age_sex_20220716.rds") %>%
column_to_rownames(var = "probeID")
residByAgeSex_df[1:5, 1:5]3.2 Step 2: Detecting Regions of Co-Methylation after Removing for Age and Sex Effects
Now we can search through the pre-defined regions to find those regions which have co-methylation.
Finding all the regions of concurrent methylation in chromosome 16 takes 7-8 minutes on a 2016 MacBook Pro laptop with 16Gb RAM and 15 minutes on a 2020 Dell Latitude 5501 with 16 RAM. NOTE: in the second section of code, we are selecting the option output = "dataframe" to show how the function marks regions of co-methylation; in practice, we use the first section of code and set output = "CpGs" for use in downstream analysis.
### For Downstream Analysis ###
system.time(
coMethAgeSex_ls <- CoMethAllRegions(
dnam = residByAgeSex_df,
betaToM = FALSE,
method = "spearman",
arrayType = "EPIC",
CpGs_ls = closeByGeneChr16_ls,
output = "CpGs",
nCores_int = 3L
)
)
# saveRDS(
# coMethAgeSex_ls,
# file = "results/chr16_coMethylated_regions_age_sex_20220716.rds"
# )
### For Inspection and Troubleshooting ONLY ###
system.time(
coMethAgeSexDF_ls <- CoMethAllRegions(
dnam = residByAgeSex_df,
betaToM = FALSE,
method = "spearman",
arrayType = "EPIC",
CpGs_ls = closeByGeneChr16_ls,
output = "dataframe",
nCores_int = 3L
)
)
# # We save ONLY the first result (as an example)
# saveRDS(
# coMethAgeSexDF_ls[[1]],
# file = "results/chr16_coMethylated_region1_df_age_sex_20220716.rds"
# )As previously, we have also calculated these co-methylated regions.
# For Analysis
coMethAgeSex_ls <- readRDS(
file = "results/chr16_coMethylated_regions_age_sex_20220716.rds"
)
# Inspection ONLY
coMethAgeSex_df <- readRDS(
file = "results/chr16_coMethylated_region1_df_age_sex_20220716.rds"
)
coMethAgeSex_dfThis table of results includes the probes in this region ordered by their position on the genome. The r_drop column measures the correlation (note the method = "spearman" argument) between the methylation values for that probe and the average methylation for all the other probes in the region. In our previous simulation studies, we have found that an r_drop threshold of 0.4 yielded the best results. We see that the correlation between probe cg19378133 and the average of the other probes is below the threshold, so it makes sense to not retain that probe. However, probe cg03586879 passes our threshold, but we are still discarding it. Why? Because probe cg19378133 breaks the continuity of correlated probes across this subregion of the genome; this is what we mean by co-methylation: the probes are all next to each other on the genome and they are all correlated with each other.
We also explore this example visually. As we can see in the figure below, probes cg02230017, cg08545287, and cg03934713 all make up a subregion of contiguous and concurrent (co-) methylation.
DNAm_df[rownames(DNAm_df) %in% coMethAgeSex_df$CpG, ] %>%
pathwayPCA::TransposeAssay(omeNames = "rowNames") %>%
select(coMethAgeSex_df$CpG) %>%
cor() %>%
corrplot::corrplot(method = "number")
We now see that we only want to consider the last three probes in this region as a region of co-methylation. Therefore, the actual results list (from setting option output = "CpGs") looks like this:
coMethAgeSex_ls[1]## $`chr16:6069019-6069198`
## [1] "cg02230017" "cg08545287" "cg03934713"4 Identifying Differetially (Co-)Methylated Regions with lmmTestAllRegions()
In the previous step, we identified CpG-dense regions which have concurrent methylation. Now, we will fit separate linear mixed models to test the association between the phenotype of interest and each detected region.
4.1 The Input Data
We need the following pieces of data:
- the list of concurrently methylated regions (from
CoMethAllRegions()) - a data frame of DNA methylation beta values (not the residuals) in probe by sample form (with probe IDs as row names)
- phenotype / clinical data containing a column called
"Sample"(case sensitivity applies)
4.2 Mixed Model Types: Simple or Random Coefficient
Our method offers two mixed models to find associations between DNA methylation in a region with detected co-methylatioin and the phenotype of interest. In the models below, let \(M_{ij}\) be the DNA methylation M-value for CpG \(j\) in subject \(i\); let \(X_i\) be the phenotype for subject \(i\); let \(U_i\) be a random effect for sample \(i\); and let \(\textbf{V}_{i}\) be a vector of covariates for sample \(i\). For more details, please see Gomez, Odom, et al., (2019). In our example, \(X_{i}\) indicates if the subject is a current heroin user or a healthy control; \(U_{i}\) is a subject effect that allows the average methylation M-value (\(\bar{M}_i\)) for each subject to vary; and \(\textbf{V}_{i}\) is a vector of the biological sex, age, and square of age for subject \(i\).
4.2.1 Simple Linear Mixed Model
The simple linear mixed model assumes the same effect of \(X\) on the methylation for each subject and each CpG, while allowing the average baseline methylation to vary across subjects (via \(\beta_0 + U_{i}\)). It is defined as \[ M_{ij} = \beta_0 + \beta_1X_{i} + U_{i} + \textbf{bV}_i + \varepsilon_{ij}. \]
4.2.2 Random Coefficient Linear Mixed Model
The random coefficient linear mixed model allows the effect of \(X\) on methylation to be different for each CpG (via \(\beta_1 + \alpha_{1j})\), while simultaneously allowing for different baseline methylation for subject and CpG (via \(\beta_0 + \alpha_{0j} + U_{i}\)). It is defined as \[ M_{ij} = (\beta_0 + \alpha_{0j}) + (\beta_1 + \alpha_{1j})X_{i} + U_{i} + \textbf{bV}_i + \varepsilon_{ij}, \] where \([\alpha_{0j}, \alpha_{1j}]^T \sim N(\textbf{0}, \textbf{G})\) are random coefficients for intercepts and slopes, and \(\textbf{G}\) is an unstructured covariance matrix.
4.3 Differentially Co-Methylated Regions after Adjusting for Age and Sex
We specify the following:
betas = DNAm_df: our data frame of methylation beta valuesregion_ls = coMethAgeSex_ls: the list of regions of concurrent methylation (recall that this example only includes genic regions on chromosome 16 for computational concerns)pheno_df = pheno_df: clinical data with subject ID column (which matches the columns ofDNAm_df)contPheno_char = "is_case": the name of the column with the outcome of interest; in this example, heroin use or not. NOTE: at this time, we only support binary and continuous outcomes; we are working to add support for categorical and survival outcomes in the future.covariates_char = c("is_female", "Age", "Age_sq"): which covariates should be used as moderators of the response’s effect on DNA methylation? Often sex, age, and plate/batch effects should be included here.modelType = "randCoef": use a random coefficient mixed model, as described abovearrayType = "EPIC": which Illumina manifest matches your data? We currently support 450k and EPIC, but we are researching appropriate statistical methods for whole genome methylation analysis.outLogFile = paste0("results/lmmLog_", Sys.Date(), ".txt"): where should the log file go? This file will contain messages about the convergence of the linear mixed model’s optimization routines. If you are running a quick test, this can beNULL.nCores_int = 8L: how many cores should you use? Note that fitting these mixed models can be memory intensive, so more cores are not always better. Try a few options with smaller subsets of your data to find the “sweet spot”.
This takes 4 minutes for ~1100 regions using 8 cores on a 2016 MacBook Pro laptop with 16Gb RAM and 4.7 minutes on a 2020 Dell Latitude 5501 with 16 RAM (also using 8 cores).
system.time(
lmmHeroinOut_df <- lmmTestAllRegions(
betas = DNAm_df,
region_ls = coMethAgeSex_ls,
pheno_df = pheno_df,
contPheno_char = "is_case",
covariates_char = c("is_female", "Age", "Age_sq"),
modelType = "randCoef",
arrayType = "EPIC",
# outLogFile = paste0("results/lmmLog_", Sys.Date(), ".txt"),
nCores_int = 8L
)
)
write_csv(
x = lmmHeroinOut_df,
file = "results/chr16_coMethylated_results_heroin_use_20220716.csv"
)Again, we have calculated this output for you:
lmmHeroinOut_df <- read_csv(
file = "results/chr16_coMethylated_results_heroin_use_20220716.csv"
)
lmmHeroinOut_df5 Annotate the Significant Results with AnnotateResults() and Inspect
Now we would like to see which regions of the genome have co-methylation that is associated with heroin use, after adjusting for age and sex. Because this is synthetic data, we will use a more relaxed false-discovery rate cutoff.
annotatedRes_df <-
# Filter to rows with significant regions
lmmHeroinOut_df %>%
filter(FDR < 0.2) %>%
# AnnotateResults() uses 450k by default, so change the array type to EPIC
AnnotateResults(arrayType = "EPIC") %>%
# Wrangle
mutate(range = paste0(chrom, ":", start, "-", end)) %>%
select(-chrom, -start, -end, -UCSC_RefGene_Accession) %>%
select(range, everything())annotatedRes_df %>%
arrange(FDR) %>%
kable("html") %>%
kable_styling("striped") %>%
scroll_box(width = "100%")| range | nCpGs | Estimate | StdErr | Stat | pValue | FDR | UCSC_RefGene_Group | UCSC_RefGene_Name | Relation_to_Island |
|---|---|---|---|---|---|---|---|---|---|
| chr16:8806966-8807043 | 5 | 0.3249334 | 0.0980320 | 3.314565 | 0.0009179 | 0.1442744 | 1stExon;5’UTR | ABAT | OpenSea |
| chr16:4102293-4102366 | 3 | 0.2260326 | 0.0677940 | 3.334109 | 0.0008557 | 0.1442744 | Body | ADCY9 | OpenSea |
| chr16:367990-368066 | 3 | -0.3680423 | 0.1023106 | -3.597302 | 0.0003215 | 0.1442744 | Body | AXIN1 | OpenSea |
| chr16:66835860-66836152 | 4 | 0.1450601 | 0.0442113 | 3.281068 | 0.0010341 | 0.1442744 | TSS1500 | CCDC79 | S_Shore |
| chr16:67562141-67562437 | 3 | 0.1613868 | 0.0498772 | 3.235682 | 0.0012135 | 0.1442744 | Body;TSS1500 | FAM65A;LOC100505942 | N_Shore |
| chr16:12241554-12241816 | 3 | 0.3453184 | 0.1038393 | 3.325509 | 0.0008826 | 0.1442744 | Body | SNX29 | OpenSea |
| chr16:89600003-89600081 | 3 | 0.1419718 | 0.0435319 | 3.261329 | 0.0011089 | 0.1442744 | Body | SPG7 | N_Shelf |
| chr16:31106703-31106788 | 3 | 0.2374078 | 0.0708141 | 3.352551 | 0.0008007 | 0.1442744 | TSS1500 | VKORC1 | S_Shore |
| chr16:30543920-30544242 | 3 | -0.2234801 | 0.0656618 | -3.403503 | 0.0006653 | 0.1442744 | 3’UTR | ZNF747 | Island |
| chr16:89363811-89364225 | 7 | 0.1850613 | 0.0583705 | 3.170460 | 0.0015220 | 0.1557890 | Body | ANKRD11 | OpenSea |
| chr16:67562299-67562517 | 4 | 0.1514994 | 0.0480094 | 3.155621 | 0.0016016 | 0.1557890 | Body;TSS1500 | FAM65A;LOC100505942 | N_Shore |
| chr16:77232089-77232522 | 4 | 0.1224214 | 0.0400089 | 3.059858 | 0.0022144 | 0.1974522 | 3’UTR;Body;TSS1500 | MON1B;SYCE1L | S_Shelf |
Note that these results are based on synthetic data, so they should not be interpreted biologically. In the real data, we identified regions in the bodies of CASKIN1 and RHOT2 as having both statistically significant and potentially biological meaningful implications (both with established implications for alcohol use disorder).
6 Conclusion
There were two main components of our workflow:
- We first detected genomic regions of co-methylation, after adjusting the DNA methylation for subject age and sex.
- We tested if methylation in these co-methylated regions was related to long-term heroin use, and then mapped regions with heroin-use-associated methylation to genes.
To our knowledge, this is the first time this technique has been applied to methylation data from individuals with opioid use disorder. We have provided preliminary support for this novel investigatory methodology’s ability to identify global changes in DNA methylation among individuals who misuse addictive substances. Previous studies have suggested some of the regions we identified may be involved with some of the underlying mechanisms within the opioid system. We are currently conducting additional analysis and modeling, and a manuscript with our findings is in preparation. Future studies may be able to integrate this tool to better understand the epigenetic effects of opioid use and substnance misuse as a whole.
sessionInfo()## R version 4.2.0 (2022-04-22)
## Platform: x86_64-apple-darwin17.0 (64-bit)
## Running under: macOS Big Sur/Monterey 10.16
##
## Matrix products: default
## BLAS: /Library/Frameworks/R.framework/Versions/4.2/Resources/lib/libRblas.0.dylib
## LAPACK: /Library/Frameworks/R.framework/Versions/4.2/Resources/lib/libRlapack.dylib
##
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] forcats_0.5.1 stringr_1.4.0 dplyr_1.0.9 purrr_0.3.4
## [5] readr_2.1.2 tidyr_1.2.0 tibble_3.1.7 ggplot2_3.3.6
## [9] tidyverse_1.3.1 kableExtra_1.3.4 pathwayPCA_1.12.0 coMethDMR_1.0.0
## [13] corrplot_0.92 BiocStyle_2.24.0
##
## loaded via a namespace (and not attached):
## [1] utf8_1.2.2
## [2] tidyselect_1.1.2
## [3] lme4_1.1-29
## [4] RSQLite_2.2.14
## [5] AnnotationDbi_1.58.0
## [6] grid_4.2.0
## [7] BiocParallel_1.30.2
## [8] munsell_0.5.0
## [9] preprocessCore_1.58.0
## [10] codetools_0.2-18
## [11] withr_2.5.0
## [12] colorspace_2.0-3
## [13] Biobase_2.56.0
## [14] filelock_1.0.2
## [15] highr_0.9
## [16] knitr_1.39
## [17] rstudioapi_0.13
## [18] stats4_4.2.0
## [19] MatrixGenerics_1.8.0
## [20] labeling_0.4.2
## [21] lars_1.3
## [22] GenomeInfoDbData_1.2.8
## [23] bit64_4.0.5
## [24] farver_2.1.0
## [25] rhdf5_2.40.0
## [26] vctrs_0.4.1
## [27] generics_0.1.2
## [28] xfun_0.31
## [29] BiocFileCache_2.4.0
## [30] R6_2.5.1
## [31] GenomeInfoDb_1.32.2
## [32] illuminaio_0.38.0
## [33] locfit_1.5-9.5
## [34] bitops_1.0-7
## [35] rhdf5filters_1.8.0
## [36] cachem_1.0.6
## [37] reshape_0.8.9
## [38] DelayedArray_0.22.0
## [39] assertthat_0.2.1
## [40] promises_1.2.0.1
## [41] BiocIO_1.6.0
## [42] scales_1.2.0
## [43] vroom_1.5.7
## [44] gtable_0.3.0
## [45] rlang_1.0.2
## [46] genefilter_1.78.0
## [47] systemfonts_1.0.4
## [48] splines_4.2.0
## [49] rtracklayer_1.56.0
## [50] GEOquery_2.64.2
## [51] broom_0.8.0
## [52] BiocManager_1.30.18
## [53] yaml_2.3.5
## [54] modelr_0.1.8
## [55] GenomicFeatures_1.48.1
## [56] backports_1.4.1
## [57] httpuv_1.6.5
## [58] tools_4.2.0
## [59] bookdown_0.26
## [60] nor1mix_1.3-0
## [61] ellipsis_0.3.2
## [62] jquerylib_0.1.4
## [63] RColorBrewer_1.1-3
## [64] siggenes_1.70.0
## [65] BiocGenerics_0.42.0
## [66] Rcpp_1.0.8.3
## [67] plyr_1.8.7
## [68] sparseMatrixStats_1.8.0
## [69] progress_1.2.2
## [70] zlibbioc_1.42.0
## [71] RCurl_1.98-1.6
## [72] prettyunits_1.1.1
## [73] openssl_2.0.1
## [74] bumphunter_1.38.0
## [75] S4Vectors_0.34.0
## [76] SummarizedExperiment_1.26.1
## [77] haven_2.5.0
## [78] fs_1.5.2
## [79] magrittr_2.0.3
## [80] data.table_1.14.2
## [81] lmerTest_3.1-3
## [82] reprex_2.0.1
## [83] IlluminaHumanMethylationEPICanno.ilm10b4.hg19_0.6.0
## [84] matrixStats_0.62.0
## [85] hms_1.1.1
## [86] mime_0.12
## [87] evaluate_0.15
## [88] xtable_1.8-4
## [89] XML_3.99-0.9
## [90] mclust_5.4.10
## [91] readxl_1.4.0
## [92] IRanges_2.30.0
## [93] compiler_4.2.0
## [94] biomaRt_2.52.0
## [95] minfi_1.42.0
## [96] crayon_1.5.1
## [97] minqa_1.2.4
## [98] htmltools_0.5.2
## [99] later_1.3.0
## [100] tzdb_0.3.0
## [101] lubridate_1.8.0
## [102] DBI_1.1.2
## [103] ExperimentHub_2.4.0
## [104] dbplyr_2.1.1
## [105] MASS_7.3-57
## [106] rappdirs_0.3.3
## [107] boot_1.3-28
## [108] Matrix_1.4-1
## [109] cli_3.3.0
## [110] quadprog_1.5-8
## [111] parallel_4.2.0
## [112] GenomicRanges_1.48.0
## [113] pkgconfig_2.0.3
## [114] GenomicAlignments_1.32.0
## [115] numDeriv_2016.8-1.1
## [116] xml2_1.3.3
## [117] foreach_1.5.2
## [118] annotate_1.74.0
## [119] svglite_2.1.0
## [120] bslib_0.3.1
## [121] rngtools_1.5.2
## [122] multtest_2.52.0
## [123] beanplot_1.3.1
## [124] webshot_0.5.3
## [125] XVector_0.36.0
## [126] rvest_1.0.2
## [127] scrime_1.3.5
## [128] doRNG_1.8.2
## [129] digest_0.6.29
## [130] Biostrings_2.64.0
## [131] base64_2.0
## [132] rmarkdown_2.14
## [133] cellranger_1.1.0
## [134] DelayedMatrixStats_1.18.0
## [135] restfulr_0.0.13
## [136] curl_4.3.2
## [137] shiny_1.7.1
## [138] Rsamtools_2.12.0
## [139] rjson_0.2.21
## [140] nloptr_2.0.1
## [141] lifecycle_1.0.1
## [142] nlme_3.1-157
## [143] jsonlite_1.8.0
## [144] Rhdf5lib_1.18.2
## [145] askpass_1.1
## [146] limma_3.52.1
## [147] viridisLite_0.4.0
## [148] fansi_1.0.3
## [149] pillar_1.7.0
## [150] lattice_0.20-45
## [151] KEGGREST_1.36.0
## [152] fastmap_1.1.0
## [153] httr_1.4.3
## [154] survival_3.3-1
## [155] interactiveDisplayBase_1.34.0
## [156] glue_1.6.2
## [157] png_0.1-7
## [158] iterators_1.0.14
## [159] BiocVersion_3.15.2
## [160] bit_4.0.4
## [161] stringi_1.7.6
## [162] sass_0.4.1
## [163] HDF5Array_1.24.0
## [164] blob_1.2.3
## [165] AnnotationHub_3.4.0
## [166] memoise_2.0.1