Table of Contents
- Introduction
- Installation
- Common Workflow
- Further Reading
- Contact
- Session Information
1. Introduction
markeR is an R package designed to provide a modular, flexible framework for evaluating gene sets as phenotypic markers using transcriptomic data. It supports researchers in both quantitative analysis and visual exploration of gene set behaviour across experimental or clinical phenotypes.
In this vignette, we will illustrate markeR functionalities using a pre-processed gene expression dataset and metadata derived from the Marthandan et al. (2016) study (GSE63577). This dataset contains human fibroblast samples cultured under two conditions: replicative senescence and proliferative control.
1.1. Primary Objectives
The primary objectives of markeR are to:
- Facilitate reproducible and standardized quantification of gene set activity.
- Provide flexible options for both score-based and enrichment-based analyses.
- Enable comparison of user-defined signatures against reference gene sets (e.g., MSigDB) for biological interpretation.
- Offer modular visualization and evaluation tools to explore gene set behaviour at both the set and individual gene level.
1.2. Features and Capabilities
The package integrates multiple approaches for characterizing phenotypes:
- Score-based methods: log2-median, ranking, and single-sample GSEA (ssGSEA) to quantify the coordinated expression of genes in a set.
- Enrichment-based methods: GSEA using moderated t- or B-statistics, with options to handle unidirectional and bidirectional gene sets.
- Gene-level exploration: expression heatmaps, violin plots, ROC curves, AUC calculations, effect size measures, and PCA analysis to assess individual gene contributions.
The package is designed to be fully customizable, supporting diverse
visualization strategies via ggplot2
,
ComplexHeatmap
, ggpubr
, cowplot
,
and grid.
Its modular structure allows easy integration of
new functionalities while providing a robust framework for reproducible,
standardized phenotypic characterization of gene sets.
2. Installation
Install the latest development release of markeR from GitHub with:
devtools::install_github("DiseaseTranscriptomicsLab/markeR@*release")
3. Common Workflow
markeR
provides a modular pipeline to quantify
transcriptomic signatures and assess their association with phenotypic
or clinical variables.
3.1. Input Requirements
Gene sets
A named list where each element is a gene set.
- Use a character vector when direction is unknown (unidirectional).
- Use a data frame with columns
gene
anddirection
(values+1
for up,-1
for down) when direction is known (bidirectional).
# Example
gene_sets
## $Set1
## [1] "GeneA" "GeneB" "GeneC" "GeneD"
##
## $Set2
## gene direction
## 1 GeneX 1
## 2 GeneY -1
## 3 GeneZ 1
For this vignette, we will use three senescence-related gene sets
that ship with the package: LiteratureMarkers
(bidirectional), REACTOME_CELLULAR_SENESCENCE
(undirected), and HernandezSegura (bidirectional). See
?genesets_example
.
# Load example gene sets
data(genesets_example)
Expression data frame
A filtered and normalised gene expression matrix (genes × samples). Row names must be gene identifiers; column names must match sample IDs in the metadata.
In this vignette, we use a pre‑processed dataset from Marthandan et
al. (2016, GSE63577) with human fibroblasts under replicative
senescence and proliferative control. See
?counts_example
for structure.
# Load example expression data
data(counts_example)
counts_example[1:5, 1:5]
## SRR1660534 SRR1660535 SRR1660536 SRR1660537 SRR1660538
## A1BG 9.94566 9.476768 8.229231 8.515083 7.806479
## A1BG-AS1 12.08655 11.550303 12.283976 7.580694 7.312666
## A2M 77.50289 56.612839 58.860268 8.997624 6.981857
## A4GALT 14.74183 15.226083 14.815891 14.675780 15.222488
## AAAS 47.92755 46.292377 43.965972 47.109493 47.213739
Sample metadata
A data frame with annotations per sample. Row names must match the expression matrix column names; the first column contains sample IDs.
We use the accompanying metadata from the Marthandan et al. (2016)
study. See ?metadata_example
.
## sampleID DatasetID CellType Condition SenescentType
## 252 SRR1660534 Marthandan2016 Fibroblast Senescent Telomere shortening
## 253 SRR1660535 Marthandan2016 Fibroblast Senescent Telomere shortening
## 254 SRR1660536 Marthandan2016 Fibroblast Senescent Telomere shortening
## 255 SRR1660537 Marthandan2016 Fibroblast Proliferative none
## 256 SRR1660538 Marthandan2016 Fibroblast Proliferative none
## 257 SRR1660539 Marthandan2016 Fibroblast Proliferative none
## Treatment
## 252 PD72 (Replicative senescence)
## 253 PD72 (Replicative senescence)
## 254 PD72 (Replicative senescence)
## 255 young
## 256 young
## 257 young
3.2. Select Mode of Analysis
Benchmarking Mode: Evaluate one or more gene sets against a metadata variable using a standardised scoring and effect size framework. This mode provides comprehensive visualisations and comparisons across methods.
Discovery Mode: Explore how a single, well-characterised gene set relates to a specific variable of interest. Suitable for hypothesis generation or signature projection.
3.3. Choose a Quantification Approach
markeR
supports two complementary strategies for
quantifying the association between gene sets and phenotypes:
3.3.1 Score-Based Approach
This strategy generates a single numeric score per sample, reflecting the activity of a gene set. It enables flexible downstream analyses, including comparisons across phenotypic groups.
Three scoring methods are available:
Log2-median: Calculates the median log2 expression of the genes in the set. Sensitive to absolute shifts in expression.
Ranking: Ranks all genes within each sample and averages the ranks of gene set members. Captures relative ordering rather than magnitude.
ssGSEA: Computes a single-sample gene set enrichment score using the ssGSEA algorithm. Reflects the coordinated up- or down-regulation of the set in each sample.
Robust gene sets are expected to perform consistently across all three.
3.3.2 Enrichment-Based Approach
This approach uses a classical gene set enrichment analysis (GSEA) framework to evaluate whether the gene set is significantly overrepresented at the top or bottom of a ranked list of genes (e.g., ranked by fold change or correlation with phenotype).
- GSEA: Computes a Normalised Enrichment Score (NES) for each contrast or variable of interest, adjusting for gene set size and multiple testing.
Use this approach when interested in collective behaviour of gene sets in relation to ranked differential signals.
3.4. Visualisation and Evaluation
3.4.1 Benchmarking Mode
In Benchmarking Mode, markeR
offers a
range of visual summaries:
- Violin or scatter plots showing score distributions by phenotype
- Volcano plots and heatmaps based on effect sizes (Cohen’s d or f)
- ROC curves and AUC values
- Null distribution testing using random gene sets matched for size and directionality
- Lollipop plots summarising enrichment scores (NES) with adjusted p-values
- Enrichment plots showing running enrichment scores across ranked gene lists
- Volcano plots
Example of common workflows in Benchmarking mode (full tutorial here):
Score-based approach
We begin by quantifying gene set activity using score-based methods. These approaches generate numeric scores per sample that reflect the coordinated expression of genes in each set. This allows easy comparison of gene set behavior across phenotypic groups, while giving a single score per sample, which can be useful in certain contexts.
First, we compute and plot scores using the log2-median method, which calculates median-centered expression for each gene and averages across the gene set to provide a robust summary of coordinated activity. We then assess whether the available gene sets distinguish between Senescent and Proliferative conditions. These results suggest that the HernandezSegura and LiteratureMarkers gene sets show very strong effect sizes (|Cohen’s d|), while the REACTOME_CELLULAR_SENESCENCE gene set shows a more modest effect. The distributions of scores also overlap more for the latter.
# Quantification approach: scores
# Method: log2-median
PlotScores(data = counts_example,
metadata = metadata_example,
gene_sets = genesets_example,
Variable="Condition",
method ="logmedian",
nrow = 1,
pointSize=4,
title="Marthandan et al. 2016",
widthTitle = 24,
labsize=12,
titlesize = 12)
Next, we calculate scores using multiple methods (log2-median, ranking, and ssGSEA) to compare results across scoring strategies. The output includes heatmaps and volcano plots to visualize effect sizes (|Cohen’s d|) between conditions. These analyses confirm that the HernandezSegura and LiteratureMarkers gene sets consistently discriminate Senescence from Proliferative samples across all methods, whereas REACTOME_CELLULAR_SENESCENCE does not.
Overall_Scores <- PlotScores(data = counts_example,
metadata = metadata_example,
gene_sets=genesets_example,
Variable="Condition",
method ="all",
ncol = NULL,
nrow = 1,
widthTitle=30,
limits = c(0,3.5),
title="Marthandan et al. 2016",
titlesize = 10,
ColorValues = list(heatmap=c("#F9F4AE", "#B44141"),
volcano=signature_colors <- c(
HernandezSegura = "#A07395",
REACTOME_CELLULAR_SENESCENCE = "#6B8E9E",
LiteratureMarkers = "#CA7E45"
)
),
mode="simple",
widthlegend=30,
sig_threshold=0.05,
cohen_threshold=0.6,
pointSize=6,
colorPalette="Paired")
Overall_Scores$heatmap
Overall_Scores$volcano
ROC curves assess the discriminatory power of gene set scores, indicating how well a signature separates different experimental or clinical groups. When discriminating between Senescent and Proliferative samples, the REACTOME_CELLULAR_SENESCENCE gene set shows more heterogeneous ROC curves and AUC values across scoring methods, reflecting less consistent performance compared to the HernandezSegura and LiteratureMarkers gene sets.
ROC_Scores(data = counts_example,
metadata = metadata_example,
gene_sets=genesets_example,
method = "all",
variable ="Condition",
colors = c(logmedian = "#3E5587", ssGSEA = "#B65285", ranking = "#B68C52"),
grid = TRUE,
spacing_annotation=0.3,
ncol=NULL,
nrow=1,
mode = "simple",
widthTitle = 28,
titlesize = 10,
title="Marthandan et al. 2016")
Finally, we perform simulations using random gene sets to estimate the false positive rate. This step informs whether observed gene set signals exceed what would be expected by chance, improving confidence in the results. By simulating 10 random gene sets of the same size, we see that the LiteratureMarkers gene set performs best, with only two of the random sets achieving higher effect sizes than the original one using the ranking method. While increasing the number of simulated sets would yield finer resolution, this comes at the cost of additional computational time.
FPR_Simulation(data = counts_example,
metadata = metadata_example,
original_signatures = genesets_example,
gene_list = row.names(counts_example),
number_of_sims = 10,
title = "Marthandan et al. 2016",
widthTitle = 30,
Variable = "Condition",
titlesize = 12,
pointSize = 5,
labsize = 10,
mode = "simple",
ColorValues=NULL,
ncol=NULL,
nrow=1 )
Enrichment-based approach
We now shift focus to enrichment-based methods, which rely on differential expression statistics to evaluate whether gene sets show coordinated changes between conditions. This approach complements score-based methods by explicitly leveraging contrasts between groups and quantifying statistical significance.
The first step is to compute differential expression. In this
example, we build the design matrix automatically from the
Condition
variable in the metadata and define the contrast
Senescent - Proliferative
. Internally, this fits a linear
model without an intercept, enabling quick comparison between levels of
a categorical variable.
# Build design matrix from variables (Condition) and apply a contrast.
# The design matrix is constructed automatically using the variable "Condition".
DEGs <- calculateDE(data = counts_example,
metadata = metadata_example,
variables = "Condition",
contrasts = c("Senescent - Proliferative"))
DEGs$`Senescent-Proliferative`[1:5,]
## logFC AveExpr t P.Value adj.P.Val B
## CCND2 3.816674 4.406721 12.379571 2.944841e-15 2.349757e-12 24.64395
## MKI67 -3.581174 6.605339 -9.187514 2.110510e-11 4.769568e-10 15.91258
## PTCHD4 3.398914 3.556007 10.729126 2.461327e-13 2.868691e-11 20.30116
## KIF20A -3.365481 5.934893 -9.718104 4.406643e-12 1.771710e-10 17.45858
## CDC20 -3.304602 6.104079 -9.791031 3.562991e-12 1.592176e-10 17.66821
Once differential expression is calculated, we can visualize the results with a volcano plot. This plot highlights genes in our predefined gene sets, showing the magnitude of differential expression (log fold-change) versus statistical significance (adjusted p-value). Up- (green) and downregulated (red) genes are color-coded. For the HernandezSegura gene set, green genes mostly appear in the positive logFC range and red genes in the negative range, though these are not the genes with the most extreme logFC values. In the LiteratureMarkers gene set, two red genes (LMNB1 and MKI67) show very negative logFC, which is consistent with their roles as proliferation markers since senescent cells are non-proliferative. Genes from the REACTOME_CELLULAR_SENESCENCE gene set are shown in blue, as this set lacks information on directionality; these genes are scattered across the full range of logFC values.
# Change order: signatures in columns, contrast in rows
plotVolcano(DEGs, genes = genesets_example,
N = NULL,
x = "logFC", y = "-log10(adj.P.Val)", pointSize = 2,
color = "#6489B4", highlightcolor = "#05254A", highlightcolor_upreg = "#038C65", highlightcolor_downreg = "#8C0303",nointerestcolor = "#B7B7B7",
threshold_y = NULL, threshold_x = NULL,
xlab = NULL, ylab = NULL, ncol = NULL, nrow = NULL, title = "Marthandan et al. 2016",
labsize = 10, widthlabs = 24, invert = TRUE)
Next, we perform Gene Set Enrichment Analysis (GSEA) using the differential expression results. This approach evaluates whether members of each gene set are non-randomly distributed toward the top or bottom of the ranked gene list, providing normalized enrichment scores (NES) and adjusted p-values. The ranking of genes can be based on different metrics, reflecting the expected directionality of the gene set. This is also reflected in the plot labels below: gene sets are marked as altered when genes are ordered by the B statistic, or enriched/depleted when ordered by the t-statistic. For a full example and explanation, see the tutorial here.
GSEAresults <- runGSEA(DEGList = DEGs,
gene_sets = genesets_example,
stat = NULL,
ContrastCorrection = FALSE)
GSEAresults
## $`Senescent-Proliferative`
## pathway pval padj log2err ES
## <char> <num> <num> <num> <num>
## 1: HernandezSegura 4.875921e-10 1.462776e-09 0.8012156 0.6263553
## 2: REACTOME_CELLULAR_SENESCENCE 1.207999e-02 1.811999e-02 0.1317424 0.3687292
## 3: LiteratureMarkers 2.020412e-02 2.020412e-02 0.1461515 0.6948634
## NES size leadingEdge stat_used
## <num> <int> <list> <char>
## 1: 2.636116 52 CNTLN, D.... t
## 2: 1.446005 128 H1-2, H2.... B
## 3: 1.650282 7 LMNB1, M.... t
To visualize the enrichment of each gene set along the ranked gene
list, we use enrichment curves from fgsea.
These plots show
the running enrichment score across the ranked list, highlighting where
members of each gene set are concentrated.
plotGSEAenrichment(GSEA_results=GSEAresults,
DEGList=DEGs,
gene_sets=genesets_example,
widthTitle=40, grid = TRUE, titlesize = 10, nrow=1, ncol=3)
We can also summarize enrichment results in a lollipop plot, which compactly shows the normalized enrichment score and highlights statistically significant gene sets. This provides a quick overview of which pathways are most strongly altered. Here, we can see that the HernandezSegura gene set clearly exhibits the strongest enrichment signal.
plotNESlollipop(GSEA_results=GSEAresults,
saturation_value=NULL,
nonsignif_color = "#F4F4F4",
signif_color = "red",
sig_threshold = 0.05,
grid = FALSE,
nrow = NULL, ncol = NULL,
widthlabels=13,
title=NULL, titlesize=12)
## $`Senescent-Proliferative`
Finally, a volcano-style scatter plot combines NES and significance for all gene sets, making it easy to identify which sets show the strongest and most statistically robust alteration.
plotCombinedGSEA(GSEAresults, sig_threshold = 0.05, PointSize=6, widthlegend = 26 )
From this exercise in benchmarking mode, we can see that two gene sets clearly perform best at discriminating between Senescent and Proliferative conditions: HernandezSegura and LiteratureMarkers. The REACTOME_CELLULAR_SENESCENCE gene set does not show a strong signal; since it is undirected and lacks information on up- or downregulation, this also dilutes the signal, highlighting the importance of providing directionality when available.
Across methods, scoring and enrichment approaches provide complementary insights. Score-based methods offer sample-level resolution, capturing strong contributions from individual genes, while enrichment-based methods evaluate coordinated behaviour across the set, typically at the group level. Even among scoring methods, rank-based approaches (e.g., ssGSEA or rank scoring) are generally more robust to technical noise, whereas magnitude-based methods (e.g., log2-median) better detect shifts in well-controlled data. Performance also depends on sample size and gene set size, explaining why HernandezSegura may appear stronger in enrichment analyses and LiteratureMarkers in scoring analyses. Integrating both approaches provides the most comprehensive view of gene set behaviour.
3.4.2 Discovery Mode
In Discovery Mode (full tutorial here), the output focuses on a single gene set:
- Score distributions by phenotype
- Pairwise contrasts (Cohen’s d) and overall effect sizes (Cohen’s f)
- Enrichment score summaries (NES) with adjusted p-values (e.g., lollipop plots)
We start by selecting the gene set of interest. Here, we use the HernandezSegura gene set from our example collection.
HernandezSegura_GeneSet <- list(HernandezSegura=genesets_example$HernandezSegura)
For illustration, we add two hypothetical variables to the metadata:
DaysToSequencing
, representing the number of days between
sample preparation and sequencing, and Researcher
,
identifying the person who processed each sample. This allows us to
explore associations between gene set activity and different types of
variables (categorical or continuous).
data(metadata_example)
set.seed("123456")
metadata_example$Researcher <- sample(c("John","Ana","Francisca"),39, replace = TRUE)
metadata_example$DaysToSequencing <- sample(c(1:20),39, replace = TRUE)
head(metadata_example)
## sampleID DatasetID CellType Condition SenescentType
## 252 SRR1660534 Marthandan2016 Fibroblast Senescent Telomere shortening
## 253 SRR1660535 Marthandan2016 Fibroblast Senescent Telomere shortening
## 254 SRR1660536 Marthandan2016 Fibroblast Senescent Telomere shortening
## 255 SRR1660537 Marthandan2016 Fibroblast Proliferative none
## 256 SRR1660538 Marthandan2016 Fibroblast Proliferative none
## 257 SRR1660539 Marthandan2016 Fibroblast Proliferative none
## Treatment Researcher DaysToSequencing
## 252 PD72 (Replicative senescence) Ana 6
## 253 PD72 (Replicative senescence) Ana 18
## 254 PD72 (Replicative senescence) John 19
## 255 young Ana 2
## 256 young Francisca 9
## 257 young John 10
We can then examine how the selected gene set associates with these
variables using enrichment-based approaches in Discovery Mode. The
resulting plot highlights significant associations across variables and
visually summarizes the direction and strength of the effect. The
“simple” mode provides comparison of effect sizes for pairwise contrasts
between only two levels of the variable, but can be changed to more
levels of comparison (see ?VariableAssociation
). Here, we
can see that the HernandezSegura gene set is significantly enriched in
samples sequenced by Francisca, when compared to those processed by Ana
or John, suggesting that this gene set may be particularly relevant to
her experimental conditions or sample processing. This gene set has also
a strong enrichment for the comparison between proliferative and
senescent, which is expected given the nature of the gene set and the
results from the Benchmarking mode.
VariableAssociation(
data = counts_example,
metadata = metadata_example,
method = "GSEA",
cols = c("Condition","Researcher","DaysToSequencing"),
mode = "simple",
gene_set = HernandezSegura_GeneSet,
saturation_value = NULL,
nonsignif_color = "white",
signif_color = "red",
sig_threshold = 0.05,
widthlabels = 30,
labsize = 10,
titlesize = 14,
pointSize = 5
) $plot
Next, we evaluate the association using a score-based method (here, log2-median). This approach calculates per-sample scores for the gene set and summarizes them across variables. In this analysis, Condition (representing whether the sample is senescent or proliferative) shows a strong effect, reflected by a large Cohen’s f. This effect size metric is directly comparable across different variable types (categorical, numeric, etc.), making it particularly versatile. In contrast, the variable Researcher does not show a significant effect here, unlike in the enrichment analysis. This divergence illustrates the value of applying both enrichment- and score-based approaches in a complementary manner.
VariableAssociation(data = counts_example,
metadata = metadata_example,
method = "logmedian",
cols = c("Condition","Researcher","DaysToSequencing"),
gene_set = HernandezSegura_GeneSet,
mode="simple",
nonsignif_color = "white", signif_color = "red", saturation_value=NULL,sig_threshold = 0.05,
widthlabels=30, labsize=10, titlesize=14, pointSize=5, discrete_colors=NULL,
continuous_color = "#8C6D03", color_palette = "Set2")$Overall
## Variable Cohen_f P_Value
## 1 Condition 1.19407625 1.267423e-08
## 2 Researcher 0.12816159 7.458318e-01
## 3 DaysToSequencing 0.00792836 9.617953e-01
Benchmarking mode offers the most comprehensive set of features and allows users to seamlessly move from discovery to benchmarking mode once a variable of interest has been identified and further testing is required. The main difference from Discovery mode is that Benchmarking is designed to evaluate multiple gene sets simultaneously, whereas Discovery mode focuses on quantifying a single, robust gene set.
3.5. Individual Gene Exploration
To better understand the contribution of individual genes within a
gene set and identify whether specific genes drive the overall signal,
markeR
offers a suite of gene-level exploratory analyses
(via VisualiseIndividualGenes
), including:
- Expression heatmaps of genes across samples and groups
- Violin plots showing expression distributions of individual genes
- Correlation heatmaps to reveal co-expression patterns among genes in the set
- ROC curves and AUC values for individual genes to evaluate their discriminatory power
- Effect size calculations (Cohen’s d) per gene to quantify differential expression
- Principal Component Analysis (PCA) on gene set genes to assess variance explained and sample clustering
See ?VisualiseIndividualGenes
for more details on the
available options and the extended online tutorial here.
3.6. Compare with Reference Gene Sets
markeR
allows comparison of user-defined gene sets to
reference sets (e.g., from MSigDB) using:
Jaccard Index: Measures gene overlap relative to union size.
Log Odds Ratio (logOR): Computes enrichment using a user-defined gene universe and Fisher’s exact test.
Filters can be applied based on similarity thresholds (e.g., minimum Jaccard, OR, or p-value).
Example of Gene Set Similarity (full tutorial here). Here, we are seeing how two user-defined signatures compare to a set of other user-defined signatures, as well as to a collection of reference gene sets from MSigDB (C2:CP:KEGG_LEGACY). The results are visualized in a heatmap, showing the similarity between the signatures based on log odds ratio (at least one with log10OR > 2).
# Example data
signature1 <- c("TP53", "BRCA1", "MYC", "EGFR", "CDK2")
signature2 <- c("ATXN2", "FUS", "MTOR", "CASP3")
signature_list <- list(
"User_Apoptosis" = c("TP53", "CASP3", "BAX"),
"User_CellCycle" = c("CDK2", "CDK4", "CCNB1", "MYC"),
"User_DNARepair" = c("BRCA1", "RAD51", "ATM"),
"User_MTOR" = c("MTOR", "AKT1", "RPS6KB1")
)
geneset_similarity(
signatures = list(Sig1 = signature1, Sig2 = signature2),
other_user_signatures = signature_list,
collection = "C2",
subcollection = "CP:KEGG_LEGACY",
metric = "odds_ratio",
# Define gene universe (e.g., genes from HPA or your dataset)
universe = unique(c(
signature1, signature2,
unlist(signature_list),
msigdbr::msigdbr(species = "Homo sapiens", category = "C2")$gene_symbol
)),
or_threshold = 100, #log10OR = 2
pval_threshold = 0.05,
width_text=50
)$plot
6. Session Information
## R version 4.5.1 (2025-06-13)
## Platform: x86_64-pc-linux-gnu
## Running under: Ubuntu 24.04.2 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 LC_TIME=C.UTF-8
## [4] LC_COLLATE=C.UTF-8 LC_MONETARY=C.UTF-8 LC_MESSAGES=C.UTF-8
## [7] LC_PAPER=C.UTF-8 LC_NAME=C LC_ADDRESS=C
## [10] LC_TELEPHONE=C LC_MEASUREMENT=C.UTF-8 LC_IDENTIFICATION=C
##
## time zone: UTC
## tzcode source: system (glibc)
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] markeR_0.99.3
##
## loaded via a namespace (and not attached):
## [1] pROC_1.19.0.1 gridExtra_2.3 rlang_1.1.6
## [4] magrittr_2.0.3 clue_0.3-66 GetoptLong_1.0.5
## [7] msigdbr_25.1.1 matrixStats_1.5.0 compiler_4.5.1
## [10] mgcv_1.9-3 png_0.1-8 systemfonts_1.2.3
## [13] vctrs_0.6.5 reshape2_1.4.4 stringr_1.5.1
## [16] pkgconfig_2.0.3 shape_1.4.6.1 crayon_1.5.3
## [19] fastmap_1.2.0 backports_1.5.0 labeling_0.4.3
## [22] effectsize_1.0.1 rmarkdown_2.29 ragg_1.4.0
## [25] purrr_1.1.0 xfun_0.53 cachem_1.1.0
## [28] jsonlite_2.0.0 BiocParallel_1.42.1 broom_1.0.9
## [31] parallel_4.5.1 cluster_2.1.8.1 R6_2.6.1
## [34] bslib_0.9.0 stringi_1.8.7 RColorBrewer_1.1-3
## [37] limma_3.64.3 car_3.1-3 jquerylib_0.1.4
## [40] Rcpp_1.1.0 assertthat_0.2.1 iterators_1.0.14
## [43] knitr_1.50 parameters_0.28.0 IRanges_2.42.0
## [46] splines_4.5.1 Matrix_1.7-3 tidyselect_1.2.1
## [49] abind_1.4-8 yaml_2.3.10 doParallel_1.0.17
## [52] codetools_0.2-20 curl_7.0.0 lattice_0.22-7
## [55] tibble_3.3.0 plyr_1.8.9 withr_3.0.2
## [58] bayestestR_0.16.1 evaluate_1.0.4 desc_1.4.3
## [61] circlize_0.4.16 pillar_1.11.0 ggpubr_0.6.1
## [64] carData_3.0-5 foreach_1.5.2 stats4_4.5.1
## [67] insight_1.4.0 generics_0.1.4 S4Vectors_0.46.0
## [70] ggplot2_3.5.2 scales_1.4.0 glue_1.8.0
## [73] tools_4.5.1 data.table_1.17.8 fgsea_1.34.2
## [76] locfit_1.5-9.12 ggsignif_0.6.4 babelgene_22.9
## [79] fs_1.6.6 fastmatch_1.1-6 cowplot_1.2.0
## [82] grid_4.5.1 tidyr_1.3.1 datawizard_1.2.0
## [85] edgeR_4.6.3 colorspace_2.1-1 nlme_3.1-168
## [88] Formula_1.2-5 cli_3.6.5 textshaping_1.0.1
## [91] ComplexHeatmap_2.24.1 dplyr_1.1.4 gtable_0.3.6
## [94] ggh4x_0.3.1 rstatix_0.7.2 sass_0.4.10
## [97] digest_0.6.37 BiocGenerics_0.54.0 ggrepel_0.9.6
## [100] rjson_0.2.23 htmlwidgets_1.6.4 farver_2.1.2
## [103] htmltools_0.5.8.1 pkgdown_2.1.3 lifecycle_1.0.4
## [106] GlobalOptions_0.1.2 statmod_1.5.0