
Accessing component databases of VEuPathDB
Rohit Satyam
King Abdullah University of Science & Technology, Saudi ArabiaAlberto Maillo
King Abdullah University of Science & Technology, Saudi ArabiaDavid Gomez-Cabrero
King Abdullah University of Science & Technology, Saudi ArabiaArnab Pain
King Abdullah University of Science & Technology, Saudi Arabia13 June, 2025
Gene_ID_Conversion.Rmd
Abstract
This section covers the capability of plasmoRUtils
to access various VEuPathDb database and carry out ID conversion
tasks, fetch preconfigured data tables and few other tasks.
Introduction
The plasmoRUtils package streamlines access to VEuPathDB’s family of 12 specialized databases such as ToxoDB, PlasmoDB, and PiroplasmaDB (Amos et al. 2022). It provides direct data retrieval capabilities through VEuPathDB’s RESTful API, enabling seamless integration of biological data into R workflows. The package supports downloading both standard and customized data tables, making it particularly valuable for researchers needing to combine data from multiple sources for downstream analysis.
# Load package and some other useful packages by using
suppressPackageStartupMessages(
suppressWarnings({
library(plasmoRUtils)
library(dplyr)
library(plyr)}))
Gene ID Conversion
A common challenge in bioinformatics involves mapping between
different identifier systems across databases. For apicomplexan
research, this might include converting between UniProt IDs, legacy gene
identifiers, and current Ensembl gene IDs. The toGeneid()
function addresses this need by retrieving up-to-date annotations from
VEuPathDB databases, supporting bidirectional conversion between various
ID types through flexible parameter specification.
Retrieving Gene Annotations and Alternative IDs
The toGeneid()
function enables annotation retrieval
when provided with Ensembl gene IDs. By default, it returns essential
information including gene names, symbols, outdated gene IDs, and
protein Uniprot IDs. The function’s versatility extends to supporting
custom field requests through its customFields
parameter,
with available options documented in the getTable()
help
section.
## Get annotations for list of geneIDs for PF3D7
toGeneid(c("PF3D7_0420300", "PF3D7_0621000"), from="ensembl")
#> # A tibble: 2 × 10
#> `Gene ID` `Product Description` `Gene Strand` `Gene Name or Symbol`
#> <chr> <chr> <chr> <chr>
#> 1 PF3D7_0420300 AP2 domain transcription fa… forward ApiAP2
#> 2 PF3D7_0621000 RNA polymerase subunit sigm… forward ApSigma
#> # ℹ 6 more variables: `Previous ID(s)` <chr>, `Entrez Gene ID` <chr>,
#> # `UniProt ID(s)` <chr>, `Protein Length` <chr>, `# TM Domains` <chr>,
#> # `SignalP Peptide` <chr>
## Get annotations for list of geneIDs for organisms other than PF3D7
toGeneid(inputid = c("TGME49_304740","TGME49_208030"),from="ensembl",org="Toxoplasma gondii ME49", db="toxodb")
#> # A tibble: 2 × 10
#> `Gene ID` `Product Description` `Gene Strand` `Gene Name or Symbol`
#> <chr> <chr> <chr> <chr>
#> 1 TGME49_208030 microneme protein MIC4 forward MIC4
#> 2 TGME49_304740 rhoptry kinase family prote… reverse ROP35
#> # ℹ 6 more variables: `Previous ID(s)` <chr>, `Entrez Gene ID` <chr>,
#> # `UniProt ID(s)` <chr>, `Protein Length` <chr>, `# TM Domains` <chr>,
#> # `SignalP Peptide` <chr>
## Convert uniprot IDs back to gene IDs. It will also provide Product description and Gene Symbol
toGeneid(inputid = c("Q8I1N6","C6KT48"),from="uniprot",to="ensembl" )
#> # A tibble: 2 × 4
#> `Gene ID` `Product Description` `Gene Name or Symbol` `UniProt ID(s)`
#> <chr> <chr> <chr> <chr>
#> 1 PF3D7_0420300 AP2 domain transcription … ApiAP2 Q8I1N6
#> 2 PF3D7_0621000 RNA polymerase subunit si… ApSigma C6KT48
## Using customFields to get only columns of interest
toGeneid(inputid = c("TGME49_304740","TGME49_208030"),
from="ensembl",org="Toxoplasma gondii ME49",
db="toxodb",
customFields=c("primary_key","predicted_go_component","annotated_go_function"))
#> # A tibble: 2 × 3
#> `Gene ID` `Computed GO Components` `Curated GO Functions`
#> <chr> <chr> <chr>
#> 1 TGME49_208030 extracellular region N/A
#> 2 TGME49_304740 N/A N/A
Note: Successful ID conversion requires precise organism nomenclature matching VEuPathDB’s conventions. For example,
Toxoplasma gondii ME49
must include proper spacing and special characters. Invalid query would beToxoplasma gondiiME49
or short formsTgME49
.
This functionality proves particularly valuable for enhancing differential expression analysis results with comprehensive annotations, enabling complete workflow automation without leaving command-line interfaces on HPC systems.
Accessing Preconfigured Data Tables from VEuPathDB’s component sites
Since VEuPathDB API documentation specifically encourages to use
specific organism database as quoted below, we developed
getTable()
function to fetch fields of interests from each
database separately.
There are 12 component sites and one portal: VEuPathDB.org. The component sites are: AmoebaDB, CryptoDB, FungiDB, GiardiaDB, HostDB, MicrosporidiaDB, PiroplasmaDB, PlasmoDB, ToxoDB, TrichDB, TriTrypDB and VectorBase. For most record types (all but dataset and organism), when running a search, the portal reaches out to component sites to get the search results. That means it will be faster to use a component site directly when you can.
Some frequently required fields have been provided in the help
section of getTable()
.
## To fetch table for all the genes present in an organism
getTable(org="Plasmodium falciparum 3D7", db="plasmodb") %>% head()
#> # A tibble: 6 × 10
#> `Gene ID` `Product Description` `Gene Strand` `Gene Name or Symbol`
#> <chr> <chr> <chr> <chr>
#> 1 PF3D7_0100100 erythrocyte membrane protei… forward VAR
#> 2 PF3D7_0100200 rifin reverse RIF
#> 3 PF3D7_0100300 erythrocyte membrane protei… reverse VAR
#> 4 PF3D7_0100400 rifin forward RIF
#> 5 PF3D7_0100500 erythrocyte membrane protei… reverse N/A
#> 6 PF3D7_0100600 rifin reverse RIF
#> # ℹ 6 more variables: `Previous ID(s)` <chr>, `Entrez Gene ID` <chr>,
#> # `UniProt ID(s)` <chr>, `Protein Length` <chr>, `# TM Domains` <chr>,
#> # `SignalP Peptide` <chr>
## User can also provide custom fields. For example we wish to download the P. falciparum 3D7 Proteome and phosphoproteome data during intraerythrocytic development (Quantitative) (Pease et al.)
getTable(org="Plasmodium falciparum 3D7", db="plasmodb", customFields = c("primary_key","pan_6365","pan_6366","pan_6367")) %>% head()
#> # A tibble: 6 × 4
#> `Gene ID` Ring Ave (Global pro…¹ Troph Ave (Global pr…² Schizont Ave (Global…³
#> <chr> <chr> <chr> <chr>
#> 1 PF3D7_01… N/A N/A N/A
#> 2 PF3D7_01… N/A N/A N/A
#> 3 PF3D7_01… N/A N/A N/A
#> 4 PF3D7_01… N/A N/A N/A
#> 5 PF3D7_01… N/A N/A N/A
#> 6 PF3D7_01… N/A N/A N/A
#> # ℹ abbreviated names: ¹`Ring Ave (Global proteome and phosphoproteome)`,
#> # ²`Troph Ave (Global proteome and phosphoproteome)`,
#> # ³`Schizont Ave (Global proteome and phosphoproteome)`
For more information on fields that can be supplied to
getTable()
, use the following steps:
Go to database of your interest (Say “PlasmoDB”)
Click on
Annotation, curation and identifiers
tab on your left and selectList of IDs
Scroll down and click on
Build a Web Services URL from this Search >>
hyperlink.Under the section
Choose Columns:
choose fields of your interest. Most of these fields are included in help section ofgetTable()
. However, fields specific to a particular database such as dataset related fields (starts with “pan_”) are excluded.Once you select fields, they are updated in POST section of the webpage query builder.
The constituent databases of VEuPathDB also provide some
preconfigured tables which can not be fetched via
getTable()
function. To enable users to fetch such tables,
we wrote another function called getPreconfiguredTable()
the usage of which has been shown below.
## Fetch pathway table for all the genes from MPMP database
getPreconfiguredTable(org = "Plasmodium falciparum 3D7",db = "plasmodb",customField = "MetabolicPathwaysMPMP") %>% head()
#> # A tibble: 6 × 4
#> `Gene ID` pathway_id Pathway Activity
#> <chr> <chr> <chr> <chr>
#> 1 PF3D7_0100100 lys_met Peptides with confirmed methylated lysine r… erythro…
#> 2 PF3D7_0100100 virulence Candidate genes related to virulence erythro…
#> 3 PF3D7_0100100 Par_RBC Protein-Protein Interactions between Human … erythro…
#> 4 PF3D7_0100100 PfEMP1 PfEMP1 domain architectures erythro…
#> 5 PF3D7_0100100 gene_lumef Gene expression affected by lumefantrine erythro…
#> 6 PF3D7_0100100 PQS P. falciparum genes harboring G-quadruplexes erythro…
Please note that the MPMP pathway version provided by PlasmoDB is outdated
(03-2019
). Some pathways have been revised or removed in
its entirety. If you wish to access the latest MPMP version, you can use
data("mpmp.28Aug2024")
for you analysis which was scraped
by us. If you wish to use this geneset for MPMP pathway enrichment
analysis using pathfindR,
you can do so by using data("pathfindrMPMP")
.
Similarly, predictions like TMHMM
and
SignalP
and InterPro
have not been updated
given the recent funding crunch and should be therefore used with
caution. We will discuss more about it in a separate tutorial.
Note: We urge the users to cite the original articles of the related datasets alongside plasmoRUtils.
Fetching genome metadata and strain names
In above examples, we saw the importance of passing exact name to the
org
argument for toGeneid()
to function
properly. A helper function is provided to achieve this called
listVeupathdb()
. By default, 11 columns are returned
including organism name, the respective database present in “VEuPathDB
Project” column and some more additional information. However, you can
limit the search to columns of your interests as shown below. As stated
by VEuPathDB:
The best use of the VEuPathDB portal is to get a table with all organisms in our sites, and for each organism: the component site, and the urls to access their fasta and gff files.
listVeupathdb() %>% head()
#> # A tibble: 6 × 11
#> Organism Species Genome Fasta Downloa…¹ CDS Fasta Download L…²
#> <chr> <chr> <chr> <chr>
#> 1 Edhazardia aedis USNM 4… Edhaza… http://MicrosporidiaD… http://MicrosporidiaD…
#> 2 Kluyveromyces marxianus… Kluyve… http://FungiDB.org/co… http://FungiDB.org/co…
#> 3 Aspergillus luchuensis … Asperg… http://FungiDB.org/co… http://FungiDB.org/co…
#> 4 Epichloe glyceriae E277 Epichl… http://FungiDB.org/co… http://FungiDB.org/co…
#> 5 Aspergillus versicolor … Asperg… http://FungiDB.org/co… http://FungiDB.org/co…
#> 6 Aspergillus sydowii CBS… Asperg… http://FungiDB.org/co… http://FungiDB.org/co…
#> # ℹ abbreviated names: ¹`Genome Fasta Download Link`,
#> # ²`CDS Fasta Download Link`
#> # ℹ 7 more variables: `Transcript Fasta Download Link` <chr>,
#> # `Protein Fasta Download Link` <chr>, `VEuPathDB Project` <chr>,
#> # Genes <chr>, `GFF Download Link` <chr>, `Genome Source` <chr>,
#> # `Structural Annotation Source` <chr>
listVeupathdb(customFields=c("species", "project_id")) %>% head()
#> # A tibble: 6 × 2
#> Species `VEuPathDB Project`
#> <chr> <chr>
#> 1 Entamoeba nuttalli AmoebaDB
#> 2 Acanthamoeba castellanii AmoebaDB
#> 3 Mastigamoeba balamuthi AmoebaDB
#> 4 Acanthamoeba sp. AmoebaDB
#> 5 Acanthamoeba sp. AmoebaDB
#> 6 Acanthamoeba sp. AmoebaDB
Since this function also provide URLs of FASTA and GFF files, you can use it to find the URLs of the files you are interested in and import them in R directly without leaving the console.
listVeupathdb() %>%
subset(.,Organism =="Edhazardia aedis USNM 41457") %>%
select(`GFF Download Link`) %>% as.character() %>%
rtracklayer::import.gff3() %>% head()
#> GRanges object with 6 ranges and 12 metadata columns:
#> seqnames ranges strand | source type score
#> <Rle> <IRanges> <Rle> | <factor> <factor> <numeric>
#> [1] AFBI030000... 5287-6076 + | VEuPathDB protein_coding_gene NA
#> [2] AFBI030000... 5287-6076 + | VEuPathDB mRNA NA
#> [3] AFBI030000... 5287-6076 + | VEuPathDB exon NA
#> [4] AFBI030000... 5447-6043 + | VEuPathDB CDS NA
#> [5] AFBI030000... 5287-5446 + | VEuPathDB five_prime_UTR NA
#> [6] AFBI030000... 6044-6076 + | VEuPathDB three_prime_UTR NA
#> phase ID description ebi_biotype Parent
#> <integer> <character> <character> <character> <CharacterList>
#> [1] <NA> EDEG_00001 hypothetic... protein_co...
#> [2] <NA> EDEG_00001... hypothetic... <NA> EDEG_00001
#> [3] <NA> exon_EDEG_... <NA> <NA> EDEG_00001...
#> [4] 0 EDEG_00001... <NA> <NA> EDEG_00001...
#> [5] <NA> utr_EDEG_0... <NA> <NA> EDEG_00001...
#> [6] <NA> utr_EDEG_0... <NA> <NA> EDEG_00001...
#> gene_ebi_biotype gene_id protein_source_id Note
#> <character> <character> <character> <CharacterList>
#> [1] <NA> <NA> <NA>
#> [2] protein_co... <NA> <NA>
#> [3] <NA> EDEG_00001 <NA>
#> [4] <NA> EDEG_00001 EDEG_00001...
#> [5] <NA> <NA> <NA>
#> [6] <NA> <NA> <NA>
#> -------
#> seqinfo: 342 sequences from an unspecified genome; no seqlengths
Mapping PDB IDs to Gene IDs
There is currently no facility in VEuPathDB to convert the PDB IDs to
respective gene IDs. If PDB ID corresponds to a multimer complex and if
you have multiple such PDB ids, it becomes arduous to map them to Gene
IDs manually. To provide solution to this issue, we can first convert
the PDB chains to Uniprot IDs using our own pdb2uniprot()
function and then can use toGeneid()
function to obtain
gene IDs.
pdbids <- c("7D2W","4U5A","6E10")
df <- lapply(pdbids, pdb2uniprot) %>% plyr::ldply()
geneids <- toGeneid(inputid = unique(df$attribute),from = "uniprot",to = "ensembl")
## Combining the geneIDs with df
S4Vectors::merge(df,geneids,all=TRUE, by.x="attribute", by.y="UniProt ID(s)") %>% head()
#> attribute entity_id chain_id struct_asym_id unp_start unp_end
#> 1 A0A143ZZR8 1 A A 27 206
#> 2 A0A143ZZR8 1 B B 27 206
#> 3 Q75UY1 1 C C 42 241
#> 4 Q75UY1 1 D D 42 241
#> 5 Q75UY1 1 A A 42 241
#> 6 Q75UY1 1 B B 42 241
#> start.residue_number start.author_residue_number start.author_insertion_code
#> 1 3 NA
#> 2 3 NA
#> 3 2 NA
#> 4 2 NA
#> 5 2 NA
#> 6 2 NA
#> end.residue_number end.author_residue_number end.author_insertion_code query
#> 1 182 181 7D2W
#> 2 182 181 7D2W
#> 3 201 NA 4U5A
#> 4 201 NA 4U5A
#> 5 201 NA 4U5A
#> 6 201 NA 4U5A
#> Gene ID Product Description Gene Name or Symbol
#> 1 PF3D7_1372300 Plasmodium exported protein (PHISTa-like) N/A
#> 2 PF3D7_1372300 Plasmodium exported protein (PHISTa-like) N/A
#> 3 <NA> <NA> <NA>
#> 4 <NA> <NA> <NA>
#> 5 <NA> <NA> <NA>
#> 6 <NA> <NA> <NA>
Session
sessionInfo()
#> R version 4.4.1 (2024-06-14 ucrt)
#> Platform: x86_64-w64-mingw32/x64
#> Running under: Windows 11 x64 (build 26100)
#>
#> Matrix products: default
#>
#>
#> locale:
#> [1] LC_COLLATE=English_India.utf8 LC_CTYPE=English_India.utf8
#> [3] LC_MONETARY=English_India.utf8 LC_NUMERIC=C
#> [5] LC_TIME=English_India.utf8
#>
#> time zone: Asia/Riyadh
#> tzcode source: internal
#>
#> attached base packages:
#> [1] stats graphics grDevices utils datasets methods base
#>
#> other attached packages:
#> [1] plyr_1.8.9 dplyr_1.1.4 plasmoRUtils_1.0.0
#> [4] rlang_1.1.6 readr_2.1.5 janitor_2.2.1
#> [7] randomcoloR_1.1.0.1 BiocStyle_2.32.1
#>
#> loaded via a namespace (and not attached):
#> [1] segmented_2.1-4 fs_1.6.6
#> [3] ProtGenerics_1.36.0 matrixStats_1.5.0
#> [5] bitops_1.0-9 lubridate_1.9.4
#> [7] pRoloc_1.44.1 httr_1.4.7
#> [9] RColorBrewer_1.1-3 doParallel_1.0.17
#> [11] ggsci_3.2.0 tools_4.4.1
#> [13] MSnbase_2.30.1 backports_1.5.0
#> [15] utf8_1.2.5 R6_2.6.1
#> [17] lazyeval_0.2.2 withr_3.0.2
#> [19] prettyunits_1.2.0 gridExtra_2.3
#> [21] preprocessCore_1.66.0 cli_3.6.5
#> [23] Biobase_2.64.0 textshaping_1.0.1
#> [25] gt_1.0.0 sass_0.4.10
#> [27] topGO_2.56.0 mvtnorm_1.3-3
#> [29] randomForest_4.7-1.2 proxy_0.4-27
#> [31] pkgdown_2.1.3 Rsamtools_2.20.0
#> [33] systemfonts_1.2.3 txdbmaker_1.0.1
#> [35] AnnotationForge_1.46.0 dichromat_2.0-0.1
#> [37] parallelly_1.44.0 limma_3.60.6
#> [39] rstudioapi_0.17.1 impute_1.78.0
#> [41] RSQLite_2.4.0 FNN_1.1.4.1
#> [43] generics_0.1.4 BiocIO_1.14.0
#> [45] vroom_1.6.5 gtools_3.9.5
#> [47] car_3.1-3 dendextend_1.19.0
#> [49] GO.db_3.19.1 Matrix_1.7-1
#> [51] MALDIquant_1.22.3 drawProteins_1.24.0
#> [53] S4Vectors_0.42.1 abind_1.4-8
#> [55] lifecycle_1.0.4 yaml_2.3.10
#> [57] snakecase_0.11.1 carData_3.0-5
#> [59] SummarizedExperiment_1.34.0 recipes_1.3.1
#> [61] SparseArray_1.4.8 BiocFileCache_2.12.0
#> [63] Rtsne_0.17 grid_4.4.1
#> [65] blob_1.2.4 promises_1.3.3
#> [67] crayon_1.5.3 PSMatch_1.8.0
#> [69] lattice_0.22-6 beachmat_2.20.0
#> [71] echarts4r_0.4.5.9000 annotate_1.82.0
#> [73] GenomicFeatures_1.56.0 chromote_0.5.1
#> [75] mzR_2.38.0 KEGGREST_1.44.1
#> [77] pillar_1.10.2 knitr_1.50
#> [79] GenomicRanges_1.56.2 rjson_0.2.23
#> [81] lpSolve_5.6.23 future.apply_1.11.3
#> [83] codetools_0.2-20 mgsub_1.7.3
#> [85] glue_1.8.0 V8_6.0.3
#> [87] pcaMethods_1.96.0 data.table_1.17.4
#> [89] MultiAssayExperiment_1.30.3 vctrs_0.6.5
#> [91] png_0.1-8 gtable_0.3.6
#> [93] kernlab_0.9-33 cachem_1.1.0
#> [95] gower_1.0.2 xfun_0.52
#> [97] prodlim_2025.04.28 S4Arrays_1.4.1
#> [99] mime_0.13 coda_0.19-4.1
#> [101] survival_3.8-3 ncdf4_1.24
#> [103] timeDate_4041.110 SingleCellExperiment_1.26.0
#> [105] iterators_1.0.14 hardhat_1.4.1
#> [107] lava_1.8.1 statmod_1.5.0
#> [109] MLInterfaces_1.84.0 ipred_0.9-15
#> [111] nlme_3.1-166 bit64_4.6.0-1
#> [113] progress_1.2.3 filelock_1.0.3
#> [115] LaplacesDemon_16.1.6 GenomeInfoDb_1.40.1
#> [117] bslib_0.9.0 affyio_1.74.0
#> [119] irlba_2.3.5.1 rpart_4.1.23
#> [121] colorspace_2.1-1 BiocGenerics_0.50.0
#> [123] DBI_1.2.3 nnet_7.3-19
#> [125] tidyselect_1.2.1 processx_3.8.6
#> [127] bit_4.6.0 compiler_4.4.1
#> [129] curl_6.2.3 rvest_1.0.4
#> [131] httr2_1.1.2 graph_1.82.0
#> [133] SparseM_1.84-2 xml2_1.3.8
#> [135] plotly_4.10.4 desc_1.4.3
#> [137] DelayedArray_0.30.1 bookdown_0.43
#> [139] rtracklayer_1.64.0 scales_1.4.0
#> [141] hexbin_1.28.5 affy_1.82.0
#> [143] rappdirs_0.3.3 stringr_1.5.1
#> [145] digest_0.6.37 mixtools_2.0.0.1
#> [147] rmarkdown_2.29 XVector_0.44.0
#> [149] htmltools_0.5.8.1 pkgconfig_2.0.3
#> [151] SingleR_2.6.0 sparseMatrixStats_1.16.0
#> [153] MatrixGenerics_1.16.0 dbplyr_2.5.0
#> [155] fastmap_1.2.0 htmlwidgets_1.6.4
#> [157] UCSC.utils_1.0.0 shiny_1.10.0
#> [159] DelayedMatrixStats_1.26.0 farver_2.1.2
#> [161] jquerylib_0.1.4 jsonlite_2.0.0
#> [163] mclust_6.1.1 BiocParallel_1.38.0
#> [165] mzID_1.42.0 ModelMetrics_1.2.2.2
#> [167] BiocSingular_1.20.0 RCurl_1.98-1.17
#> [169] magrittr_2.0.3 scuttle_1.14.0
#> [171] Formula_1.2-5 GenomeInfoDbData_1.2.12
#> [173] Rcpp_1.0.14 viridis_0.6.5
#> [175] MsCoreUtils_1.16.1 vsn_3.72.0
#> [177] pROC_1.18.5 stringi_1.8.7
#> [179] zlibbioc_1.50.0 MASS_7.3-61
#> [181] listenv_0.9.1 parallel_4.4.1
#> [183] splines_4.4.1 Biostrings_2.72.1
#> [185] hms_1.1.3 ps_1.9.1
#> [187] igraph_2.1.4 ggpubr_0.6.0
#> [189] QFeatures_1.14.2 ggsignif_0.6.4
#> [191] reshape2_1.4.4 biomaRt_2.60.1
#> [193] stats4_4.4.1 ScaledMatrix_1.12.0
#> [195] XML_3.99-0.18 evaluate_1.0.3
#> [197] BiocManager_1.30.25 tzdb_0.5.0
#> [199] foreach_1.5.2 httpuv_1.6.16
#> [201] tidyr_1.3.1 purrr_1.0.4
#> [203] future_1.49.0 clue_0.3-66
#> [205] bio3d_2.4-5 ggplot2_3.5.2
#> [207] rsvd_1.0.5 broom_1.0.8
#> [209] xtable_1.8-4 restfulr_0.0.15
#> [211] AnnotationFilter_1.28.0 easyPubMed_2.13
#> [213] e1071_1.7-16 rstatix_0.7.2
#> [215] later_1.4.2 class_7.3-22
#> [217] viridisLite_0.4.2 ragg_1.4.0
#> [219] tibble_3.2.1 websocket_1.4.4
#> [221] memoise_2.0.1 AnnotationDbi_1.66.0
#> [223] GenomicAlignments_1.40.0 IRanges_2.38.1
#> [225] cluster_2.1.8 globals_0.18.0
#> [227] timechange_0.3.0 caret_7.0-1
#> [229] sampling_2.10