The methods of group similarity implemented in simona are mainly from the supplementary file of the paper “Mazandu et al., Gene Ontology semantic similarity tools: survey on features and challenges for biological knowledge discovery. Briefings in Bioinformatics 2017”. Original denotations have been slightly modified to make them more consistent. Also more explanations have been added in this vignette.

There are two groups of terms denoted as TpT_p and TqT_q represented as two sets:

$$ T_p = \{ a_1, a_2, ...\} \\ T_q = \{ b_1, b_2, ... \} $$

where aia_i is a term in set TpT_p and bjb_j is a term in set TqT_q.

The wrapper function group_sim() calculates semantic similarities between two groups of terms with a specific method. Note the method name can be partially matched.

group_sim(dag, group1, group2, method = ..., control = list(...))

Some of the group similarity methods have no assumption of which similarity measure between single terms to use. If there are annotation already provided in the DAG object, by default Sim_Lin_1998 is used, or else Sim_WP_1994 is used. The term similarity method can be set via the term_sim_method parameter in control. Additionally parameters for a specific term_sim_method can also be set in control.

group_sim(dag, group1, group2, method = ..., 
    control = list(term_sim_method = ...))

All supported group similarity methods are:

##  [1] "GroupSim_pairwise_avg"            "GroupSim_pairwise_max"           
##  [3] "GroupSim_pairwise_BMA"            "GroupSim_pairwise_BMM"           
##  [5] "GroupSim_pairwise_ABM"            "GroupSim_pairwise_HDF"           
##  [7] "GroupSim_pairwise_MHDF"           "GroupSim_pairwise_VHDF"          
##  [9] "GroupSim_pairwise_Froehlich_2007" "GroupSim_pairwise_Joeng_2014"    
## [11] "GroupSim_SimALN"                  "GroupSim_SimGIC"                 
## [13] "GroupSim_SimDIC"                  "GroupSim_SimUIC"                 
## [15] "GroupSim_SimUI"                   "GroupSim_SimDB"                  
## [17] "GroupSim_SimUB"                   "GroupSim_SimNTO"                 
## [19] "GroupSim_SimCOU"                  "GroupSim_SimCOT"                 
## [21] "GroupSim_SimLP"                   "GroupSim_Ye_2005"                
## [23] "GroupSim_SimCHO"                  "GroupSim_SimALD"                 
## [25] "GroupSim_Jaccard"                 "GroupSim_Dice"                   
## [27] "GroupSim_Overlap"                 "GroupSim_Kappa"

Pairwise term similarity-based methods

GroupSim_pairwise_avg

Denote S(a,b)S(a, b) as the semantic similarity between term aa and bb where aa is from group pp and bb is from group qq, The similarity between group pp and group qq is the average similarity of every pair of individual terms in the two groups:

GroupSim(p,q)=1|Tp|*|Tq|aTp,bTqS(a,b) \mathrm{GroupSim}(p, q) = \frac{1}{|T_p|*|T_q|} \sum_{a \in T_p, b \in T_q}S(a, b)

The term semantic similarity method and the IC method can be set via control argument, for example:

