Full Analysis with epiRomics
Alex M. Mawla & Mark O. Huising
2026-04-20
Source:vignettes/articles/full-analysis-with-epiRomics.Rmd
full-analysis-with-epiRomics.RmdAbstract
Summary epiRomics is an R package designed to integrate multi-omics data in order to identify and visualize enhancer regions alongside gene expression and other epigenomic modifications. Regulatory network analysis can be done using combinatory approaches to infer regions of significance such as enhancers, when combining ChIP and histone data. Downstream analysis can identify co-occurrence of these regions of interest with other user-supplied data, such as chromatin availability or gene expression. Finally, this package allows for results to be visualized at high resolution in a stand-alone browser.
Availability and Implementation epiRomics is released under the Artistic-2.0 License. The source code and documents are freely available through GitHub (https://github.com/Huising-Lab/epiRomics).
Contact ammawla@ucdavis.edu or mhuising@ucdavis.edu
Supplementary information Supplementary data and methods are available online on bioRxiv or GitHub.
Competing Interest Statement The authors have declared no competing interest.
Citation
If you use epiRomics in published research, please cite:
Mawla AM, van der Meulen T, Huising MO. Chromatin accessibility differences between alpha, beta, and delta cells identifies common and cell type-specific enhancers. BMC Genomics. 2023;24:202. https://doi.org/10.1186/s12864-023-09293-6
Mawla AM & Huising MO. epiRomics: a multi-omics R package to identify and visualize enhancers. bioRxiv. 2021. https://doi.org/10.1101/2021.08.19.456732
Interactive showcases
Two companion Shiny applications render published enhanceosome browsers built with epiRomics:
Shiny epiRomics Mouse Islet Enhanceosome Browser — mouse pancreatic alpha, beta, and delta cell results from Mawla et al. 2023.
Shiny epiRomics Human Islet Enhanceosome Browser — human alpha and beta cell results used as the example dataset for this package (Mawla and Huising 2021).
Installing and loading the epiRomics package
Install epiRomics and its annotation dependencies from
Bioconductor (these are declared as Suggests, so install
them explicitly the first time you render the vignette):
if (!requireNamespace("BiocManager", quietly = TRUE))
install.packages("BiocManager")
BiocManager::install(c(
"epiRomics",
"TxDb.Hsapiens.UCSC.hg38.knownGene",
"org.Hs.eg.db"
))With the package and its annotation companions installed, attach them at the top of the session:
epiRomics provides the analysis API; the hg38
TxDb supplies gene models and the org.Hs.eg.db
package resolves gene symbols.
Locating the example data
This vignette runs end-to-end against the full 1.3 GB curated archive
on Zenodo (see the Zenodo record at https://zenodo.org/records/19189987 and the package
README for curation methodology). Because the archive is too large to
ship inside a Bioconductor tarball, we download it once via
cache_data() and then resolve its on-disk location at
render time.
Run this once (interactive session) to populate the
BiocFileCache entry used by every chunk below:
epiRomics::cache_data()The block below probes four sources in priority order — the installed
package extdata/, the development source tree, the on-disk
BiocFileCache, and (only in an interactive session) an
automatic cache_data() fallback — then reports which source
was picked and whether the BigWig signal tracks are present:
## Priority 1: installed package extdata (full install or R CMD check)
extdata_dir <- system.file("extdata", package = "epiRomics")
has_data <- nzchar(extdata_dir) &&
file.exists(file.path(extdata_dir, "example_epiRomics_Db_sheet.csv")) &&
dir.exists(file.path(extdata_dir, "BigWigs"))
## Priority 2: development source tree (rendering from package root)
if (!has_data) {
dev_path <- normalizePath(file.path("..", "inst", "extdata"),
mustWork = FALSE)
if (dir.exists(dev_path) &&
file.exists(file.path(dev_path,
"example_epiRomics_Db_sheet.csv")) &&
dir.exists(file.path(dev_path, "BigWigs"))) {
extdata_dir <- dev_path
has_data <- TRUE
}
}
## Priority 3: BiocFileCache (lightweight install + cached data)
if (!has_data && isTRUE(epiRomics::has_cache())) {
cache_dir <- epiRomics::get_cache_path()
if (!is.null(cache_dir) &&
file.exists(file.path(cache_dir,
"example_epiRomics_Db_sheet.csv")) &&
dir.exists(file.path(cache_dir, "BigWigs"))) {
extdata_dir <- cache_dir
has_data <- TRUE
}
}
## Priority 4: interactive download as last resort
if (!has_data && interactive()) {
message("Downloading epiRomics example data (~1.3 GB, one-time)...")
tryCatch(
{ epiRomics::cache_data(); has_data <- TRUE
extdata_dir <- epiRomics::get_cache_path() },
error = function(e) message("Download failed: ", e$message,
"\nRun cache_data() to retry.")
)
}
if (!has_data) extdata_dir <- NULL
## BigWig plots require the signal files
has_bigwigs <- has_data && file.exists(
file.path(extdata_dir, "BigWigs",
"Kaestner_Alpha.atac.signal.subset.multi_chr.bigwig")
)
## Report a machine-independent discovery descriptor rather than the
## absolute path, so the rendered vignette is identical across users.
data_source <- if (!has_data) {
"none (run epiRomics::cache_data() in an interactive session)"
} else if (identical(extdata_dir,
system.file("extdata", package = "epiRomics"))) {
"installed package extdata/"
} else if (grepl("inst/extdata$", extdata_dir)) {
"development source tree (inst/extdata/)"
} else {
"BiocFileCache (epiRomics::cache_data())"
}
has_data
#> [1] FALSE
has_bigwigs
#> [1] FALSE
data_source
#> [1] "none (run epiRomics::cache_data() in an interactive session)"has_data and has_bigwigs drive the
per-chunk eval = guards below: chunks that read BED/CSV
files are gated on has_data; chunks that plot BigWig signal
are gated on has_bigwigs. If the archive is not yet cached,
the vignette still renders — the gated chunks simply skip — but running
epiRomics::cache_data() once in an interactive session
unlocks the full walkthrough.
Data status: Example data is not cached. Code chunks are shown but not executed.
Run epiRomics::cache_data() in an R session to download
the data, then rebuild this vignette.
Brief explanation of example data
This package includes example data delineating human pancreatic islet enhancers between alpha and beta cells.
Human pancreatic islet alpha and beta ATAC- and companion RNA-Seq data were retrieved from GEO accession GSE76268 (Ackermann et al. 2016).
ATAC samples were processed using the ENCODE-DCC ATAC sequencing pipeline, aligning to the hg38 build of the human genome (Harrow et al. 2012; ENCODE Project Consortium 2012; Davis et al. 2018).
Peak calls generated through the pipeline using MACS2 (Zhang et al. 2008) were analyzed downstream through the Bioconductor package DiffBind (Ross-Innes et al. 2012) to identify differentially enriched chromatin regions between the two cell types.
RNA samples were quality controlled using fastp (Chen et al. 2018), and aligned using STAR (Dobin et al. 2013) to hg38. Wiggle files produced by the STAR aligner were merged by cell type using UCSC command line tools.
BigWig files were subsetted using UCSC command line tools (Kent et al. 2010).
ChIP-seq peak calls for human pancreatic islet transcription factors Foxa2, MafB, Nkx2.2, Nkx6.1, and Pdx1 were retrieved from EMBL-EBI E-MTAB-1919 (Pasquali et al. 2014).
Histone modification peak calls for H3K27ac, H3K4me1, H3K27me3, H3K9me3, H3K4me3, H3K36me3 were retrieved from GEO GSE16256 (Bernstein et al. 2010), and H2A.Z from EMBL-EBI E-MTAB-1919 (Pasquali et al. 2014).
All peak calls retrieved as hg19 were lifted over to hg38 using the UCSC liftOver (Kent et al. 2002).
The FANTOM5 human enhancer database (Lizio et al. 2015), ultra-conserved non-coding elements from the UCNE database (Dimitrieva and Bucher 2013), and the human islet regulome (Miguel-Escalada et al. 2019) were also retrieved and lifted to hg38.
Building the database
The database sheet CSV defines the input annotations. Each row specifies a BED file with its name, relative path, genome build, format, and type.
Required columns: name, path,
genome, format, type.
genome is a free-form string (hg38,
mm10, rn6, …) that is cross-checked against
the TxDb you supply to build_database() via
GenomeInfoDb::genome(txdb), so a mismatch raises an error
rather than proceeding silently. Any organism available in
Bioconductor’s annotation catalogue works; the code path does not branch
on a whitelist. Type can be: histone, methyl,
SNP, chip, or functional.
example_epiRomics_Db_sheet <- read.csv(
file.path(extdata_dir,
"example_epiRomics_Db_sheet.csv"))
head(example_epiRomics_Db_sheet, 5)build_database() constructs the S4 database from the
data sheet. The CSV uses relative paths; the data_dir
parameter resolves them against the cache directory automatically.
database <- build_database(
db_file = file.path(extdata_dir,
"example_epiRomics_Db_sheet.csv"),
txdb_organism = paste0(
"TxDb.Hsapiens.UCSC.hg38.knownGene::",
"TxDb.Hsapiens.UCSC.hg38.knownGene"),
genome = "hg38",
organism = "org.Hs.eg.db",
data_dir = extdata_dir)database is the S4 epiRomicsS4 object used
by every downstream step. It holds the overlap-annotated
GRanges and metadata; all pipeline functions accept it
directly.
Quick visualization with plot_quick_view
Before running the full enhancer pipeline, epiRomics provides a zero-configuration entry point for the ultra-fast visualization of any gene locus with BigWig signal tracks. No database setup is required—just pass a gene symbol and BigWig file paths.
## Load BigWig track connection sheet and resolve paths
track_connection <- read.csv(
file.path(extdata_dir,
"example_epiRomics_BW_sheet.csv"))
track_connection$path <- file.path(
extdata_dir,
track_connection$path)
## Reorder: Beta first, then Alpha
beta_idx <- grep("Beta", track_connection$name)
alpha_idx <- grep("Alpha", track_connection$name)
track_connection <- track_connection[
c(beta_idx, alpha_idx), ]
## Named vector of BigWig paths and colors
bw_paths <- setNames(
track_connection$path,
track_connection$name)
bw_colors <- track_connection$colorThe mirror = TRUE option displays ATAC and RNA signals
in a mirrored layout for direct comparison between cell types. The
genome = "hg38" argument below is a string label that
plot_quick_view() uses to pick the matching
TxDb.* / org.*.db pair loaded in the session;
any genome label that corresponds to an installed Bioconductor
TxDb.* package (mm10, rn6,
dm6, …) works identically.
Beta cell genes
INS
The INS locus encodes insulin and is the canonical beta-cell marker. We expect sharp beta-specific ATAC accessibility and RNA signal here; mirrored against the alpha-cell tracks, the asymmetry should be visible at a glance.
plot_quick_view("INS",
bw_paths = bw_paths,
colors = bw_colors,
mirror = TRUE,
genome = "hg38")The beta-cell ATAC and RNA tracks (above the axis) clearly dominate the INS body while alpha-cell signal (below the axis) is nearly flat, confirming the expected cell-type-specific activity.
PDX1
PDX1 is the beta-cell master transcription factor. Unlike INS, much of its regulatory activity sits in distal enhancers rather than directly over the gene body, so we expect strong beta-specific accessibility spread across a broader locus.
plot_quick_view("PDX1",
bw_paths = bw_paths,
colors = bw_colors,
mirror = TRUE,
genome = "hg38")The beta-cell ATAC signal spans the gene body and its flanking regulatory region; RNA expression is beta-biased but less extreme than INS because PDX1 is expressed at lower steady-state levels.
MAFA
MAFA is a late-stage beta-cell identity factor.
plot_quick_view("MAFA",
bw_paths = bw_paths,
colors = bw_colors,
mirror = TRUE,
genome = "hg38")As with INS and PDX1, the beta-cell tracks dominate; alpha-cell signal is sparse.
Alpha cell genes
GCG
GCG encodes glucagon, the alpha-cell counterpart to insulin. The expected pattern here is a mirror image of INS: strong alpha-cell ATAC and RNA signal, minimal beta-cell activity.
plot_quick_view("GCG",
bw_paths = bw_paths,
colors = bw_colors,
mirror = TRUE,
genome = "hg38")The alpha-cell tracks now dominate the gene body; beta-cell tracks are close to baseline, exactly mirroring the INS pattern above.
MAFB
MAFB is the alpha-cell identity factor (mouse Mafa-like role).
plot_quick_view("MAFB",
bw_paths = bw_paths,
colors = bw_colors,
mirror = TRUE,
genome = "hg38")The alpha-cell ATAC and RNA tracks show the expected alpha-biased signal; MAFB is expressed in both cell types developmentally but remains more abundant in adult alpha cells.
Both Alpha and Beta Cell Genes
UCN3
UCN3 (urocortin-3) is expressed in both alpha and beta cells, though most strongly in mature beta cells. We therefore expect some signal in both tracks, with a beta-cell bias.
plot_quick_view("UCN3",
bw_paths = bw_paths,
colors = bw_colors,
mirror = TRUE,
genome = "hg38")Both cell types show signal over the gene body, with the beta-cell side typically stronger — a useful contrast to the strictly cell-type-specific INS and GCG loci above.
Identifying putative enhancers
There is a lot of flexibility for data exploration here. In this example, we search for putative enhancers using two histone marks known to co-occur at enhancer regions: H3K4me1 and H3K27ac.
epiRomicsS4 objects expose their five slots through
public getter functions: annotations(),
meta(), txdb(), organism(), and
genome(). Users should prefer these accessors over reaching
into slots directly with @ or slot() — they
are the stable public API for reading and (via <-
assignment) writing slot contents.
putative_enhancers <-
find_enhancers_by_comarks(
database,
histone_mark_1 = "h3k4me1",
histone_mark_2 = "h3k27ac")
## First 5 putative enhancers
head(as.data.frame(
annotations(putative_enhancers)), 5)Cross-referencing enhancer calls to databases
FANTOM Enhancer Database
Cross-referencing putative enhancers against the FANTOM5 enhancer database increases confidence in the calls.
putative_enhancers_filtered_fantom <-
filter_enhancers(
putative_enhancers,
database,
type = "hg38_custom_fantom")
head(as.data.frame(
annotations(putative_enhancers_filtered_fantom)), 5)Human Islet Regulome Active Enhancers
putative_enhancers_filtered_regulome_active <-
filter_enhancers(
putative_enhancers,
database,
type = "hg38_custom_regulome_active")
head(as.data.frame(
annotations(putative_enhancers_filtered_regulome_active)), 5)Human Islet Regulome Super Enhancers
putative_enhancers_filtered_regulome_super <-
filter_enhancers(
putative_enhancers,
database,
type = "hg38_custom_regulome_super")
head(as.data.frame(
annotations(putative_enhancers_filtered_regulome_super)), 5)Ultra-Conserved Non-Coding Elements (UCNEs)
putative_enhancers_filtered_ucnes <-
filter_enhancers(
putative_enhancers,
database,
type = "hg38_custom_ucnes")
head(as.data.frame(
annotations(putative_enhancers_filtered_ucnes)), 5)Screening for high TF co-binding sites (enhanceosomes)
Enhancers can be redundant and not all play an active regulatory role. By intersecting putative enhancer calls with all available ChIP-seq data, we identify enhanceosome regions with high transcription factor co-binding. Results are sorted by co-binding count.
putative_enhanceosome_fantom <-
find_enhanceosomes(
putative_enhancers_filtered_fantom,
database)
head(as.data.frame(
annotations(putative_enhanceosome_fantom)), 5)
## Evaluate calls on chromosome 1
fantom_df <- as.data.frame(
annotations(putative_enhanceosome_fantom))
head(fantom_df[fantom_df$seqnames == "chr1", ], 5)
## Total enhanceosome regions
length(annotations(putative_enhanceosome_fantom))Transcription factor co-binding analysis
analyze_tf_cobinding() uses Fisher’s exact test to
identify statistically significant TF co-binding pairs within
enhanceosome regions, reporting odds ratios and pointwise mutual
information (PMI).
tf_cobinding <- analyze_tf_cobinding(
putative_enhanceosome_fantom,
database)
tf_cobinding$pairwiseVisualizing enhanceosome regions with signal tracks
plot_tracks(
putative_enhanceosome_fantom,
index = which(names(
annotations(putative_enhanceosome_fantom)) == "14"),
database = database,
track_connection = track_connection)Regulome active-filtered enhanceosomes
putative_enhanceosome_regulome_active <-
find_enhanceosomes(
putative_enhancers_filtered_regulome_active,
database)
head(as.data.frame(
annotations(putative_enhanceosome_regulome_active)), 5)
reg_active_df <- as.data.frame(
annotations(putative_enhanceosome_regulome_active))
head(reg_active_df[reg_active_df$seqnames == "chr1", ], 5)
length(
annotations(putative_enhanceosome_regulome_active))The head of the data frame shows the coordinates and per-TF hit
counts of the first five enhanceosomes that overlap active Islet
Regulome elements, and the final length() prints the total
number of such regions. We pick one for track visualization below.
plot_tracks(
putative_enhanceosome_regulome_active,
index = which(
names(annotations(
putative_enhanceosome_regulome_active)) == "21"),
database = database,
track_connection = track_connection)Regulome super-enhancer filtered enhanceosomes
putative_enhanceosome_regulome_super <-
find_enhanceosomes(
putative_enhancers_filtered_regulome_super,
database)
head(as.data.frame(
annotations(putative_enhanceosome_regulome_super)), 5)
reg_super_df <- as.data.frame(
annotations(putative_enhanceosome_regulome_super))
head(reg_super_df[reg_super_df$seqnames == "chr1", ], 5)
length(
annotations(putative_enhanceosome_regulome_super))As before, the printed data frame and total count summarise the enhanceosome set narrowed to super-enhancer calls; the track plot below renders a representative locus.
plot_tracks(
putative_enhanceosome_regulome_super,
index = which(
names(annotations(
putative_enhanceosome_regulome_super)) == "13"),
database = database,
track_connection = track_connection)UCNE-filtered enhanceosomes
putative_enhanceosome_ucnes <-
find_enhanceosomes(
putative_enhancers_filtered_ucnes,
database)
head(as.data.frame(
annotations(putative_enhanceosome_ucnes)), 5)
## Evaluate calls on chromosome 1
ucne_df <- as.data.frame(
annotations(putative_enhanceosome_ucnes))
head(ucne_df[ucne_df$seqnames == "chr1", ], 5)
## Total enhanceosome regions
length(annotations(putative_enhanceosome_ucnes))The UCNE-filtered result isolates enhanceosomes that sit inside ultra-conserved non-coding elements — a small, high-confidence subset. The track plot below visualises one of these regions.
plot_tracks(
putative_enhanceosome_ucnes,
index = 4,
database = database,
track_connection = track_connection)Integrating differential chromatin accessibility
How can we use putative enhanceosome regions to infer biology between cell states? Here we integrate DiffBind ATAC-seq data showing differential chromatin accessibility between alpha and beta cells.
b.v.a <- read.csv(
file.path(extdata_dir,
"DBA_Beta_Versus_Alpha.csv"))
b.v.a <- GRanges(b.v.a)
## Filter for beta-enriched chromatin regions
beta.enriched <- b.v.a[b.v.a$Fold >= 1 & b.v.a$Conc>=3, ]
## Connect to putative enhanceosomes
beta_enhancer_regions <-
get_regions_of_interest(
putative_enhanceosome_fantom,
beta.enriched)
head(as.data.frame(
annotations(beta_enhancer_regions)), 5)These rows are enhanceosome regions that also overlap beta-enriched differential chromatin accessibility calls, combining the TF co-binding signal with cell-type-specific accessibility evidence. One such region is rendered below.
plot_tracks(
beta_enhancer_regions,
index = which(names(
annotations(beta_enhancer_regions)) == "305"),
database = database,
track_connection = track_connection)Gene-centered visualization
plot_gene_tracks() provides a gene-focused view
combining all database tracks (histone marks, ChIP, annotations) with
BigWig signal—no enhanceosome pipeline required.
plot_gene_tracks("INS",
database,
track_connection)The INS view centres on the beta-cell insulin locus so the ATAC and RNA tracks peak sharply relative to histone-mark context.
plot_gene_tracks("GCG",
database,
track_connection)The GCG view shows the alpha-cell glucagon locus for comparison; the signal pattern should mirror-image INS, with alpha-specific accessibility and transcription dominating.
## Wider padding captures distal regulatory elements
plot_gene_tracks("PDX1",
database,
track_connection,
padding = 5000)The PDX1 view uses wider padding
(padding = 5000) to expose distal regulatory elements
around this master transcription factor.
Chromatin state classification
classify_chromatin_states() classifies genomic regions
by their histone mark combinations into six biologically meaningful
states: active, bivalent, poised, primed, repressed, and unmarked.
Classification follows ChromHMM/Roadmap Epigenomics conventions (Ernst
& Kellis 2012; Kundaje et al. 2015), with TSS proximity
refinement.
chromatin_states <-
classify_chromatin_states(database)
table(chromatin_states$chromatin_state)
table(chromatin_states$genomic_context)
pdx1_states <- chromatin_states[
chromatin_states$seqnames == "chr13" &
chromatin_states$start >= 27900000 &
chromatin_states$end <= 28000000, ]
head(pdx1_states, 10)Visualize the PDX1 locus with chromatin state tracks overlaid:
plot_gene_tracks("PDX1",
database,
track_connection,
chromatin_states = chromatin_states,
show_chromatin = TRUE,
show_enhancer_highlight = TRUE,
padding = 5000)Session Information
sessionInfo()
#> R version 4.5.3 (2026-03-11)
#> Platform: x86_64-pc-linux-gnu
#> Running under: Ubuntu 24.04.4 LTS
#>
#> Matrix products: default
#> BLAS: /usr/lib/x86_64-linux-gnu/openblas-pthread/libblas.so.3
#> LAPACK: /usr/lib/x86_64-linux-gnu/openblas-pthread/libopenblasp-r0.3.26.so; LAPACK version 3.12.0
#>
#> locale:
#> [1] LC_CTYPE=C.UTF-8 LC_NUMERIC=C
#> [3] LC_TIME=C.UTF-8 LC_COLLATE=C.UTF-8
#> [5] LC_MONETARY=C.UTF-8 LC_MESSAGES=C.UTF-8
#> [7] LC_PAPER=C.UTF-8 LC_NAME=C
#> [9] LC_ADDRESS=C LC_TELEPHONE=C
#> [11] LC_MEASUREMENT=C.UTF-8 LC_IDENTIFICATION=C
#>
#> time zone: UTC
#> tzcode source: system (glibc)
#>
#> attached base packages:
#> [1] stats4 stats graphics grDevices utils datasets
#> [7] methods base
#>
#> other attached packages:
#> [1] org.Hs.eg.db_3.22.0
#> [2] TxDb.Hsapiens.UCSC.hg38.knownGene_3.22.0
#> [3] GenomicFeatures_1.62.0
#> [4] AnnotationDbi_1.72.0
#> [5] Biobase_2.70.0
#> [6] GenomicRanges_1.62.1
#> [7] Seqinfo_1.0.0
#> [8] IRanges_2.44.0
#> [9] S4Vectors_0.48.1
#> [10] BiocGenerics_0.56.0
#> [11] generics_0.1.4
#> [12] epiRomics_0.99.5
#>
#> loaded via a namespace (and not attached):
#> [1] splines_4.5.3
#> [2] BiocIO_1.20.0
#> [3] bitops_1.0-9
#> [4] ggplotify_0.1.3
#> [5] filelock_1.0.3
#> [6] tibble_3.3.1
#> [7] R.oo_1.27.1
#> [8] polyclip_1.10-7
#> [9] XML_3.99-0.23
#> [10] lifecycle_1.0.5
#> [11] httr2_1.2.2
#> [12] lattice_0.22-9
#> [13] MASS_7.3-65
#> [14] magrittr_2.0.5
#> [15] sass_0.4.10
#> [16] rmarkdown_2.31
#> [17] jquerylib_0.1.4
#> [18] yaml_2.3.12
#> [19] plotrix_3.8-14
#> [20] ggtangle_0.1.1
#> [21] cowplot_1.2.0
#> [22] DBI_1.3.0
#> [23] RColorBrewer_1.1-3
#> [24] abind_1.4-8
#> [25] purrr_1.2.2
#> [26] R.utils_2.13.0
#> [27] RCurl_1.98-1.18
#> [28] yulab.utils_0.2.4
#> [29] tweenr_2.0.3
#> [30] rappdirs_0.3.4
#> [31] gdtools_0.5.0
#> [32] enrichplot_1.30.5
#> [33] ggrepel_0.9.8
#> [34] tidytree_0.4.7
#> [35] pkgdown_2.2.0
#> [36] ChIPseeker_1.46.1
#> [37] codetools_0.2-20
#> [38] DelayedArray_0.36.1
#> [39] DOSE_4.4.0
#> [40] ggforce_0.5.0
#> [41] tidyselect_1.2.1
#> [42] aplot_0.2.9
#> [43] UCSC.utils_1.6.1
#> [44] farver_2.1.2
#> [45] matrixStats_1.5.0
#> [46] BiocFileCache_3.0.0
#> [47] GenomicAlignments_1.46.0
#> [48] jsonlite_2.0.0
#> [49] annotatr_1.36.0
#> [50] systemfonts_1.3.2
#> [51] tools_4.5.3
#> [52] ggnewscale_0.5.2
#> [53] treeio_1.34.0
#> [54] TxDb.Hsapiens.UCSC.hg19.knownGene_3.22.1
#> [55] ragg_1.5.2
#> [56] Rcpp_1.1.1-1
#> [57] glue_1.8.1
#> [58] SparseArray_1.10.10
#> [59] xfun_0.57
#> [60] qvalue_2.42.0
#> [61] MatrixGenerics_1.22.0
#> [62] GenomeInfoDb_1.46.2
#> [63] dplyr_1.2.1
#> [64] withr_3.0.2
#> [65] BiocManager_1.30.27
#> [66] fastmap_1.2.0
#> [67] boot_1.3-32
#> [68] caTools_1.18.3
#> [69] digest_0.6.39
#> [70] R6_2.6.1
#> [71] gridGraphics_0.5-1
#> [72] textshaping_1.0.5
#> [73] GO.db_3.22.0
#> [74] gtools_3.9.5
#> [75] RSQLite_2.4.6
#> [76] cigarillo_1.0.0
#> [77] R.methodsS3_1.8.2
#> [78] tidyr_1.3.2
#> [79] fontLiberation_0.1.0
#> [80] data.table_1.18.2.1
#> [81] rtracklayer_1.70.1
#> [82] httr_1.4.8
#> [83] htmlwidgets_1.6.4
#> [84] S4Arrays_1.10.1
#> [85] scatterpie_0.2.6
#> [86] regioneR_1.42.0
#> [87] pkgconfig_2.0.3
#> [88] gtable_0.3.6
#> [89] blob_1.3.0
#> [90] S7_0.2.1-1
#> [91] XVector_0.50.0
#> [92] htmltools_0.5.9
#> [93] fontBitstreamVera_0.1.1
#> [94] fgsea_1.36.2
#> [95] scales_1.4.0
#> [96] png_0.1-9
#> [97] ggfun_0.2.0
#> [98] knitr_1.51
#> [99] tzdb_0.5.0
#> [100] reshape2_1.4.5
#> [101] rjson_0.2.23
#> [102] nlme_3.1-168
#> [103] curl_7.0.0
#> [104] cachem_1.1.0
#> [105] stringr_1.6.0
#> [106] BiocVersion_3.22.0
#> [107] KernSmooth_2.23-26
#> [108] parallel_4.5.3
#> [109] restfulr_0.0.16
#> [110] desc_1.4.3
#> [111] pillar_1.11.1
#> [112] grid_4.5.3
#> [113] vctrs_0.7.3
#> [114] gplots_3.3.0
#> [115] tidydr_0.0.6
#> [116] dbplyr_2.5.2
#> [117] cluster_2.1.8.2
#> [118] evaluate_1.0.5
#> [119] readr_2.2.0
#> [120] cli_3.6.6
#> [121] compiler_4.5.3
#> [122] Rsamtools_2.26.0
#> [123] rlang_1.2.0
#> [124] crayon_1.5.3
#> [125] plyr_1.8.9
#> [126] fs_2.1.0
#> [127] ggiraph_0.9.6
#> [128] stringi_1.8.7
#> [129] BiocParallel_1.44.0
#> [130] Biostrings_2.78.0
#> [131] lazyeval_0.2.3
#> [132] GOSemSim_2.36.0
#> [133] fontquiver_0.2.1
#> [134] Matrix_1.7-4
#> [135] BSgenome_1.78.0
#> [136] hms_1.1.4
#> [137] patchwork_1.3.2
#> [138] bit64_4.6.0-1
#> [139] ggplot2_4.0.2
#> [140] KEGGREST_1.50.0
#> [141] SummarizedExperiment_1.40.0
#> [142] AnnotationHub_4.0.0
#> [143] igraph_2.2.3
#> [144] memoise_2.0.1
#> [145] bslib_0.10.0
#> [146] ggtree_4.0.5
#> [147] fastmatch_1.1-8
#> [148] bit_4.6.0
#> [149] ape_5.8-1