Advanced usage of BioCor
Lluís Revilla
August Pi i Sunyer Biomedical Research Institute (IDIBAPS); Liver Unit, Hospital Cliniclrevilla@clinic.cat
Juan José Lozano
Centro de Investigación Biomédica en Red de Enfermedades Hepaticas y Digestivas (CIBEREHD); Barcelona, SpainPau Sancho-Bru
August Pi i Sunyer Biomedical Research Institute (IDIBAPS); Liver Unit, Hospital Clinic23 April 2024
Source:vignettes/BioCor_2_advanced.Rmd
BioCor_2_advanced.RmdAbstract
Describes how to use the BioCor package to answer several biological questions and how to use functional similarities with other measures.
Introduction
In this vignette we assume that the reader is already familiarized with the introduction vignette, but wants to know how can it help to answer other questions in other situations
We follow the same convention about names of the objects used
genesReact and genesKegg are the list with
information about the pathways the genes are involved in.
Merging similarities
If one calculates similarities with KEGG data and Reactome or other input for the same genes or clusters BioCor provides a couple of functions to merge them.
We can set a weight to each similarity input with
weighted.sum, multiply them also using a weight for each
similarity (with weighted.prod), doing the mean or just
adding them up. Similarities allow us to apply a function to combine the
matrices of a list. Here we use some of the genes used in the first
vignette:
kegg <- mgeneSim(c("672", "675", "10"), genesKegg)
react <- mgeneSim(c("672", "675", "10"), genesReact)
## We can sum it adding a weight to each origin
weighted.sum(c(kegg["672", "675"], react["672", "675"]), w = c(0.3, 0.7))
## [1] 0.7247289
## Or if we want to perform for all the matrix
## A list of matrices to merge
sim <- list("kegg" = kegg, "react" = react)
similarities(sim, weighted.sum, w = c(0.3, 0.7))
## 672 675 10
## 672 1.0000000 0.7247289 0.1428635
## 675 0.7247289 1.0000000 0.1453360
## 10 0.1428635 0.1453360 1.0000000
similarities(sim, weighted.prod, w = c(0.3, 0.7))
## 672 675 10
## 672 0.2100000 0.0173101952 0.0000000000
## 675 0.0173102 0.2100000000 0.0003532339
## 10 0.0000000 0.0003532339 0.2100000000
similarities(sim, prod)
## 672 675 10
## 672 1.0000000 0.082429501 0.000000000
## 675 0.0824295 1.000000000 0.001682066
## 10 0.0000000 0.001682066 1.000000000
similarities(sim, mean)
## 672 675 10
## 672 1.0000000 0.5412148 0.1020454
## 675 0.5412148 1.0000000 0.1061662
## 10 0.1020454 0.1061662 1.0000000This functions are similar to weighted.mean, except that
first the multiplication by the weights is done and then the
NAs are removed:
weighted.mean(c(1, NA), w = c(0.5, 0.5), na.rm = TRUE)
## [1] 1
weighted.mean(c(1, 0.5, NA), w = c(0.5, 0.25, 0.25), na.rm = TRUE)
## [1] 0.8333333
weighted.sum(c(1, NA), w = c(0.5, 0.5))
## [1] 0.5
weighted.sum(c(1, 0.5, NA), w = c(0.5, 0.25, 0.25))
## [1] 0.625
weighted.prod(c(1, NA), w = c(0.5, 0.5))
## [1] 0.5
weighted.prod(c(1, 0.5, NA), w = c(0.5, 0.25, 0.25))
## [1] 0.0625Assessing a differential study
In this example we will use the functional similarities in a classical differential study.
Obtaining data
We start using data from the RNAseq workflow and following the analysis comparing treated and untreated:
suppressPackageStartupMessages(library("airway"))
data("airway")
library("DESeq2")
dds <- DESeqDataSet(airway, design = ~ cell + dex)
dds$dex <- relevel(dds$dex, "untrt")
dds <- DESeq(dds)
## estimating size factors
## estimating dispersions
## gene-wise dispersion estimates
## mean-dispersion relationship
## final dispersion estimates
## fitting model and testing
res <- results(dds, alpha = 0.05)
summary(res)
##
## out of 33469 with nonzero total read count
## adjusted p-value < 0.05
## LFC > 0 (up) : 2211, 6.6%
## LFC < 0 (down) : 1817, 5.4%
## outliers [1] : 0, 0%
## low counts [2] : 16687, 50%
## (mean count < 7)
## [1] see 'cooksCutoff' argument of ?results
## [2] see 'independentFiltering' argument of ?results
plot(res$log2FoldChange, -log10(res$padj),
pch = 16, xlab = "log2FC",
ylab = "-log10(p.ajd)", main = "Untreated vs treated"
)
logFC <- 2.5
abline(v = c(-logFC, logFC), h = -log10(0.05), col = "red")
Volcano plot. The airway data
As we can see here there are around 4000 genes differentially expressed genes, some of them differentially expressed above 2^2.5.
Selecting differentially expressed genes
Usually in such a study one selects genes above certain logFC or fold change threshold, here we use the absolute value of 2.5:
fc <- res[abs(res$log2FoldChange) >= logFC & !is.na(res$padj), ]
fc <- fc[fc$padj < 0.05, ]
# Convert Ids (used later)
genes <- select(org.Hs.eg.db,
keys = rownames(res), keytype = "ENSEMBL",
column = c("ENTREZID", "SYMBOL")
)
## 'select()' returned 1:many mapping between keys and columns
genesFC <- genes[genes$ENSEMBL %in% rownames(fc), ]
genesFC <- genesFC[!is.na(genesFC$ENTREZID), ]
genesSim <- genesFC$ENTREZID
names(genesSim) <- genesFC$SYMBOL
genesSim <- genesSim[!duplicated(genesSim)]
# Calculate the functional similarity
sim <- mgeneSim(genes = genesSim, info = genesReact, method = "BMA")
## Warning in mgeneSim(genes = genesSim, info = genesReact, method = "BMA"): Some
## genes are not in the list provided.Once the similarities for the selected genes are calculated we can now visualize the effect of each method:
nas <- apply(sim, 1, function(x) {
all(is.na(x))
})
sim <- sim[!nas, !nas]
MDSs <- cmdscale(1 - sim)
plot(MDSs, type = "n", main = "BMA similarity", xlab = "MDS1", ylab = "MDS2")
up <- genes[genes$ENSEMBL %in% rownames(fc)[fc$log2FoldChange >= logFC], "SYMBOL"]
text(MDSs, labels = rownames(MDSs), col = ifelse(rownames(MDSs) %in% up, "black", "red"))
abline(h = 0, v = 0)
legend("top", legend = c("Up-regulated", "Down-regulated"), fill = c("black", "red"))
Functional similarity of genes with logFC above 2,5. Similar genes cluster together.
This plot illustrate that some differentially expressed genes are quite similar according to their pathways. Suggesting that there might be a relationship between them. Furthermore, some up-regulated genes seem functionally related to down-regulated genes indicating a precise regulation of the pathways where they are involved.
Note that here we are only using 74 genes from the original 127.
Are differentially expressed genes selected by their functionality?
In the previous section we have seen that some differentially expressed genes are functionally related and that they have a high logFC value. Are genes differentially expressed more functional related than those which aren’t differential expressed?
For simplicity we will use a subset of 400 genes represented again in a volcano plot and we will look for the functional similarities between those genes:
set.seed(250)
subRes <- res[!is.na(res$log2FoldChange), ]
subs <- sample.int(nrow(subRes), size = 400)
subRes <- subRes[subs, ]
g <- genes[genes$ENSEMBL %in% rownames(subRes), "ENTREZID"]
gS <- mgeneSim(g[g %in% names(genesReact)], genesReact, "BMA")
deg <- rownames(subRes[subRes$padj < 0.05 & !is.na(subRes$padj), ])
keep <- rownames(gS) %in% genes[genes$ENSEMBL %in% deg, "ENTREZID"]
plot(subRes$log2FoldChange, -log10(subRes$padj),
pch = 16, xlab = "log2FC",
ylab = "-log10(p.ajd)", main = "Untreated vs treated"
)
abline(v = c(-logFC, logFC), h = -log10(0.05), col = "red")
Volcano plot of the subset of 400 genes. This subset will be used in the following sections
We can answer this by testing it empirically:
library("boot")
# The mean of genes differentially expressed
(scoreDEG <- mean(gS[!keep, keep], na.rm = TRUE))
## [1] 0.1226459
b <- boot(data = gS, R = 1000, statistic = function(x, i) {
g <- !rownames(x) %in% rownames(x)[i]
mean(x[g, i], na.rm = TRUE)
})
(p.val <- (1 + sum(b$t > scoreDEG)) / 1001)
## [1] 0.6633367
hist(b$t, main = "Distribution of scores", xlab = "Similarity score")
abline(v = scoreDEG, col = "red")
Distribution of scores between differentially expressed genes and those who aren’t. The line indicates the mean score of the similarity between differentially expressed genes and those which aren’t differentially expressed.
Comparing the genes differentially expressed and those who aren’t doesn’t show that they are non-randomly selected (p-value 0.6633367). The genes with a p-value below the threshold are not more closely functionally related than all the other genes1.
Are functionally related the selected differentially expressed genes?
We have seen that the genes differentially expressed aren’t selected by their functionality. However they could be more functionally related than the other genes. Are the differentially expressed genes more functionally similar than it would be expected ?
(scoreW <- combineScores(gS[keep, keep], "avg"))
## [1] 0.158997
b <- boot(data = gS, R = 1000, statistic = function(x, i) {
mean(x[i, i], na.rm = TRUE)
})
(p.val <- (1 + sum(b$t > scoreW)) / 1001) # P-value
## [1] 0.003996004
hist(b$t, main = "Distribution of scores", xlab = "Similarity score")
abline(v = scoreW, col = "red")
Distribution of the similarity within differentially expressed genes. The line indicates the mean funtional similarity whitin them.
If we selected randomly the genes from our pool we would expect a score around 0.158997 with a probability of 0.003996. That means that the differentially expressed genes is highly different compared to the other genes if we use a significance threshold of 0.05 2.
Influence of the fold change in the functionally similarity of the genes
We have seen that the genes differentially expressed are not selected by functional similarity but they are functionally related. Now we would like to know if selecting a fold change threshold affects the functional similarity between them.
To know the relationship between the fold change and the similarity between genes we have several methods:
s <- seq(0, max(abs(subRes$log2FoldChange)) + 0.05, by = 0.05)
l <- sapply(s, function(x) {
deg <- rownames(subRes[abs(subRes$log2FoldChange) >= x, ])
keep <- rownames(gS) %in% genes[genes$ENSEMBL %in% deg, "ENTREZID"]
BetweenAbove <- mean(gS[keep, keep], na.rm = TRUE)
AboveBelow <- mean(gS[keep, !keep], na.rm = TRUE)
BetweenBelow <- mean(gS[!keep, !keep], na.rm = TRUE)
c(
"BetweenAbove" = BetweenAbove, "AboveBelow" = AboveBelow,
"BetweenBelow" = BetweenBelow
)
})
L <- as.data.frame(cbind(logfc = s, t(l)))
plot(L$logfc, L$BetweenAbove,
type = "l", xlab = "abs(log2) fold change",
ylab = "Similarity score",
main = "Similarity scores along logFC", col = "darkred"
)
lines(L$logfc, L$AboveBelow, col = "darkgreen")
lines(L$logfc, L$BetweenBelow, col = "black")
legend("topleft",
legend = c(
"Between genes above and below threshold",
"Whitin genes above threshold",
"Whitin genes below threshold"
),
fill = c("darkgreen", "darkred", "black")
)
Similarity of genes along abs(logFC). Assessing the similarity of genes according to their absolute log2 fold change.
The functional similarity of the genes above the threshold increases with a more restrictive threshold, indicating that a logFC threshold selects genes by their functionality, or in other words that genes differentially expressed tend to be of related pathways. The similarity between those genes above the threshold and below remains constant as well as within genes below the threshold.
l <- sapply(s, function(x) {
# Names of genes up and down regulated
degUp <- rownames(subRes[subRes$log2FoldChange >= x, ])
degDown <- rownames(subRes[subRes$log2FoldChange <= -x, ])
# Translate to ids in gS
keepUp <- rownames(gS) %in% genes[genes$ENSEMBL %in% degUp, "ENTREZID"]
keepDown <- rownames(gS) %in% genes[genes$ENSEMBL %in% degDown, "ENTREZID"]
# Calculate the mean similarity between each subgrup
between <- mean(gS[keepUp, keepDown], na.rm = TRUE)
c("UpVsDown" = between)
})
L <- as.data.frame(cbind("logfc" = s, "UpVsDown" = l))
plot(L$logfc, L$UpVsDown,
type = "l",
xlab = "abs(log2) fold change threshold",
ylab = "Similarity score",
main = "Similarity scores along logFC"
)
legend("topright",
legend = "Up vs down regulated genes",
fill = "black"
)
Functional similarity between the up-regulated and down-regulated genes.
The maximal functional similarity between genes up-regulated and down-regulated are at 1.15 log2 fold change.
Assessing a new pathway
Sometimes the top differentially expressed genes or some other key
genes are selected as a signature or a potential new group of related
genes. In those cases we can test how does the network of genes change
if we add them. Here we create a new pathway named deg, and
we see the effect on the functional similarity score for all the
genes:
# Adding a new pathway "deg" to those genes
genesReact2 <- genesReact
diffGenes <- genes[genes$ENSEMBL %in% deg, "ENTREZID"]
genesReact2[diffGenes] <- sapply(genesReact[diffGenes], function(x) {
c(x, "deg")
})
plot(ecdf(mgeneSim(names(genesReact), genesReact)))
curve(ecdf(mgeneSim(names(genesReact2), genesReact2)), color = "red")This would take lot of time, for a illustration purpose we reduce to some genes:
library("Hmisc")
genesReact2 <- genesReact
diffGenes <- genes[genes$ENSEMBL %in% deg, "ENTREZID"]
# Create the new pathway called deg
genesReact2[diffGenes] <- sapply(genesReact[diffGenes], function(x) {
c(x, "deg")
})
ids <- unique(genes[genes$ENSEMBL %in% rownames(subRes), "ENTREZID"])
Ecdf(c(mgeneSim(ids, genesReact, method = "BMA"),
mgeneSim(ids, genesReact2, method = "BMA")
),
group = c(rep("Reactome", length(ids)^2), rep("Modified", length(ids)^2)),
col = c("black", "red"), xlab = "Functional similarities",
main = "Empirical cumulative distribution")
The effect of adding a new pathway to a functional similarity. In red the same network as in black but with the added pathway.
Merging sources of information
Sometimes we have several origin of information, either several
databases, or information from other programs… We can merge this in the
single object required by the function in BioCor using the
function combineSources3.
This functions helps to evaluate what happens when we add more pathway information. For instance here we add the information in Kegg and the information in Reactome and we visualize it using the same procedure as previously:
genesKegg <- as.list(org.Hs.egPATH)
gSK <- mgeneSim(rownames(gS), genesKegg)
mix <- combineSources(genesKegg, genesReact)
## Warning in combineSources(genesKegg, genesReact): More than 50% of genes identifiers of a source are unique
## Check the identifiers of the genes
gSMix <- mgeneSim(rownames(gS), mix)
Ecdf(c(gS, gSK, gSMix),
group = c(
rep("Reactome", length(gS)), rep("Kegg", length(gSK)),
rep("Mix", length(gSMix))
),
col = c("black", "red", "blue"), xlab = "Functional similarities",
main = "Empirical cumulative distribution."
)
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
Comparison of functional similarity in different databases.
When mixed, there is a huge increase in the genes that share a
pathway using the max method. Observe in next figure
(@ref(fig:combineSource2)) how does the method affects to the
results:
gSK2 <- mgeneSim(rownames(gS), genesKegg, method = "BMA")
gS2 <- mgeneSim(rownames(gS), genesReact, method = "BMA")
gSMix2 <- mgeneSim(rownames(gS), mix, method = "BMA")
Ecdf(c(gS2, gSK2, gSMix2),
group = c(
rep("Reactome", length(gS)), rep("Kegg", length(gSK)),
rep("Mix", length(gSMix))
),
col = c("black", "red", "blue"), xlab = "Functional similarities (BMA)", main = "Empirical cumulative distribution."
)
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
Comparison of functional similarity in different gene sets.
Now we can appreciate that most of the functional similarity is brought by Reactome database
miRNA analysis
miRNAs are RNAs that interact with many genes and transcripts and is subject to change with location and time, thus defining the effect of an miRNA is difficult. In this section we try to answer how functionally similar are miRNA between them to provide a help for potentially closely related miRNA.
First we look for human miRNAs and prepare them for using them as the input for the cluster functions (Restricting to 10 miRNA with less than 150 genes).4
library("targetscan.Hs.eg.db")
## select human mirna
humanMirnaIdx <- grep("hsa", mappedkeys(targetscan.Hs.egMIRNA))
## select seed-based families for human mirna
humanMirna <- mappedkeys(targetscan.Hs.egMIRNA)[humanMirnaIdx]
## select targets of families
humanMirnaFamilies <- unlist(mget(humanMirna, targetscan.Hs.egMIRBASE2FAMILY))
humanMirnaTargets <- mget(humanMirnaFamilies, revmap(targetscan.Hs.egTARGETS))
names(humanMirnaTargets) <- humanMirna
# Restrict to miRNA with more than one target and less than 150
miRNAs <- sample(humanMirnaTargets[lengths(humanMirnaTargets) > 1 &
lengths(humanMirnaTargets) < 150], 10)
lengths(miRNAs)
## hsa-miR-615-3p hsa-miR-487b hsa-miR-210 hsa-miR-450a hsa-miR-191
## 28 23 70 24 118
## hsa-miR-187 hsa-miR-551a hsa-miR-99a hsa-miR-184 hsa-miR-99b
## 13 17 92 54 92Now we calculate the functional similarity of those miRNAs using a cluster approach.
cluster1 <- mclusterSim(miRNAs, genesReact, method = "BMA")
## Warning in mclusterSim(miRNAs, genesReact, method = "BMA"): Some genes are not
## in the list provided.
knitr::kable(round(cluster1, 2), caption = "The similarity between miRNA", format = "html")| hsa-miR-615-3p | hsa-miR-487b | hsa-miR-210 | hsa-miR-450a | hsa-miR-191 | hsa-miR-187 | hsa-miR-551a | hsa-miR-99a | hsa-miR-184 | hsa-miR-99b | |
|---|---|---|---|---|---|---|---|---|---|---|
| hsa-miR-615-3p | 1.00 | 0.72 | 0.67 | 0.71 | 0.60 | 0.55 | 0.76 | 0.57 | 0.48 | 0.57 |
| hsa-miR-487b | 0.72 | 1.00 | 0.64 | 0.74 | 0.61 | 0.60 | 0.75 | 0.63 | 0.53 | 0.63 |
| hsa-miR-210 | 0.67 | 0.64 | 1.00 | 0.70 | 0.80 | 0.57 | 0.71 | 0.66 | 0.57 | 0.66 |
| hsa-miR-450a | 0.71 | 0.74 | 0.70 | 1.00 | 0.62 | 0.58 | 0.79 | 0.68 | 0.62 | 0.68 |
| hsa-miR-191 | 0.60 | 0.61 | 0.80 | 0.62 | 1.00 | 0.57 | 0.59 | 0.76 | 0.66 | 0.76 |
| hsa-miR-187 | 0.55 | 0.60 | 0.57 | 0.58 | 0.57 | 1.00 | 0.59 | 0.69 | 0.75 | 0.69 |
| hsa-miR-551a | 0.76 | 0.75 | 0.71 | 0.79 | 0.59 | 0.59 | 1.00 | 0.58 | 0.55 | 0.58 |
| hsa-miR-99a | 0.57 | 0.63 | 0.66 | 0.68 | 0.76 | 0.69 | 0.58 | 1.00 | 0.78 | 1.00 |
| hsa-miR-184 | 0.48 | 0.53 | 0.57 | 0.62 | 0.66 | 0.75 | 0.55 | 0.78 | 1.00 | 0.78 |
| hsa-miR-99b | 0.57 | 0.63 | 0.66 | 0.68 | 0.76 | 0.69 | 0.58 | 1.00 | 0.78 | 1.00 |
So for instance hsa-miR-99b, hsa-miR-99a are functionally related despite being from different families.
Comparing with GO similarities
As suggested in the main vignette functional similarities can be compared to semantic similarities such as those based on GO. Here a comparison using the biological process from the gene ontologies is done:
library("GOSemSim")
##
## GOSemSim v2.28.1 For help: https://yulab-smu.top/biomedical-knowledge-mining-book/
##
## If you use GOSemSim in published research, please cite:
## - Guangchuang Yu. Gene Ontology Semantic Similarity Analysis Using GOSemSim. In: Kidder B. (eds) Stem Cell Transcriptional Networks. Methods in Molecular Biology, 2020, 2117:207-215. Humana, New York, NY. doi:10.1007/978-1-0716-0301-7_11
## - Guangchuang Yu, Fei Li, Yide Qin, Xiaochen Bo, Yibo Wu, Shengqi Wang. GOSemSim: an R package for measuring semantic similarity among GO terms and gene products Bioinformatics 2010, 26(7):976-978. doi:10.1093/bioinformatics/btq064
##
## Attaching package: 'GOSemSim'
## The following objects are masked from 'package:BioCor':
##
## clusterSim, combineScores, geneSim, mclusterSim, mgeneSim
BP <- godata("org.Hs.eg.db", ont = "BP", computeIC = TRUE)
## preparing gene to GO mapping data...
## preparing IC data...
gsGO <- GOSemSim::mgeneSim(rownames(gS), semData = BP, measure = "Resnik", verbose = FALSE)
keep <- rownames(gS) %in% rownames(gsGO)
hist(as.dist(gS[keep, keep] - gsGO),
main = "Difference between functional similarity and biological process",
xlab = "Functional similarity - biological process similarity"
)
Comparison of similarities. Functional similarities compared to biological process semantic similarity.
On this graphic we can observe that some genes have a large functional similarity and few biological similarity. They are present together in several pathways while they share few biological process, indicating that they might be key elements of the pathways they are in. On the other hand some other pairs of genes show higher biological process similarity than functional similarity, indicating specialization or compartmentalization of said genes.
Session Info
sessionInfo()
## R version 4.3.3 (2024-02-29)
## Platform: x86_64-pc-linux-gnu (64-bit)
## Running under: Ubuntu 22.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.20.so; LAPACK version 3.10.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] stats4 stats graphics grDevices utils datasets methods
## [8] base
##
## other attached packages:
## [1] GOSemSim_2.28.1 targetscan.Hs.eg.db_0.6.1
## [3] Hmisc_5.1-2 boot_1.3-29
## [5] DESeq2_1.42.1 airway_1.22.0
## [7] SummarizedExperiment_1.32.0 GenomicRanges_1.54.1
## [9] GenomeInfoDb_1.38.8 MatrixGenerics_1.14.0
## [11] matrixStats_1.3.0 BioCor_1.27.2
## [13] org.Hs.eg.db_3.18.0 AnnotationDbi_1.64.1
## [15] IRanges_2.36.0 S4Vectors_0.40.2
## [17] Biobase_2.62.0 BiocGenerics_0.48.1
## [19] BiocStyle_2.30.0
##
## loaded via a namespace (and not attached):
## [1] DBI_1.2.2 bitops_1.0-7 gridExtra_2.3
## [4] GSEABase_1.64.0 rlang_1.1.3 magrittr_2.0.3
## [7] compiler_4.3.3 RSQLite_2.3.6 png_0.1-8
## [10] systemfonts_1.0.6 vctrs_0.6.5 reactome.db_1.86.2
## [13] stringr_1.5.1 pkgconfig_2.0.3 crayon_1.5.2
## [16] fastmap_1.1.1 backports_1.4.1 XVector_0.42.0
## [19] utf8_1.2.4 rmarkdown_2.26 graph_1.80.0
## [22] ragg_1.3.0 purrr_1.0.2 bit_4.0.5
## [25] xfun_0.43 zlibbioc_1.48.2 cachem_1.0.8
## [28] jsonlite_1.8.8 blob_1.2.4 highr_0.10
## [31] DelayedArray_0.28.0 BiocParallel_1.36.0 parallel_4.3.3
## [34] cluster_2.1.6 R6_2.5.1 stringi_1.8.3
## [37] bslib_0.7.0 rpart_4.1.23 jquerylib_0.1.4
## [40] Rcpp_1.0.12 bookdown_0.39 knitr_1.46
## [43] base64enc_0.1-3 Matrix_1.6-5 nnet_7.3-19
## [46] rstudioapi_0.16.0 abind_1.4-5 yaml_2.3.8
## [49] codetools_0.2-19 lattice_0.22-5 tibble_3.2.1
## [52] KEGGREST_1.42.0 evaluate_0.23 foreign_0.8-86
## [55] desc_1.4.3 Biostrings_2.70.3 pillar_1.9.0
## [58] BiocManager_1.30.22 checkmate_2.3.1 RCurl_1.98-1.14
## [61] ggplot2_3.5.0 munsell_0.5.1 scales_1.3.0
## [64] xtable_1.8-4 glue_1.7.0 tools_4.3.3
## [67] data.table_1.15.4 annotate_1.80.0 locfit_1.5-9.9
## [70] fs_1.6.3 XML_3.99-0.16.1 grid_4.3.3
## [73] colorspace_2.1-0 GenomeInfoDbData_1.2.11 htmlTable_2.4.2
## [76] Formula_1.2-5 cli_3.6.2 textshaping_0.3.7
## [79] fansi_1.0.6 S4Arrays_1.2.1 gtable_0.3.5
## [82] yulab.utils_0.1.4 sass_0.4.9 digest_0.6.35
## [85] SparseArray_1.2.4 htmlwidgets_1.6.4 memoise_2.0.1
## [88] htmltools_0.5.8.1 pkgdown_2.0.9 lifecycle_1.0.4
## [91] httr_1.4.7 GO.db_3.18.0 bit64_4.0.5