group_sim(dag, group1, group2, method = "GroupSim_pairwise_avg"
    control = list(term_sim_method = "Sim_Lin_1998", IC_method = "IC_annotation")`.

Other parameters for the term_sim_method can also be set in the control list.

Paper link: https://doi.org/10.1093/bioinformatics/btg153.

GroupSim_pairwise_max

The similarity is defined as the maximal S(a,b)S(a, b) among all pairs of terms in group pp and qq:

GroupSim(p,q)=maxaTp,bTqS(a,b) \mathrm{GroupSim}(p, q) = \max_{a \in T_p, b \in T_q}S(a, b)

The term semantic similarity method and the IC method can be set via control argument, for example:

group_sim(dag, group1, group2, method = "GroupSim_pairwise_max"
    control = list(term_sim_method = "Sim_Lin_1998", IC_method = "IC_annotation")`.

Other parameters for the term_sim_method can also be set in the control list.

Paper link: https://doi.org/10.1109/TCBB.2005.50.

GroupSim_pairwise_BMA

BMA stands for “best-match average”. First define similarity of a term xx to a group of terms TT as

S(x,T)=maxyTS(x,y) S(x, T) = \max_{y \in T} S(x, y)

which corresponds to the most similar term in TT to xx. Then the BMA similarity is calculated as:

GroupSim(p,q)=12(1|Tp|aTpS(a,Tq)+1|Tq|bTqS(b,Tp)) \mathrm{GroupSim}(p, q) = \frac{1}{2}\left( \frac{1}{|T_p|}\sum_{a \in T_p} S(a, T_q) + \frac{1}{|T_q|}\sum_{b \in T_q} S(b, T_p) \right)

The term semantic similarity method and the IC method can be set via control argument, for example:

group_sim(dag, group1, group2, method = "GroupSim_pairwise_BMA"
    control = list(term_sim_method = "Sim_Lin_1998", IC_method = "IC_annotation")`.

Other parameters for the term_sim_method can also be set in the control list.

Paper link: https://doi.org/10.1155/2012/975783.

GroupSim_pairwise_BMM

BMM stands for “best-match max”. It is defined as:

GroupSim(p,q)=max{1|Tp|aTpS(a,Tq),1|Tq|bTqS(b,Tp)} \mathrm{GroupSim}(p, q) = \max \left \{ \frac{1}{|T_p|}\sum_{a \in T_p} S(a, T_q), \frac{1}{|T_q|}\sum_{b \in T_q} S(b, T_p) \right \}

The term semantic similarity method and the IC method can be set via control argument, for example:

group_sim(dag, group1, group2, method = "GroupSim_pairwise_BMM"
    control = list(term_sim_method = "Sim_Lin_1998", IC_method = "IC_annotation")`.

Other parameters for the term_sim_method can also be set in the control list.

Paper link: https://doi.org/10.1186/1471-2105-7-302.

GroupSim_pairwise_ABM

ABM stands for “average best-match”. It is defined as:

GroupSim(p,q)=1|Tq|+|Tq|(aTpS(a,Tq)+bTqS(b,Tp)) \mathrm{GroupSim}(p, q) = \frac{1}{|T_q| + |T_q|} \left( \sum_{a \in T_p} S(a, T_q) + \sum_{b \in T_q} S(b, T_p) \right)

The term semantic similarity method and the IC method can be set via control argument, for example:

group_sim(dag, group1, group2, method = "GroupSim_pairwise_ABM"
    control = list(term_sim_method = "Sim_Lin_1998", IC_method = "IC_annotation")`.

Other parameters for the term_sim_method can also be set in the control list.

Paper link: https://doi.org/10.1186/1471-2105-14-284.

GroupSim_pairwise_HDF

First define the distance of a term xx to a group of terms TT:

D(x,T)=1S(x,T)D(x, T) = 1 - S(x, T)

Then the Hausdorff distance between two groups are:

HDF(p,q)=max{maxaTpD(a,Tq),maxbTqD(b,Tq)} \mathrm{HDF}(p, q) = \max \left\{ \max_{a \in T_p} D(a, T_q), \max_{b \in T_q} D(b, T_q) \right\}

This final similarity is:

GroupSim(p,q)=1HDF(p,q) \mathrm{GroupSim}(p, q) = 1 - \mathrm{HDF}(p, q)

The term semantic similarity method and the IC method can be set via control argument, for example:

group_sim(dag, group1, group2, method = "GroupSim_pairwise_HDF"
    control = list(term_sim_method = "Sim_Lin_1998", IC_method = "IC_annotation")`.

Other parameters for the term_sim_method can also be set in the control list.

GroupSim_pairwise_MHDF

Instead of using the maximal distance from a group to the other group, MHDF uses mean distance:

MHDF(p,q)=max{1|Tp|aTpD(a,Tq),1|Tq|bTqD(b,Tq)} \mathrm{MHDF}(p, q) = \max \left\{ \frac{1}{|T_p|} \sum_{a \in T_p} D(a, T_q), \frac{1}{|T_q|} \sum_{b \in T_q} D(b, T_q) \right\}

This final similarity is:

GroupSim(p,q)=1MHDF(p,q) \mathrm{GroupSim}(p, q) = 1 - \mathrm{MHDF}(p, q)

The term semantic similarity method and the IC method can be set via control argument, for example:

group_sim(dag, group1, group2, method = "GroupSim_pairwise_MHDF"
    control = list(term_sim_method = "Sim_Lin_1998", IC_method = "IC_annotation")`.

Other parameters for the term_sim_method can also be set in the control list.

Paper link: https://doi.org/10.1109/ICPR.1994.576361.

GroupSim_pairwise_VHDF

It is defined as:

VHDF(p,q)=12(1|Tp|aTpD2(a,Tq)+1|Tq|bTqD2(b,Tq)) \mathrm{VHDF}(p, q) = \frac{1}{2} \left( \sqrt{\frac{1}{|T_p|} \sum_{a \in T_p} D^2(a, T_q)} + \sqrt{\frac{1}{|T_q|} \sum_{b \in T_q} D^2(b, T_q)} \right)

This final similarity is:

GroupSim(p,q)=1VHDF(p,q) \mathrm{GroupSim}(p, q) = 1 - \mathrm{VHDF}(p, q)

The term semantic similarity method and the IC method can be set via control argument, for example:

group_sim(dag, group1, group2, method = "GroupSim_pairwise_VHDF"
    control = list(term_sim_method = "Sim_Lin_1998", IC_method = "IC_annotation")`.

Other parameters for the term_sim_method can also be set in the control list.

Paper link: https://doi.org/10.1073/pnas.0702965104.

GroupSim_pairwise_Froehlich_2007

The similarity is:

GroupSim(p,q)=exp(HDF(p,q)) \mathrm{GroupSim}(p, q) = \exp(-\mathrm{HDF}(p, q))

The term semantic similarity method and the IC method can be set via control argument, for example:

group_sim(dag, group1, group2, method = "GroupSim_pairwise_Froehlich_2007"
    control = list(term_sim_method = "Sim_Lin_1998", IC_method = "IC_annotation")`.

Other parameters for the term_sim_method can also be set in the control list.

Paper link: https://doi.org/10.1186/1471-2105-8-166.

GroupSim_pairwise_Joeng_2014

Similar to VHDF, but it directly uses the similarity:

GroupSim(p,q)=12(1|Tp|aTpS2(a,Tq)+1|Tq|bTqS2(b,Tq)) \mathrm{GroupSim}(p, q) = \frac{1}{2} \left( \sqrt{\frac{1}{|T_p|} \sum_{a \in T_p} S^2(a, T_q)} + \sqrt{\frac{1}{|T_q|} \sum_{b \in T_q} S^2(b, T_q)} \right)

The term semantic similarity method and the IC method can be set via control argument, for example:

group_sim(dag, group1, group2, method = "GroupSim_pairwise_Joeng_2014"
    control = list(term_sim_method = "Sim_Lin_1998", IC_method = "IC_annotation")`.

Other parameters for the term_sim_method can also be set in the control list.

Paper link: https://doi.org/10.1109/tcbb.2014.2343963.

Pairwise edge-based methods

GroupSim_SimALN

It is based on the average distance between every pair of terms in the two groups:

GroupSim(p,q)=exp(1|Tp|*|Tq|aTp,bTqDsp(a,b)) \mathrm{GroupSim}(p, q) = \exp\left(-\frac{1}{|T_p|*|T_q|} \sum_{a \in T_p, b \in T_q} D_\mathrm{sp}(a, b)\right)

Or use the longest distance between two terms:

GroupSim(p,q)=exp(1|Tp|*|Tq|aTp,bTqlen(a,b)) \mathrm{GroupSim}(p, q) = \exp\left(-\frac{1}{|T_p|*|T_q|} \sum_{a \in T_p, b \in T_q} \mathrm{len}(a, b)\right)

There is a parameter distance which takes value of "longest_distances_via_LCA" (the default) or "shortest_distances_via_NCA":

group_sim(dag, group1, group2, method = "GroupSim_SimALN",
    control = list(distance = "shortest_distances_via_NCA"))

Paper link: https://doi.org/10.1109/CBMS.2008.27.

Groupwise IC-based methods

This category of methods depend on the IC of terms in the two groups as well as their ancestor terms.

GroupSim_SimGIC, GroupSim_SimDIC and GroupSim_SimUIC,

Denote AA and BB as the two sets of ancestors of terms in group pp and qq respectively:

𝒜p=aTp𝒜a𝒜q=bTq𝒜b \begin{align*} \mathcal{A}_p &= \bigcup_{a \in T_p} \mathcal{A}_a \\ \mathcal{A}_q &= \bigcup_{b \in T_q} \mathcal{A}_b \\ \end{align*}

The GroupSim_SimGIC, GroupSim_SimDIC and GroupSim_SimUIC are very similar. They are based on the IC of the ancestor terms, defined as:

GroupSimSimGIC(p,q)=x𝒜p𝒜qIC(x)x𝒜p𝒜qIC(x)GroupSimSimDIC(p,q)=2*x𝒜p𝒜qIC(x)x𝒜pIC(x)+x𝒜qIC(x)GroupSimSimUIC(p,q)=2*x𝒜p𝒜qIC(x)max{x𝒜pIC(x),x𝒜qIC(x)} \begin{align*} \mathrm{GroupSim}_\mathrm{SimGIC}(p, q) &= \frac{\sum\limits_{x \in \mathcal{A}_p \cap \mathcal{A}_q} \mathrm{IC}(x)}{\sum\limits_{x \in \mathcal{A}_p \cup \mathcal{A}_q} \mathrm{IC}(x)} \\ \mathrm{GroupSim}_\mathrm{SimDIC}(p, q) &= \frac{2 * \sum\limits_{x \in \mathcal{A}_p \cap \mathcal{A}_q} \mathrm{IC}(x)}{\sum\limits_{x \in \mathcal{A}_p} \mathrm{IC}(x) + \sum\limits_{x \in \mathcal{A}_q} \mathrm{IC}(x)} \\ \mathrm{GroupSim}_\mathrm{SimUIC}(p, q) &= \frac{2 * \sum\limits_{x \in \mathcal{A}_p \cap \mathcal{A}_q} \mathrm{IC}(x)}{\max\left\{\sum\limits_{x \in \mathcal{A}_p} \mathrm{IC}(x), \sum\limits_{x \in \mathcal{A}_q} \mathrm{IC}(x) \right\}} \\ \end{align*}

IC method can be set via the control argument. By default if there is annotation associated, IC_annotation is used, or else IC_offspring is used.

group_sim(dag, group1, group2, method = "GroupSim_SimGIC",
    control = list(IC_method = ...))

GroupSim_SimUI, GroupSim_SimDB, GroupSim_SimUB and GroupSim_SimNTO

These four methods are based on the counts of ancestor terms:

GroupSimSimUI(p,q)=|𝒜p𝒜q||𝒜p𝒜q|GroupSimSimDB(p,q)=2*|𝒜p𝒜q||𝒜p|+|𝒜q|GroupSimSimUB(p,q)=|𝒜p𝒜q|max{|𝒜p|,|𝒜q|}GroupSimSimNTO(p,q)=|𝒜p𝒜q|min{|𝒜p|,|𝒜q|} \begin{align*} \mathrm{GroupSim}_\mathrm{SimUI}(p, q) &= \frac{|\mathcal{A}_p \cap \mathcal{A}_q|}{|\mathcal{A}_p \cup \mathcal{A}_q|} \\ \mathrm{GroupSim}_\mathrm{SimDB}(p, q) &= \frac{2*|\mathcal{A}_p \cap \mathcal{A}_q|}{|\mathcal{A}_p| + |\mathcal{A}_q|} \\ \mathrm{GroupSim}_\mathrm{SimUB}(p, q) &= \frac{|\mathcal{A}_p \cap \mathcal{A}_q|}{\max\{|\mathcal{A}_p|, |\mathcal{A}_q|\}} \\ \mathrm{GroupSim}_\mathrm{SimNTO}(p, q) &= \frac{|\mathcal{A}_p \cap \mathcal{A}_q|}{\min\{|\mathcal{A}_p|, |\mathcal{A}_q|\}} \end{align*}

group_sim(dag, group1, group2, method = "GroupSim_SimUI")

GroupSim_SimCOU

Let’s write 𝒜p\mathcal{A}_p and 𝒜q\mathcal{A}_q as two vectors 𝐯𝐩\mathbf{v_p} and 𝐯𝐪\mathbf{v_q}. Taking 𝐯𝐩\mathbf{v_p} as an example, it is 𝐯𝐩=(w1,...,wn)\mathbf{v_p} = (w_1, ..., w_n) where nn is the number of total terms in the DAG. The value wiw_i is assigned to the corresponding term tit_i and is defined as:

𝓌i={IC(ti)ifti𝒜p0otherwise \mathcal{w}_{i} = \left\{ \begin{array}{ll} \mathrm{IC}(t_i) & \textrm{if} t_i \in \mathcal{A}_p \\ 0 & \textrm{otherwise} \end{array} \right.

The semantic similarity is defined as the cosine similarity between the two vectors:

GroupSim(a,b)=𝐯𝐩𝐯𝐪𝐯𝐩𝐯𝐪 \mathrm{GroupSim}(a, b) = \frac{ \mathbf{v_p} \cdot \mathbf{v_q} }{\left \| \mathbf{v_p} \right \| \cdot \left \| \mathbf{v_q} \right \|}

It can also be written as:

GroupSim(a,b)=x𝒜p𝒜qIC(x)2x𝒜pIC(x)2x𝒜qIC(x)2 \mathrm{GroupSim}(a, b) = \frac{\sum\limits_{x \in \mathcal{A}_p \cap \mathcal{A}_q}\mathrm{IC}(x)^2}{\sqrt{\sum\limits_{x \in \mathcal{A}_p}\mathrm{IC}(x)^2} \cdot \sqrt{\sum\limits_{x \in \mathcal{A}_q}\mathrm{IC}(x)^2}}

IC method can be set via the control argument. By default if there is annotation associated, IC_annotation is used, or else IC_offspring is used.

group_sim(dag, group1, group2, method = "GroupSim_SimCOU",
    control = list(IC_method = ...))

GroupSim_SimCOT

The semantic similarity is defined as:

GroupSim(a,b)=𝐯𝐩𝐯𝐪𝐯𝐩2+𝐯𝐪2𝐯𝐩𝐯𝐪=x𝒜p𝒜qIC(x)2x𝒜p𝒜qIC(x)2 \begin{align*} \mathrm{GroupSim}(a, b) &= \frac{ \mathbf{v_p} \cdot \mathbf{v_q} }{\left \| \mathbf{v_p} \right \|^2 + \left \| \mathbf{v_q} \right \|^2 - \mathbf{v_p} \cdot \mathbf{v_q}} \\ &= \frac{\sum\limits_{x \in \mathcal{A}_p \cap \mathcal{A}_q}\mathrm{IC}(x)^2}{\sum\limits_{x \in \mathcal{A}_p \cup \mathcal{A}_q}\mathrm{IC}(x)^2} \end{align*}

IC method can be set via the control argument. By default if there is annotation associated, IC_annotation is used, or else IC_offspring is used.

group_sim(dag, group1, group2, method = "GroupSim_SimCOT",
    control = list(IC_method = ...))

Groupwise edge-based methods

GroupSim_SimLP

It is the largest depth of terms in 𝒜p𝒜q\mathcal{A}_p \cap \mathcal{A}_q.

GroupSim(p,q)=max{δ(t):t𝒜p𝒜q} \mathrm{GroupSim}(p, q) = \max\{\delta(t): t \in \mathcal{A}_p \cap \mathcal{A}_q\}

group_sim(dag, group1, group2, method = "GroupSim_SimLP")

Link: https://bioconductor.org/packages/release/bioc/vignettes/GOstats/inst/doc/GOvis.html#go-induced-distances.

GroupSim_Ye_2005

It is a normalized version of GroupSim_SimLP:

GroupSim(p,q)=max{δ(t)δminδmaxδmin:t𝒜p𝒜q}=max{δ(t)δmax:t𝒜p𝒜q} \begin{align*} \mathrm{GroupSim}(p, q) &= \max\left\{\frac{\delta(t) - \delta_\mathrm{min}}{\delta_\mathrm{max} - \delta_\mathrm{min}}: t \in \mathcal{A}_p \cap \mathcal{A}_q\right\} \\ &= \max\left\{\frac{\delta(t) }{\delta_\mathrm{max}}: t \in \mathcal{A}_p \cap \mathcal{A}_q\right\} \end{align*}

Since the minimal depth is zero for root.

group_sim(dag, group1, group2, method = "GroupSim_Ye_2005")

Paper link: https://doi.org/10.1038/msb4100034.

Annotated items-based methods

This category of methods consider the items annotated to the two groups of terms.

GroupSim_SimCHO

It is based on the annotated items. Denote σ(t)\sigma(t) as the total number of annotated items of tt (after merging all its offspring terms). The similarity is calculated as:

GroupSim(p,q)=log(Cpq)log(Cmin/Cmax) \mathrm{GroupSim}(p, q) = \frac{\log(C_{pq})}{\log(C_\mathrm{min}/C_\mathrm{max})}

where Cpq=min{σ(t):tTpTq}C_{pq} = \min\{\sigma(t): t \in T_p \cap T_q \}, CminC_\mathrm{min} is the minimal number of annotated items in the DAG which in most cases is 1, CmaxC_\mathrm{max} is the maximal number of annotated items, which is the total number of items annotated to the complete DAG.

The similarity can also be written in form of ICanno\mathrm{IC}_\mathrm{anno}:

GroupSim(p,q)=maxxTpTqIC(x)ICmax \mathrm{GroupSim}(p, q) = \frac{\max\limits_{x \in T_p \cup T_q}\mathrm{IC}(x)}{\mathrm{IC}_\mathrm{max}}

group_sim(dag, group1, group2, method = "GroupSim_SimCHO")

GroupSim_SimALD

The similarity is calculated as:

GroupSim(p,q)=max{1σ(x)Cmax:xTpTq} \mathrm{GroupSim}(p, q) = \max\left\{ 1 - \frac{\sigma(x)}{C_\mathrm{max}}: x \in T_p \cap T_q \right\}

group_sim(dag, group1, group2, method = "GroupSim_SimALD")

Set-based methods

Since TpT_p and TqT_q are two sets, the Kappa coeffcient, Jaccard coeffcient, Dice coeffcient and overlap coeffcient can be naturally used.

group_sim(dag, group1, group2, method = "GroupSim_Jaccard", 
    control = list(universe = ...))
group_sim(dag, group1, group2, method = "GroupSim_Dice", 
    control = list(universe = ...))
group_sim(dag, group1, group2, method = "GroupSim_Overlap", 
    control = list(universe = ...))
group_sim(dag, group1, group2, method = "GroupSim_Kappa", 
    control = list(universe = ...))

Session info

## R version 4.4.1 (2024-06-14)
## Platform: x86_64-apple-darwin20
## Running under: macOS Sonoma 14.5
## 
## Matrix products: default
## BLAS:   /Library/Frameworks/R.framework/Versions/4.4-x86_64/Resources/lib/libRblas.0.dylib 
## LAPACK: /Library/Frameworks/R.framework/Versions/4.4-x86_64/Resources/lib/libRlapack.dylib;  LAPACK version 3.12.0
## 
## locale:
## [1] en_GB.UTF-8/en_GB.UTF-8/en_GB.UTF-8/C/en_GB.UTF-8/en_GB.UTF-8
## 
## time zone: Europe/Berlin
## tzcode source: internal
## 
## attached base packages:
## [1] stats     graphics  grDevices utils     datasets  methods   base     
## 
## other attached packages:
## [1] simona_1.3.12 knitr_1.48   
## 
## loaded via a namespace (and not attached):
##  [1] sass_0.4.9            xml2_1.3.6            shape_1.4.6.1        
##  [4] digest_0.6.37         magrittr_2.0.3        evaluate_0.24.0      
##  [7] grid_4.4.1            RColorBrewer_1.1-3    iterators_1.0.14     
## [10] circlize_0.4.16       fastmap_1.2.0         foreach_1.5.2        
## [13] doParallel_1.0.17     jsonlite_1.8.8        GlobalOptions_0.1.2  
## [16] promises_1.3.0        ComplexHeatmap_2.20.0 codetools_0.2-20     
## [19] textshaping_0.4.0     jquerylib_0.1.4       cli_3.6.3            
## [22] shiny_1.9.1           rlang_1.1.4           crayon_1.5.3         
## [25] scatterplot3d_0.3-44  cachem_1.1.0          yaml_2.3.10          
## [28] tools_4.4.1           parallel_4.4.1        colorspace_2.1-1     
## [31] httpuv_1.6.15         GetoptLong_1.0.5      BiocGenerics_0.50.0  
## [34] mime_0.12             R6_2.5.1              png_0.1-8            
## [37] matrixStats_1.3.0     stats4_4.4.1          lifecycle_1.0.4      
## [40] S4Vectors_0.42.1      fs_1.6.4              htmlwidgets_1.6.4    
## [43] IRanges_2.38.1        clue_0.3-65           ragg_1.3.2           
## [46] cluster_2.1.6         pkgconfig_2.0.3       desc_1.4.3           
## [49] pkgdown_2.1.0         bslib_0.8.0           later_1.3.2          
## [52] Rcpp_1.0.13           systemfonts_1.1.0     xfun_0.47            
## [55] xtable_1.8-4          rjson_0.2.22          htmltools_0.5.8.1    
## [58] igraph_2.0.3          rmarkdown_2.28        Polychrome_1.5.1     
## [61] compiler_4.4.1