Skip to contents

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 be Toxoplasma gondiiME49 or short forms TgME49.

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:

  1. Go to database of your interest (Say “PlasmoDB”)

  2. Click on Annotation, curation and identifiers tab on your left and select List of IDs

  3. Scroll down and click on Build a Web Services URL from this Search >> hyperlink.

  4. Under the section Choose Columns: choose fields of your interest. Most of these fields are included in help section of getTable(). However, fields specific to a particular database such as dataset related fields (starts with “pan_”) are excluded.

  5. 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

References

Amos, Beatrice, Cristina Aurrecoechea, Matthieu Barba, Ana Barreto, Evelina Y Basenko, Robert Belnap, Ann S Blevins, et al. 2022. “VEuPathDB: The Eukaryotic Pathogen, Vector and Host Bioinformatics Resource Center.” Nucleic Acids Research 50 (D1): D898–911.