SEARchways vignette
Gene set enrichment analyses
Kim Dill-McFarland, kadm@uw.edu
version May 04, 2025
Overview
This package performs gene set enrichment analyses including functional/pathway enrichment and fold change GSEA. Checkout BIGpicture for ways to plot these results.
Setup
Software and packages
This pipeline should be completed in R version 4.5.0 (2025-04-11) or newer. Please also download the following packages.
# CRAN packages
install.packages(c("tidyverse"))
# Bioconductor packages
## These are dependencies for kimma that are beneficial to download separately
install.packages("BiocManager")
BiocManager::install(c("limma", "edgeR"))
# GitHub packages
install.packages("devtools")
devtools::install_github("BIGslu/kimma")
devtools::install_github("BIGslu/SEARchways")And load them into your current R session.
library(tidyverse)## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.1.4 ✔ readr 2.1.5
## ✔ forcats 1.0.0 ✔ stringr 1.5.1
## ✔ ggplot2 3.5.2 ✔ tibble 3.2.1
## ✔ lubridate 1.9.4 ✔ tidyr 1.3.1
## ✔ purrr 1.0.4
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errorslibrary(kimma)
library(SEARchways)To ensure reproducibility with this document, set the random seed value.
set.seed(651)Data background
Example data were obtained from media controls and human rhinovirus (HRV) infected human plasmacytoid dendritic cells. You can learn more about these data in
Dill-McFarland KA, Schwartz JT, Zhao H, Shao B, Fulkerson PC, Altman MC, Gill MA. 2022. Eosinophil-mediated suppression and Anti-IL-5 enhancement of plasmacytoid dendritic cell interferon responses in asthma. J Allergy Clin Immunol. Epub ahead of print. doi: 10.1016/j.jaci.2022.03.025. — GitHub
Specifically, this vignette uses bulk human RNA-seq processed from fastq to counts using SEAsnake and counts to voom using edgeR/limma. This results in voom-normalized, log2 counts per million (CPM) gene expression and associated sample and gene metadata. In total, 6 donors and 1000 random genes were selected for this vignette. Actual donor identifiers have been altered for privacy.
Load data
Gene fold change estimates.
example.gene.list <- SEARchways::example.gene.list
str(example.gene.list)## List of 2
## $ HRV1: Named num [1:100] 0.1636 0.26787 0.53233 -0.09021 0.00173 ...
## ..- attr(*, "names")= chr [1:100] "ENSG00000198898" "ENSG00000121410" "ENSG00000092096" "ENSG00000076928" ...
## $ HRV2: Named num [1:100] -4.663 -0.367 -0.876 0.898 -0.601 ...
## ..- attr(*, "names")= chr [1:100] "ENSG00000174123" "ENSG00000165675" "ENSG00000103024" "ENSG00000151067" ...Extract gene names only.
example.gene.list2 <- list(
"HRV1"=names(example.gene.list$HRV1),
"HRV2"=names(SEARchways::example.gene.list$HRV2))
str(example.gene.list2)## List of 2
## $ HRV1: chr [1:100] "ENSG00000198898" "ENSG00000121410" "ENSG00000092096" "ENSG00000076928" ...
## $ HRV2: chr [1:100] "ENSG00000174123" "ENSG00000165675" "ENSG00000103024" "ENSG00000151067" ...Create data frame formats as well.
example.gene.df <- data.frame(
group = "HRV1",
gene = names(example.gene.list$HRV1),
logFC = example.gene.list$HRV1
)
head(example.gene.df)## group gene logFC
## ENSG00000198898 HRV1 ENSG00000198898 0.16360114
## ENSG00000121410 HRV1 ENSG00000121410 0.26786954
## ENSG00000092096 HRV1 ENSG00000092096 0.53232684
## ENSG00000076928 HRV1 ENSG00000076928 -0.09020896
## ENSG00000171307 HRV1 ENSG00000171307 0.00172912
## ENSG00000116863 HRV1 ENSG00000116863 0.90057251example.gene.df2 <- example.gene.df[,1:2]
head(example.gene.df2)## group gene
## ENSG00000198898 HRV1 ENSG00000198898
## ENSG00000121410 HRV1 ENSG00000121410
## ENSG00000092096 HRV1 ENSG00000092096
## ENSG00000076928 HRV1 ENSG00000076928
## ENSG00000171307 HRV1 ENSG00000171307
## ENSG00000116863 HRV1 ENSG00000116863Functional / pathway / hypergeometric enrichment
Fisher’s exact test of selected genes.
flexEnrich: Flexible enrichment
flexEnrich allows more customization than
BIGprofiler but can run slower. By default, enrichment is
performed against a protein-coding background and gene sets must be
between 10 and 1E10 genes.
flexEnrich(
gene_list = example.gene.list2,
ID = "ENSEMBL", species = "human",
collection = "H")## # A tibble: 51 × 11
## group n_query_genes n_background_genes gs_collection pathway n_pathway_genes
## <chr> <int> <int> <chr> <chr> <int>
## 1 HRV1 100 21940 H HALLMAR… 200
## 2 HRV1 100 21940 H HALLMAR… 139
## 3 HRV1 100 21940 H HALLMAR… 201
## 4 HRV1 100 21940 H HALLMAR… 200
## 5 HRV1 100 21940 H HALLMAR… 200
## 6 HRV1 100 21940 H HALLMAR… 200
## 7 HRV1 100 21940 H HALLMAR… 158
## 8 HRV1 100 21940 H HALLMAR… 200
## 9 HRV1 100 21940 H HALLMAR… 200
## 10 HRV1 100 21940 H HALLMAR… 200
## # ℹ 41 more rows
## # ℹ 5 more variables: n_query_genes_in_pathway <int>, `k/K` <dbl>, pval <dbl>,
## # FDR <dbl>, genes <list>Use a data frame input.
flexEnrich(
gene_df = example.gene.df2,
ID = "ENSEMBL", species = "human",
collection = "H")Change the background to a custom one (custom_bg) or
include all genes by turning off protein_coding.
flexEnrich(
gene_list = example.gene.list2,
ID = "ENSEMBL", species = "human",
collection = "H",
protein_coding = FALSE)Restrict the gene sets included to those with at least 3 query genes and total size between 100 and 1000 genes.
flexEnrich(
gene_list = example.gene.list2,
ID = "ENSEMBL", species = "human",
collection = "H",
minOverlap = 3,
minGeneSetSize = 100,
maxGeneSetSize = 1000)Speed up computation by not saving the query genes found in each gene set.
flexEnrich(
gene_list = example.gene.list2,
ID = "ENSEMBL", species = "human",
collection = "H",
print_genes = FALSE)BIGprofiler: Less flexible but faster enrichment
BIGprofiler performs the same enrichment as
flexEnrich but with fewer options. It does not allow you to
alter the background, always utilizing protein-coding only. Thus, when
this background is desired, BIGprofiler may be utilized
since it runs faster than flexEnrich on large queries.
BIGprofiler(
gene_list = example.gene.list2,
ID = "ENSEMBL", species = "human",
collection = "H")## [1] "HRV1"## ## [1] "HRV2"## # A tibble: 51 × 13
## group n_query_genes gs_collection gs_subcollection n_background_genes pathway
## <chr> <int> <chr> <chr> <dbl> <chr>
## 1 HRV1 100 H "" 4386 HALLMA…
## 2 HRV1 100 H "" 4386 HALLMA…
## 3 HRV1 100 H "" 4386 HALLMA…
## 4 HRV1 100 H "" 4386 HALLMA…
## 5 HRV1 100 H "" 4386 HALLMA…
## 6 HRV1 100 H "" 4386 HALLMA…
## 7 HRV1 100 H "" 4386 HALLMA…
## 8 HRV1 100 H "" 4386 HALLMA…
## 9 HRV1 100 H "" 4386 HALLMA…
## 10 HRV1 100 H "" 4386 HALLMA…
## # ℹ 41 more rows
## # ℹ 7 more variables: n_pathway_genes <dbl>, n_query_genes_in_pathway <dbl>,
## # `k/K` <dbl>, pval <dbl>, FDR <dbl>, qvalue <dbl>, genes <list>Restrict the gene sets included to those with total size between 100
and 1000 genes (similar to flexEnrich).
BIGprofiler(
gene_list = example.gene.list2,
ID = "ENSEMBL", species = "human",
collection = "H",
minGSSize = 100, maxGSSize = 1000)iterEnrich: Iterative enrichment
iterEnrich allows you to run enrichment on inputs other
than genes. For example, if you have SNPs, each SNP may annotate to >
1 gene. Regular enrichment would result in over-representation of
pathways associated these multiple annotation SNPs. Thus,
iterEnrich subsamples gene lists to 1 annotation per
feature (e.g. SNP, probe, etc) and runs enrichment across many
subsamples to determine a more representative result.
Here, we create a fake gene annotation where each SNP is associated with 2 genes.
example.gene.anno <- data.frame(
SNP = c(paste0("SNP",1:50), paste0("SNP",1:50)),
gene = sample(example.gene.list2[["HRV1"]],
replace = TRUE, size = 100)
)
head(example.gene.anno)## SNP gene
## 1 SNP1 ENSG00000076928
## 2 SNP2 ENSG00000213995
## 3 SNP3 ENSG00000073803
## 4 SNP4 ENSG00000184232
## 5 SNP5 ENSG00000107018
## 6 SNP6 ENSG00000105289Run iterative enrichment. We reduce the number of iters and increase the processors to run quickly for this vignette.
iter <- iterEnrich(
anno_df = example.gene.anno,
ID = "ENSEMBL", species = "human",
collection = "H",
niter = 10, ncores = 6)
# Summary results across iters
iter$summary## group pathway pval_min pval_max
## 1 anno_df HALLMARK_COAGULATION 0.02908510 1.00000000
## 2 anno_df HALLMARK_COMPLEMENT 0.32085006 1.00000000
## 3 anno_df HALLMARK_E2F_TARGETS 0.30054682 0.31953545
## 4 anno_df HALLMARK_ESTROGEN_RESPONSE_EARLY 0.30054682 0.31953545
## 5 anno_df HALLMARK_ESTROGEN_RESPONSE_LATE 0.30054682 0.31953545
## 6 anno_df HALLMARK_G2M_CHECKPOINT 0.04915364 0.05611768
## 7 anno_df HALLMARK_GLYCOLYSIS 0.30054682 1.00000000
## 8 anno_df HALLMARK_HEME_METABOLISM 0.31953545 1.00000000
## 9 anno_df HALLMARK_INFLAMMATORY_RESPONSE 0.30180159 0.32085006
## 10 anno_df HALLMARK_KRAS_SIGNALING_UP 0.31953545 1.00000000
## 11 anno_df HALLMARK_OXIDATIVE_PHOSPHORYLATION 0.04915364 0.31953545
## 12 anno_df HALLMARK_P53_PATHWAY 0.05611768 0.30054682
## 13 anno_df HALLMARK_PI3K_AKT_MTOR_SIGNALING 0.17076862 0.18263882
## 14 anno_df HALLMARK_REACTIVE_OXYGEN_SPECIES_PATHWAY 0.08357430 0.08971189
## 15 anno_df HALLMARK_SPERMATOGENESIS 0.21409813 1.00000000
## 16 anno_df HALLMARK_UV_RESPONSE_DN 0.22666112 0.24181594
## 17 anno_df HALLMARK_UV_RESPONSE_UP 0.24581556 0.26202152
## pval_median FDR k_median K k/K_median genes
## 1 0.02908510 0.3180002 2 139 0.014388489 ENSG0000....
## 2 0.32085006 0.3636301 1 201 0.004975124 ENSG0000....
## 3 0.31953545 0.3636301 1 200 0.005000000 ENSG0000....
## 4 0.31953545 0.3636301 1 200 0.005000000 ENSG0000....
## 5 0.31953545 0.3636301 1 200 0.005000000 ENSG0000....
## 6 0.05611768 0.3180002 2 200 0.010000000 ENSG0000....
## 7 1.00000000 1.0000000 0 200 0.000000000 ENSG0000....
## 8 0.31953545 0.3636301 1 200 0.005000000 ENSG0000....
## 9 0.32085006 0.3636301 1 201 0.004975124 ENSG0000....
## 10 0.31953545 0.3636301 1 200 0.005000000 ENSG0000....
## 11 0.31953545 0.3636301 1 200 0.005000000 ENSG0000....
## 12 0.05611768 0.3180002 2 200 0.010000000 ENSG0000....
## 13 0.18263882 0.3636301 1 105 0.009523810 ENSG0000....
## 14 0.08971189 0.3636301 1 49 0.020408163 ENSG0000....
## 15 1.00000000 1.0000000 0 135 0.000000000 ENSG0000....
## 16 0.24181594 0.3636301 1 144 0.006944444 ENSG0000....
## 17 0.26202152 0.3636301 1 158 0.006329114 ENSG0000....# All individual iter results
head(iter$`k/K_iterations`)## pathway k/K_1 k/K_2 k/K_3 k/K_4
## 10 HALLMARK_COAGULATION 0.000 0.014388489 0.014388489 0.014388489
## 11 HALLMARK_COMPLEMENT 0.000 0.004975124 0.004975124 0.004975124
## 13 HALLMARK_E2F_TARGETS 0.005 0.005000000 0.005000000 0.005000000
## 15 HALLMARK_ESTROGEN_RESPONSE_EARLY 0.005 0.005000000 0.005000000 0.005000000
## 16 HALLMARK_ESTROGEN_RESPONSE_LATE 0.005 0.005000000 0.005000000 0.005000000
## 18 HALLMARK_G2M_CHECKPOINT 0.010 0.010000000 0.010000000 0.010000000
## k/K_5 k/K_6 k/K_7 k/K_8 k/K_9 k/K_10
## 10 0.014388489 0.014388489 0.014388489 0.000 0.014388489 0.000
## 11 0.004975124 0.004975124 0.004975124 0.000 0.004975124 0.000
## 13 0.005000000 0.005000000 0.005000000 0.005 0.005000000 0.005
## 15 0.005000000 0.005000000 0.005000000 0.005 0.005000000 0.005
## 16 0.005000000 0.005000000 0.005000000 0.005 0.005000000 0.005
## 18 0.010000000 0.010000000 0.010000000 0.010 0.010000000 0.010head(iter$p_iterations)## pathway pval_1 pval_2 pval_3 pval_4
## 10 HALLMARK_COAGULATION 1.00000000 0.02908510 0.02908510 0.02908510
## 11 HALLMARK_COMPLEMENT 1.00000000 0.32085006 0.32085006 0.32085006
## 13 HALLMARK_E2F_TARGETS 0.30054682 0.31953545 0.31953545 0.31953545
## 15 HALLMARK_ESTROGEN_RESPONSE_EARLY 0.30054682 0.31953545 0.31953545 0.31953545
## 16 HALLMARK_ESTROGEN_RESPONSE_LATE 0.30054682 0.31953545 0.31953545 0.31953545
## 18 HALLMARK_G2M_CHECKPOINT 0.04915364 0.05611768 0.05611768 0.05611768
## pval_5 pval_6 pval_7 pval_8 pval_9 pval_10
## 10 0.02908510 0.02908510 0.02908510 1.00000000 0.02908510 1.00000000
## 11 0.32085006 0.32085006 0.32085006 1.00000000 0.32085006 1.00000000
## 13 0.31953545 0.31953545 0.31953545 0.30054682 0.31953545 0.30054682
## 15 0.31953545 0.31953545 0.31953545 0.30054682 0.31953545 0.30054682
## 16 0.31953545 0.31953545 0.31953545 0.30054682 0.31953545 0.30054682
## 18 0.05611768 0.05611768 0.05611768 0.04915364 0.05611768 0.04915364Change the background to a custom one (custom_bg) or
include all genes by turning off protein_coding.
iterEnrich(
anno_df = example.gene.anno,
ID = "ENSEMBL", species = "human",
collection = "H",
niter = 10, ncores = 6,
protein_coding = FALSE)Restrict the gene sets included to those with at least 3 query genes and total size between 100 and 1000 genes.
iterEnrich(
anno_df = example.gene.anno,
ID = "ENSEMBL", species = "human",
collection = "H",
niter = 10, ncores = 6,
minOverlap = 3,
minGeneSetSize = 100,
maxGeneSetSize = 1000)Speed up computation by not saving the query genes found in each gene set.
iterEnrich(
anno_df = example.gene.anno,
ID = "ENSEMBL", species = "human",
collection = "H",
niter = 10, ncores = 6,
print_genes = FALSE)Other parameters
- Change the column names to match your annotation data with
anno_featColandanno_annotationCol - Change the FDR method with
p_adjust. Seep.adjust.methodsfor options
GSEA
Gene set enrichment analysis (GSEA) as described by Broad. This method
utilizes the fgsea framework.
BIGsea(
gene_list = example.gene.list,
ID = "ENSEMBL", species = "human",
collection = "C2", subcollection = "CP")## Running HRV1## Running HRV2## Warning in preparePathwaysAndStats(pathways, stats, minSize, maxSize, gseaParam, : There are ties in the preranked stats (1% of the list).
## The order of those tied genes will be arbitrary, which may produce unexpected results.## group gs_collection gs_subcollection method
## 1 HRV1 C2 CP multi
## 2 HRV2 C2 CP multi
## 3 HRV2 C2 CP multi
## pathway pval FDR
## 1 REACTOME_TRANSPORT_OF_SMALL_MOLECULES 0.8736264 0.8736264
## 2 REACTOME_POST_TRANSLATIONAL_PROTEIN_MODIFICATION 0.5366430 0.5366430
## 3 REACTOME_RNA_POLYMERASE_II_TRANSCRIPTION 0.1465721 0.2931442
## log2err ES NES size leadingEdge
## 1 0.04782996 0.2925544 0.6594805 10 ENSG0000....
## 2 0.08431443 0.4127165 0.9346453 10 ENSG0000....
## 3 0.17821987 0.6103486 1.3822065 10 ENSG0000....Use a data frame input.
BIGsea(
gene_list = example.gene.df,
ID = "ENSEMBL", species = "human",
collection = "C2", subcollection = "CP")Restrict the gene sets included to those with total size between 100 and 1000 genes.
BIGsea(
gene_list = example.gene.list,
ID = "ENSEMBL", species = "human",
collection = "C2", subcollection = "CP",
minGeneSetSize = 100, maxGeneSetSize = 1000)Change the GSEA method. By default, rand = "multi",
which runs an adaptive multilevel splitting Monte Carlo approach.
For very small sample groups, it is recommended to run “simple” randomization, which randomizes genes in gene sets only.
BIGsea(
gene_list = example.gene.list,
ID = "ENSEMBL", species = "human",
collection = "C2", subcollection = "CP",
rand = "simple")Alternatively, choose the “label” method to randomize samples in
groups and calculate new fold changes for each randomization. This
method runs kimma under the hood to determine fold changes.
Thus, you have all the model flexibility of kimma, not just
a simple A vs B fold change.
The “label” method requires dat and
rand_var in order to calculate fold changes from the
original expression data. You also must provide
kimma::kmFit model parameters model,
run_, use_weights, etc.
Finally, we strongly recommend only running on select gene sets with
pw and reducing computational time with multiple
processors. Here, we run only 3 permutations for speed as
well.
#Force original fold change list name to match example.voom$targets variable
example.gene.list.virus <- list(
"virus" = example.gene.list$HRV1
)
gsea_label <- BIGsea(
gene_list = example.gene.list.virus,
ID = "ENSEMBL", species = "human",
collection = "C2", subcollection = "CP",
rand = "label",
#Data for fold change calculations
dat = kimma::example.voom,
rand_var = "virus",
#Parameters for kimma models
model="~virus",
run_lm=TRUE, use_weights=TRUE,
#Recommend only running on select gene sets of interest
pw = "REACTOME_TRANSPORT_OF_SMALL_MOLECULES",
nperm = 3, processors = 3)## Running virus## Loading required namespace: limma## Label randomization runs multiple kimma models. This may take a long time.## [1] "kimma permutation 1"
## [1] "kimma permutation 2"
## [1] "kimma permutation 3"
## | | | 0% | |======================= | 33% | |=============================================== | 67% | |======================================================================| 100%#GSEA results
gsea_label$gsea## group gs_collection gs_subcollection method
## 1 virus C2 CP label
## pathway pval FDR ES NES nMoreExtreme
## 1 REACTOME_TRANSPORT_OF_SMALL_MOLECULES 1 1 0.2925544 NaN 0
## size leadingEdge
## 1 10 ENSG0000....#Fold change estimates for each permutation
head(gsea_label$estimates)## V1 V2 V3
## ENSG00000000460 -0.48411710 -0.3572232 -0.22709577
## ENSG00000001460 0.69172192 -0.4946608 -0.06399949
## ENSG00000002587 -0.02961869 0.1458342 0.37694009
## ENSG00000005189 1.72116622 0.6523883 0.07678363
## ENSG00000005801 0.02720169 0.2159433 -0.19241759
## ENSG00000005884 0.89349003 0.1290022 -1.49447436Once you’ve run the “label” method once, you can reuse the fold
change permutations for additional gene sets with
rand_est.
gsea_label2 <- BIGsea(
gene_list = example.gene.list.virus,
ID = "ENSEMBL", species = "human",
collection = "C2", subcollection = "CP",
rand = "label",
#Data for fold change calculations
dat = kimma::example.voom,
rand_var = "virus",
#Pre-made fold change results
rand_est = gsea_label$estimates,
#Recommend only running on select gene sets of interest
pw = "REACTOME_POST_TRANSLATIONAL_PROTEIN_MODIFICATION",
nperm = 3, processors = 3)## Running virusR session
sessionInfo()## R version 4.5.0 (2025-04-11)
## Platform: aarch64-apple-darwin20
## Running under: macOS Sequoia 15.4.1
##
## Matrix products: default
## BLAS: /Library/Frameworks/R.framework/Versions/4.5-arm64/Resources/lib/libRblas.0.dylib
## LAPACK: /Library/Frameworks/R.framework/Versions/4.5-arm64/Resources/lib/libRlapack.dylib; LAPACK version 3.12.1
##
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
##
## time zone: America/Los_Angeles
## tzcode source: internal
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] SEARchways_1.5.0 kimma_1.6.3 lubridate_1.9.4 forcats_1.0.0
## [5] stringr_1.5.1 dplyr_1.1.4 purrr_1.0.4 readr_2.1.5
## [9] tidyr_1.3.1 tibble_3.2.1 ggplot2_3.5.2 tidyverse_2.0.0
##
## loaded via a namespace (and not attached):
## [1] DBI_1.2.3 gson_0.1.0 rlang_1.1.6
## [4] magrittr_2.0.3 DOSE_4.2.0 msigdbr_10.0.2
## [7] compiler_4.5.0 RSQLite_2.3.9 png_0.1-8
## [10] vctrs_0.6.5 reshape2_1.4.4 pkgconfig_2.0.3
## [13] crayon_1.5.3 fastmap_1.2.0 XVector_0.48.0
## [16] utf8_1.2.4 rmarkdown_2.29 tzdb_0.5.0
## [19] enrichplot_1.28.1 UCSC.utils_1.4.0 bit_4.6.0
## [22] xfun_0.52 cachem_1.1.0 aplot_0.2.5
## [25] GenomeInfoDb_1.44.0 jsonlite_2.0.0 blob_1.2.4
## [28] BiocParallel_1.42.0 parallel_4.5.0 R6_2.6.1
## [31] bslib_0.9.0 stringi_1.8.7 RColorBrewer_1.1-3
## [34] limma_3.64.0 msigdbdf_24.1.1 jquerylib_0.1.4
## [37] GOSemSim_2.34.0 Rcpp_1.0.14 assertthat_0.2.1
## [40] iterators_1.0.14 knitr_1.50 ggtangle_0.0.6
## [43] R.utils_2.13.0 IRanges_2.42.0 igraph_2.1.4
## [46] Matrix_1.7-3 splines_4.5.0 timechange_0.3.0
## [49] tidyselect_1.2.1 qvalue_2.40.0 rstudioapi_0.17.1
## [52] yaml_2.3.10 doParallel_1.0.17 codetools_0.2-20
## [55] lattice_0.22-7 plyr_1.8.9 treeio_1.32.0
## [58] Biobase_2.68.0 withr_3.0.2 KEGGREST_1.48.0
## [61] evaluate_1.0.3 gridGraphics_0.5-1 Biostrings_2.76.0
## [64] ggtree_3.16.0 pillar_1.10.2 foreach_1.5.2
## [67] stats4_4.5.0 clusterProfiler_4.16.0 ggfun_0.1.8
## [70] generics_0.1.3 S4Vectors_0.46.0 hms_1.1.3
## [73] tidytree_0.4.6 scales_1.4.0 glue_1.8.0
## [76] lazyeval_0.2.2 tools_4.5.0 data.table_1.17.0
## [79] fgsea_1.34.0 babelgene_22.9 fs_1.6.6
## [82] fastmatch_1.1-6 cowplot_1.1.3 grid_4.5.0
## [85] ape_5.8-1 AnnotationDbi_1.70.0 colorspace_2.1-1
## [88] nlme_3.1-168 patchwork_1.3.0 GenomeInfoDbData_1.2.14
## [91] cli_3.6.5 gtable_0.3.6 R.methodsS3_1.8.2
## [94] yulab.utils_0.2.0 sass_0.4.10 digest_0.6.37
## [97] BiocGenerics_0.54.0 ggrepel_0.9.6 ggplotify_0.1.2
## [100] farver_2.1.2 memoise_2.0.1 htmltools_0.5.8.1
## [103] R.oo_1.27.0 lifecycle_1.0.4 httr_1.4.7
## [106] statmod_1.5.0 GO.db_3.21.0 bit64_4.6.0-1