DOI: 10.5061/dryad.5x69p8d9w
title: Data from: Predicting undetected native vascular plant diversity at a global scale

description: Vascular plants are diverse and a major component of terrestrial
 ecosystems, yet their geographic distributions remain incomplete. Here, I
 present a global database of vascular plant distributions by integrating
 species distribution models calibrated to species’ dispersal ability and
 natural habitats to predict native range maps for 201,681 vascular plant
 species into unsurveyed areas. Using these maps, I uncover unique patterns
 of native vascular plant diversity, endemism, and phylogenetic diversity
 revealing hotspots in under-documented biodiversity-rich regions. These
 hotspots, based on detailed species-level maps, show a pronounced
 latitudinal gradient, strongly supporting the theory of increasing
 diversity towards the equator. I trained random forest models to
 extrapolate diversity patterns under unbiased global sampling and identify
 overlaps with modelled estimations but unveiled cryptic hotspots that were
 not captured by modelled estimations. Only 29-36% of extrapolated plant
 hotspots are inside protected areas, leaving more than 60% outside and
 vulnerable. However, the unprotected hotspots harbor species with unique
 attributes that make them good candidates for conservation prioritization.
description: Experimental design The
 experimental design is divided into six steps (Fig. S1) that include: (i)
 Source data, the raw occurrence data used to produce the global species’
 range polygons; (ii) Data cleaning, which describes how errors and
 inconsistencies were minimized in the input data; (iii) Geographic
 modeling of polygons, describing how alpha hull polygons were generated
 from point occurrences; (iv) Species dispersal capacity and calibration
 area, describes how I incorporated species-specific dispersal ability in
 defining a calibration area for modeling the species distributions; (v)
 Background points, describes how background points were selected for the
 species distribution modeling; and (v) Species distribution modeling,
 describes the estimation of species distributions based on environmental
 conditions that correlate with known occurrences, and calibrated to
 species’ realized niche. Step1: Source
 data. The primary source of occurrence data used in this study
 was downloaded from Global Biodiversity Information Facility (GBIF;
 https://doi.org/10.15468/dl.vgvc3z, accessed 15 May 2023) by querying
 "TaxonKey = Tracheophyta" in a full data export. The GBIF
 represents the largest effort to mobilize vascular plant diversity data at
 a global extent. The resulting occurrences consisted of 402 million
 records from 11,517 published datasets. Step
 2: Data cleaning. From the occurrences, I used a two-step
 cleaning process to reconcile taxonomic names. First, I used the World
 Checklist of Vascular Plants (WCVP) (<em>1</em>) because it
 has the most comprehensive source for species names along with their
 native ranges within broad political divisions. I then matched names from
 the GBIF occurrence records to those in WCVP and retained only species
 names that exist in WCVP (<em>1</em>). This approach ensures
 data consistency and minimizes the tendency to include inaccurate or
 outdated names. Second, I manually harmonized four frequently interchanged
 family names to be consistent with the Angiosperm Phylogeny Group
 classification system. Specifically, Xanthorrhoeaceae was changed to
 Asphodelaceae, Cephalotaxaceae to Taxaceae, Leguminosae to Fabaceae and
 Compositae to Asteraceae. Subsequently, I analyzed the validity of the
 point data by running the occurrence records for each species through the
 function <em>clean_coordinates</em> in the R package
 CoordinateCleaner v.3.0.1 (<em>2</em>) and removing
 potentially erroneous records such as duplicates, records with doubtful or
 imprecise localities (e.g., points outside ±180° or in the sea), and on
 country or province centroids. The <em>clean_coordinates</em>
 function also has a built-in global database of over 10,055 geo-referenced
 biodiversity institutions, including herbaria, botanical gardens, zoos,
 museums, GBIF contributors, and university collections. By using a
 100-meter radius filter around these institutions, the function allowed me
 to remove records that likely originate in proximity to these institutions
 indicating records from cultivation. Next, I refined the resulting point
 records to accurately reflect species distributions within their known
 native ranges. This was achieved by intersecting point records with WCVP’s
 database of native distribution maps of vascular plants within country
 borders, i.e., level 3 of the biogeographic regions recognized by the
 Biodiversity Information Standards (also known as the Taxonomic Databases
 Working Group (<em>1, 3</em>). Points that intersected with
 the boundaries of WCVP were retained for further analysis. Subsequently,
 for species with ≥5 unique localities, I used spatial thinning to reduce
 the effects of uneven, or biased sampling. This involved subsampling
 records present in each cell using the R package dismo v.1.3-14
 (<em>4</em>) with function <em>gridSample</em>.
 The resulting records from spatial thinning represent the final cleaned
 point occurrence records utilized for subsequent analyses.
 Step 3: Geographic modeling of alpha hull
 polygons. Following data cleaning, I incorporated the thinned
 point records into polygon maps by modeling with alpha hulls using the R
 package rangeBuilder v.2.1 (<em>5</em>) with options
 <em>initialAlpha</em> = 2. For each species, the alpha
 hyperparameter value was incrementally increased until a range was formed
 that encompassed at least 99% of the occurrence records. All other
 parameters were left at their default settings. From these alpha hull
 polygons, I cropped them to land areas by removing areas falling in the
 sea and inland water bodies (using vector polygons from
 https://naturalearthdata.com). I additionally finetuned the alpha hull
 polygons by clipping them using the polygons of vascular plant families
 (<em>6, 7</em>) under the assumption that a species does not
 occur outside its family’s range. Step 4:
 Species dispersal capacity and calibration area. Predicting a
 species’ habitat suitability and ability to occupy a space requires
 knowledge of inherent dispersal capacity (<em>8</em>). I
 implemented a partial-dispersal model which incorporates demographic or
 life history factors into the species distribution modeling thereby
 preventing erroneous predictions in areas that might contain suitable
 habitat but not occupied by the species (<em>8, 9</em>). I
 started by fitting a phylogeographic Spherical Brownian Motion (SBM) model
 (<em>10, 11</em>), which quantifies the dispersal of a clade
 over time as a diffusion-like process based on a single diffusion
 coefficient <em>D</em>, while accounting for Earth's
 spherical geometry (<em>9, 12</em>) using the function
 <em>fit_sbm_const </em>in the R package castor v.1.7.10
 (<em>11, 13</em>). This yields a diffusivity score
 <em>D</em>, measured in km<sup>2</sup>
 Myr<sup>–1</sup> and it is based on maximum-likelihood and
 independent contrasts between sister lineages (<em>11</em>).
 The <em>fit_sbm_con </em>function calculates diffusivity using
 a rooted phylogenetic tree and geographic coordinates (latitudes and
 longitudes) for the tips of the tree. I then computed species dispersal
 rates as the expected dispersal distance traveled by a species in a year
 (in km yr<sup>–1</sup>) using the castor function
 <em>expected_SBM_distance</em>. The phylogenetic tree for the
 SBM modeling was generated using the R package V.PhyloMaker2 v.0.1.0
 (<em>14</em>) with function <em>phylo.maker</em>
 based on the expanded megaphylogeny (GBOTB.extended.TPL) of Smith and
 Brown (<em>15</em>) as a backbone and the
 <em>build.nodes.1</em> function (<em>14</em>).
 Missing taxa from the megaphylogeny were added using V.PhyloMaker2 under
 scenario "S2” that allows generation of random trees. I generated 2
 random trees and used it for downstream analysis. The tip coordinates were
 estimated by taking the centroid of each species’ alpha hull polygon.
 Finally, I used the species dispersal rates to define a calibration area
 by buffering the dispersal rates around the alpha hull polygons of each
 species before intersecting the buffered zones with the terrestrial
 ecoregions occupied by the species (<em>16, 17</em>) to
 predict habitat suitability of each species. The SBM-derived dispersal
 rates (km/yr) allow me to account for dispersal ability and restrict
 predictions to the realized niche (suitable and accessible areas) rather
 than the fundamental niche (all potentially suitable areas).
 Step 5: Background points. Due to the inherent
 sampling biases in vascular plant occurrence records (<em>18,
 19</em>) and the lack of absence data, I generated background
 localities randomly for modeling distributions for each species but
 accounting for sampling biases in vascular plant occurrence data. This was
 achieved using a spatial kernel density estimate (<em>20</em>)
 through the R package spatialEco v.2.0-1 (<em>21</em>) to
 generate sampling density grid standardized to 0-1
 (<em>22</em>). I used the sampling density grid as estimate of
 vascular plant sampling bias in the occurrence records (Fig. S6). I then
 used the sampling bias grid to probabilistically sample 10,000 background
 points for each species within their respective calibration areas,
 representing the natural areas where vascular plants can occur.
 Step 6: Species distribution modeling.
 Previous empirical analyses (<em>18, 19</em>) have revealed
 that multidimensional gaps and biases in vascular plant occurrence records
 can influence accurate estimation of vascular plant species distributions.
 For these reasons, I modeled species distributions using alpha hull
 polygons rather than the raw biased point records, following previous
 studies (<em>23–27</em>). Here, I used alpha hull polygons in
 similar context to how IUCN expert-range polygons are generally used in
 macroecology and biogeography. IUCN expert-range polygons incorporate
 diverse information sources such as point records, species ecology and
 distribution, local inventories, atlases, and relevant literature
 (<em>26, 27</em>). However, unlike IUCN expert-range polygons,
 alpha hulls are purely data-driven. Building on the success of IUCN
 expert-range polygons in generating accurate SDMs across various taxonomic
 groups (<em>24–27</em>) and their demonstrated agreement with
 point data (<em>23</em>), I used alpha hulls for modelling
 vascular plant distributions in the absence of expert-range polygons for
 vascular plants. Following these approaches, I sampled these alpha hull
 polygons to generate 500 points per species for input into the species
 distribution model (SDM) as in previous studies
 (<em>24–27</em>) rather than using the raw and biased point
 occurrences. To achieve this, I used a regular sampling procedure that
 ensured samples were systematically aligned by transferring values
 associated with the geometries of the alpha hull polygons to a raster at a
 resolution matching the environmental data and then thinning to 500 points
 per species, as in previous studies (<em>24</em>). This
 ensures a reproducible set of points for modeling. I acknowledge that
 actual species distributions in nature are not necessarily evenly
 distributed within its range and could potentially lead to commission
 errors, previous studies suggest minimal impact for modelling species
 distributions at macroecological scales as in this study. Furthermore, as
 IUCN range polygons generally represent expected species distributions
 (<em>28</em>), and my aim is to predict plant distribution at
 a scale comparable to IUCN polygons, this approach is appropriate for the
 data and goal of predicting expected native vascular plant
 diversity. Using the occurrence data as input, I
 employed a 75% random sample for model development, while keeping the
 remaining 25% sample for model evaluation. To predict the distribution of
 species into unsampled areas, I used species distribution models to
 predict the relationship between the observed occurrences and
 environmental factors while restricted to the calibration area that was
 defined by species’ dispersal rates and intersected with terrestrial
 ecoregions. Environmental factors were obtained from the WorldClim
 database v.2.1 (<em>29</em>) which represent historical (near
 current) climate layers for a total of 19 variables at a spatial grain
 resolution of 5 arcmin (~9 km). These variables represent climate data for
 1970-2000 and were selected because they correlate with key physical and
 biological attributes of many plant species. I additionally included
 elevation from WorldClim v.2.1 (<em>29</em>), as a variable in
 the model. I then computed the variance inflation factor for the predictor
 variables to remove highly correlated variables in a stepwise manner using
 the <em>vifstep</em> function in the R package usdm v.2.1-6
 (<em>30</em>) by setting the <em>VIF</em>
 threshold = 5. From these predictor variables, I excluded pixels
 corresponding to inland water bodies (i.e., lakes) from Natural Earth
 (https://www.naturalearthdata.com/) and then masked them to each
 species' calibration area. I used maximum entropy
 (MaxEnt) to model plant species distributions because it generally
 outperforms other methods based on predictive accuracy and does not
 generate response curves that may lead to unpredictable behavior when
 applied to novel climatic conditions (<em>31</em>). MaxEnt is
 also widely used for species distribution modelling with presence-only
 data and works well for several taxonomic groups including plants. For the
 model settings, I considered important decision steps in the input data
 and parameter settings. This was achieved using the sdm function in the R
 package phyloregion v.1.0.9 (<em>32</em>) by evaluating the
 model with a range of different complexities in terms of the feature
 classes and regularization multiplier settings in MaxEnt as follows:
 <em>linear, threshold</em>, and <em>hinge
 responses</em>, and tested a set of regularization multiplier values
 (2, 5, 10, 15, 20) under a 5-folds cross-validation framework
 (<em>33</em>). To convert the continuous
 predicted probabilities into binary presence-absence maps, I used the
 equal training sensitivity (true positive rate) and specificity threshold
 (true negative rate; <em>34</em>). I also used three
 evaluation measures for quantifying the predictive accuracy of the models:
 area under the curve (AUC), True Skill Statistic (TSS), and Boyce Index.
 The AUC measures how well predictions differentiate between presence and
 non-presence sites (<em>35, 36</em>), TSS measures the ability
 of the model to discriminate between positive and negative cases while
 accounting for both correct and incorrect classifications
 (<em>37</em>), and Boyce Index assesses how well predictions
 serve as an index of the probability of presence (<em>38,
 39</em>). AUC scores vary between 0 (random model) to 1 (best
 discriminating model), whereas both TSS scores and Boyce Index vary
 between −1 and +1, with positive values indicating model predictions
 better than random and negative values corresponding to poor prediction
 (<em>40</em>). For each species, I
 generated five sets of models and reported the median to account for
 uncertainties across different model runs. The model prediction consisted
 of a range map stored in raster format at a 5 arcmin grid cell resolution.
 The suitability of the models for each species was converted to binary
 presences by using the 95% quantile of the suitability values extracted
 from the underlying occurrence records as presence threshold. I did not
 exclude models with lower test scores because while such values may mean
 poor prediction (<em>40</em>), I found that the performance of
 the models was generally high across the evaluation metrics with a median
 AUC of 0.91 (range 0.36 to 0.999), Boyce index = 0.87 (range -0.55 to
 0.99) and TSS = 0.42 (range -0.16 to 0.88) (Fig. S2). These high
 evaluation scores indicate that my models’ predictions are a reliable
 indicator of where the plant species are likely to be found. I provided an
 ODMAP standard report of my protocol (<em>41</em>) in the
 supplementary information (Dataset S1). I provided the model outputs in
 three formats: 1) raster files of raw suitability scores to allow users to
 explore alternative approaches to thresholding the predictions, 2) raster
 files of binary presences using a 95% quantile threshold for downstream
 analyses, and 3) range polygons derived by dissolving the binary presence
 rasters using the R package terra v.1.7-3 (<em>42</em>)
 followed by spline interpolation of vertices in the R package smoothr
 v.1.0.0 (<em>43</em>) to smooth the jagged edges and sharp
 corners of the polygon maps. The smoothing of the polygon boundaries was
 done to achieve a more visually appealing representation for these
 standardized range maps consistent with IUCN’s mapping standard
 (<em>44</em>). However, it is important to emphasize that this
 smoothing does not alter the underlying binary presence data, which
 maintains all the details from the model, in downstream analyses (e.g.,
 see Fig. S8). All calculations were processed using R v.4.2.0
 "Vigorous Calisthenics" (<em>45</em>) on the
 Sherlock computing clusters of Stanford University.
 Data repository. The database contains range
 maps for 201,681 species and covers the basal Eudicots (4,770 species),
 Monocots (42,370 species), Superasteridae (77,443 species), Superrosidae
 (63,372 species), Magnoliidae (7,243 species), Chloranthales (204
 species), Gymnosperms (854 species), and Ferns (5,425 species). The
 geographic coverage of the data is global spanning all continents except
 for Antarctica and focused on only the terrestrial part of a species’
 native range, and representing contemporary times. Each species’ spatial
 data is projected in the WGS84 coordinate reference system (EPSG:4326
 code) and stored as vector polygons in GEOPACKAGE format for maximum
 compression, allowing easy integration with other spatial datasets and
 data manipulation using a standard computer. I additionally stored each
 species’ binary spatial raster as well as the raw suitability scores as
 raster files to explore alternative approaches of thresholding the
 predictions. The dataset can be downloaded from Dryad repository (https://doi.org/10.5061/dryad.5x69p8d9w). Future versions of the data will be made available on GreenMaps (greenmaps.net). Technical Validation.<strong> </strong>I explored the validity of the modeled maps using three different distribution datasets: 1) traditional SDM of using raw input occurrence records, 2) extent-of-occurrence range polygons from the IUCN, and 3) range maps for trees of southern Africa. Traditional SDM of directly using raw input occurrence points. I tested the sensitivity of my results by comparing my polygon approach to a traditional species distribution model of directly using the raw and biased point records. This was achieved using range filling analysis (<em>46</em>) by calculating the proportion of each traditional SDM occupied in the realized ranges of the polygon SDM approach for 97,026 species with sufficient raw point data. I found that while models based on raw point records tend to underestimate species range, their spatial congruence with the range polygons is generally high (median range filling = 60%). In addition, range filling is negatively correlated with geographic range size meaning that large-ranged species based on the traditional SDM approach are underfilling their potential ranges (Fig. S9). These finding indicates that the impact of my method on overall patterns is likely minimal given the high congruence of 60%. The underestimation observed in the raw point record especially for large-ranged species highlights the inherent challenges of using raw point records to infer complete species distributions. Extent-of-occurrence maps from IUCN. I contrasted my modelled estimations with extent-of-occurrence polygon maps from the International Union for the Conservation of Nature (IUCN) for 9,081 plant species. The IUCN maps were derived as hand-drawn polygons based on expert knowledge of species ecology and distribution, local inventories, atlases, and relevant literature (<em>26, 27</em>). I contrasted the boundaries of my modeled maps with the range maps of 9,081 plant species from the IUCN spatial database. To achieved this, I projected each map in WGS84 coordinate reference system (EPSG:4326 code) and assessed congruence between the maps using the V-measure statistic in the R package sabre v.0.4.3 (<em>47</em>). The V-measure was initially developed for comparing spatial congruence between regionalizations (<em>47</em>) and can be used to compare the boundaries of two categorical maps. Spatial association was computed using the function <em>vmeasure_calc</em> and setting the option B &gt; 1 so that completeness is weighted more than homogeneity as my goal was to generate reliable range maps for macroecology and biogeography and not necessarily expecting a one-to-one match between my modeled maps and the hand-drawn maps from expert knowledge. V-measure scores range from 0 (incongruence) to 1 (indicating perfect similarity). I found relatively high degree of similarity between my modeled maps and range polygons from the IUCN with most of the pairwise comparisons, represented by the 25<sup>th</sup> and 75<sup>th</sup> percentiles of the V-measure scores, falling within the range of 0.03 and 0.6. This level of congruence lends support to the accuracy and reliability of my modeled maps, affirming their ability to faithfully represent the expected plant distributions as predicted by expert knowledge. Range maps for trees of southern Africa. I also compared my modeled maps with range polygons for trees of southern Africa. The dataset came from a study that mapped tree diversity hotspots in southern Africa (<em>48</em>) using polygon range maps that represent the geographic area encompassing the known extent of over 1400 tree species in southern Africa. The maps were derived from Coates Palgrave (<em>49</em>) and Van Wyk et al. (<em>50</em>) as hand-drawn polygons based on expert knowledge of the flora of southern Africa. I also contrasted the boundaries of my modeled maps with the range maps of 886 tree species of southern Africa from Coates Palgrave (<em>49</em>) and Van Wyk et al. (<em>50</em>) using the V-measure statistic. I found some degree of similarity between my modelled estimations and that of trees of southern Africa. In fact, most of the pairwise comparisons, represented by the 25<sup>th</sup> and 75<sup>th</sup> percentiles of the V-measure scores, fell within the range of 0.07 to 0.36, indicating ability of representing expected plant distributions based on expert knowledge. Quality control and data limitation. A general shortcoming of this study is the assumption of even distribution of occurrences across a species range. In reality, species distributions are not evenly distributed – some regions have more occurrences than others. This assumption of evenness may not fully capture fine-scale community processes. While I anticipate that analyses of greater sampling coverage within a plant range might refine my results and capture finer-scale community processes such as species interactions, currently, the data are not available to explore this further. The hidden hotspots uncovered in this study may eventually guide future plant collecting to realize the heterogenous distributions expected in nature. Users are encouraged to be aware of these limitations and carefully check for possible effects in their own analyses and results. Despite these limitations, the strong alignment of my modelled estimates with three different distribution datasets (traditional SDM of using raw input occurrence records, extent-of-occurrence range polygons from the IUCN, and range maps for trees of southern Africa), suggests potential predictions into unsampled regions to overcome inherent biases in vascular plant sampling. Another limitation is that I generated maps for a substantial number of plant species, 201,681 species, which accounts for ~60% of the total described plant species (<em>1, 51, 52</em>). This means that geographic information for approximately 40% of vascular plant species is still lacking. As this data is continuously updated on DRYAD, my intention is to enhance and revise the maps as new distribution data becomes available. However, these maps still hold significant value for addressing questions at regional to global scales where macroecological and biogeographical processes such as niche conservatism, speciation and dispersal operate. As additional data is acquired, I intend to refine the maps to make them more relevant to local and community-scale research questions. Data analysis The data analysis is divided into four steps. Calculations of diversity patterns. To create consistent spatial representation, I projected each plant species modeled range polygon to a Wagner IV projection and subsequently rasterized using a masked raster grid of 20 km ´ 20 km at Wagner IV projection. To contrast my modelled estimates to those of terrestrial tetrapods (amphibians, birds, reptiles and mammals), I downloaded expert range maps for each tetrapod species from the International Union of Conservation of Nature’s (IUCN) Red List of Threatened Species spatial database (<em>44</em>) for amphibians, reptiles and mammals, and BirdLife International (<em>53</em>) for birds. These range polygons were also rasterized at 20 km resolution and reprojected to the Wagner IV projection. All species rasters were converted to a community matrix using the <em>rast2comm</em> function in the R package phyloregion (<em>32</em>) resulting in a sparse community matrix of 201,681 species ´ 350,488 pixels for plants, 6,738 species ´ 293,582 pixels for amphibians, 9,469 species ´ 349,579 pixels for birds, 9,617 ´ 305,244 pixels for reptiles, and 5,394 species ´ 348,959 pixels for mammals, that are essential for subsequent biodiversity calculations. To account for shared phylogenetic relationships in my diversity calculations, I paired the plant community matrix with a phylogenetic tree generated using the expanded megaphylogeny (GBOTB.extended.TPL) of Smith and Brown (<em>15</em>) as a backbone and the <em>build.nodes.1</em> function in the R package V.PhyloMaker2 (<em>14</em>). Missing taxa from the megaphylogeny were added using scenario "S2” which allows the generation of random trees. The phylogenetic tree for amphibians was obtained from ref. (<em>54</em>), birds (<em>55</em>), reptiles (<em>56</em>), and mammals (<em>57</em>). Using the community matrix, I calculated species richness as the counts of species per pixel, phylogenetic diversity as the sum of branch lengths connecting species from the tip to the root of a phylogeny (<em>58</em>) using the function <em>PD</em> in the R package phyloregion v.1.0.9 (<em>32</em>), weighted endemism as species richness weighted by range size (<em>32</em>) using <em>weighted_endemism</em> function in phyloregion, and phylogenetic endemism as range-weighted phylogenetic diversity (<em>59</em>) using <em>phylo_endemism</em> function in the R package phyloregion v.1.0.9 (<em>32</em>). I also estimated spatial clustering or dispersion of biodiversity patterns using Moran’s <em>I</em> spatial autocorrelation measure (Monte Carlo test, 999 randomizations), with values of 1 indicating clustered patterns and 0 corresponding to even diversity patterns. To identify areas of exceptional diversity, I mapped hotspots by selecting the top 10% richest pixels for each of the calculated diversity metrics using the <em>hotspots</em> function in the R package phyloregion v.1.0.9 (<em>32</em>). Validation using machine learning random forest extrapolations. I tested my modelled estimates for the effects of sampling biases common in plant occurrence records (Fig. S6) by training a random forest model to extrapolate my modelled estimates under a globally unbiased plant sampling. I used the R package ranger (<em>60</em>) to build models predicting the key diversity metrics based on 1,000 random trees, using the following variables as predictor variables: topography (<em>29</em>), summed tetrapod diversity, plant sampling bias grid (see step 5), floristic realms of the world, and each of the 19 bioclimatic variables (<em>29</em>). I utilized tetrapod diversity derived from range polygons as a variable in training the random forest machine learning models because their geographic sampling is less biased and they represent the most reliable and readily available data on the expected distributions of species at a global scale. I resampled each predictor variable to a grain resolution of 10 km and reprojected to Wagner IV projection. I then computed the variance inflation factor for the predictor variables to remove highly correlated variables in a stepwise manner using the <em>vifstep</em> function in the R package usdm version 2.1-6 (<em>30</em>) by setting the VIF threshold = 5. My estimate of sampling density (i.e., sampling biases) was calculated using a spatial kernel density estimate (<em>20</em>) through the R package spatialEco v.2.0-1 (<em>21</em>) to generate sampling density grid standardized to 0-1 (<em>22</em>). I assessed model accuracy using 5-fold spatial block cross-validation (<em>61</em>) and by systematically dividing the data into spatial blocks using predefined floristic realms of the world (<em>16</em>). Spatial cross-validation allowed the identification of optimal combination of hyperparameters that minimize the root-mean-square error (RMSE) or "mtry", by tuning from 1 to 10, resulting in mtry = 4 for all diversity metrics except phylogenetic diversity with mtry = 5 (Fig. S7). These settings were used to train the final models and to generate extrapolated values per pixel for each metric under a scenario of globally unbiased sampling. Based on this scenario of globally unbiased sampling, I set all cells in the sampling bias grid to 1 under the assumption that all places across the globe have equal plant sampling, although I recognize that biodiversity is not randomly distributed. I then contrasted overlap of the extrapolated maps from the globally unbiased sampling to the modelled estimation maps. Spatial correlations. I conducted spatially corrected correlations between plant diversity patterns and those of the different tetrapod classes for each metric. The spatial correlations were conducted using a corrected Pearson’s correlation for spatial autocorrelation using the package SpatialPack v.0.4 (<em>62</em>). Protection levels by protected areas. To determine the level of protection of each diversity hotspot, I masked the hotspots based on modelled estimates onto existing protected area networks. I used version 1.6 of the protected areas of the IUCN categories I–IV compiled from the World Database on Protected Areas (<em>63</em>). I then calculated the proportion of the area covered by each hotspot inside and outside protected areas using the functions <em>mask</em> and <em>crop</em> in the R package terra v.1.7-3 (<em>42</em>). Last, I explored the impact of protected areas on common species-level attributes including functional traits (such as plant height and seed size), evolutionary history (evolutionary distinctiveness), and extinction risk (defined as degree of threat facing a species), for each diversity hotspot inside and outside protected areas. Majority of the trait data were compiled from the Missouri Botanical Garden’s Tropicos database (https://www.tropicos.org/), supplemented by online regional databases. A detailed list of all consulted online regional databases is included in the supplementary materials (Dataset S2) while Table S2 lists the number of species with data for each trait category. It is important to note that the trait analysis focuses only on hotspots (top 10% of highest-ranking pixels for each metric), so the traits do not represent the entire global range. The extinction risk data in this study came from published dataset of machine learning predictions of conservation status for over 150,000 land plant species (<em>64</em>). The model considers factors known to influence extinction risk, including geographic range, environmental variables, and morphological traits (<em>64</em>). It estimates the likelihood of a species facing some level of threat. The predicted probability is then scaled to match the IUCN Red List categories using the expected probability of extinction over 100 years of each taxon (<em>65</em>) as follows: Least Concern = 0.001, Near Threatened & Conservation Dependent = 0.01, Vulnerable = 0.1, Endangered = 0.67, and Critically Endangered = 0.999. I used these machine-learning predictions because they provide broader coverage (for over 150,000 species) whereas formal IUCN Red List assessments are only available for a fraction of plant species (~34,000 as of March 15, 2023). This allows me to analyze a much larger dataset and identify potential coverage within protected and unprotected hotspots. Evolutionary distinctiveness was calculated with the <em>evol_distinct</em> function in the R package phyloregion v.1.0.9 (<em>32</em>) and is defined as the partitioning of phylogenetic branch lengths by the total number of species spanned by the tree weighted by the amount of unique evolutionary history they represent (<em>66</em>). In this way, evolutionary distinctiveness is considered as a species-level attribute, whereas diversity metrics such as phylogenetic endemism quantify biodiversity across geographic sites. Thus, the analysis aims to assess whether protected hotspots harbor a higher concentration of evolutionarily distinct lineages compared to unprotected hotspots. This means that within hotspots of phylogenetic endemism, for example, one can easily assess species that overlap these protected or unprotected hotspots on the basis of their evolutionary distinctiveness. To achieve this, I first identify grid cells unique to either protected or unprotected hotspots. Then, I identified all species that intersected with each type of hotspot. For each identified species within a hotspot, I extracted its relevant trait. I compared the functional traits between protected and unprotected hotspots using a <em>t</em>-test followed by Cohen's <em>d</em> with 1000 bootstrap replicates to estimate the effect size. References 1.     R. Govaerts, E. Nic Lughadha, N. Black, R. Turner, A. Paton, The World Checklist of Vascular Plants, a continuously updated resource for exploring global plant diversity. <em>Sci. Dat</em>. 8, 215 (2021). 2.     A. D. Zizka <em>et al.,</em> CoordinateCleaner: Standardized cleaning of occurrence records from biological collection databases. <em>Methods Ecol. Evol.</em> 10, 744–751 (2019). 3.     R. K. Brummitt, <em>World geographical scheme for recording plant distributions, 2nd ed </em>(Biodiversity Information Standards (TDWG) http://www.tdwg.org/ standards/109, 2001). 4.     R. J. Hijmans, S. Phillips, J. Leathwick, J. Elith, <em>dismo: Species Distribution Modeling</em> (R package version 1.3-14, &lt;https://CRAN.R-project.org/package=dismo&gt;, 2023). 5.     A. R. D. Rabosky <em>et al.</em>, Coral snakes predict the evolution of mimicry across New World snakes. <em>Nat. Commun</em>. 7, 11484 (2016). 6.     V. H. Heywood, <em>Flowering plants of the world</em> (Batsford, 1993). 7.     APG IV, An update of the Angiosperm Phylogeny Group classification for the orders and families of flowering plants. <em>Bot. J. Linn. Soc.</em> 181, 1–20 (2016). 8.     B. L. Bateman, H. T. Murphy, A. E. Reside, K. Mokany, J. VanDerWal, Appropriateness of full-, partial- and no-dispersal scenarios in climate change impact modelling. <em>Divers. Distrib</em>. 19, 1224–1234 (2013). 9.     N. V. Barve <em>et al.,</em> The crucial role of the accessible area in ecological niche modeling and species distribution modeling. <em>Ecol. Model</em>. 222, 1810–1819 (2011). 10.  D. R. Brillinger, A particle migrating randomly on a sphere. In: <em>Selected works of David Brillinger</em>. Springer. p. 73–87 (2012). 11.  S. Louca, Phylogeographic Estimation and Simulation of Global Diffusive Dispersal. <em>Syst. Biol.</em> 70, 340–359 (2021). 12.  F. Perrin, Étude mathématique du mouvement brownien de rotation. <em>In Annales scientifiques de l’École Normale Supérieure</em> 45, 1–51 (1928). 13.  S. Louca, M. Doebeli, Efficient comparative phylogenetics on large trees. <em>Bioinformatics</em> 34, 1053–1055 (2018). 14.  Y. L. Jin, H. Qian, V.PhyloMaker: An R package that can generate very large phylogenies for vascular plants. <em>Ecography</em> 42, 1353–1359 (2019). 15.  S. A. Smith, J. W. Brown, Constructing a broadly inclusive seed plant phylogeny. <em>Am. J. Bot.</em> 105, 302–314 (2018). 16.  D. M. Olson <em>et al.,</em> Terrestrial ecoregions of the world: A new map of life on Earth. <em>Bioscience</em> 51, 933–938 (2001). 17.  Folk, R. A. C. J. Visger, P. S. Soltis, D. E. Soltis, R. P. Guralnick, Geographic range dynamics drove ancient hybridization in a lineage of Angiosperms. <em>Am. Nat</em>. 192, 171–187 (2018). 18.  C. Meyer, P. Weigelt, H. Kreft, Multidimensional biases, gaps and uncertainties in global plant occurrence information. <em>Ecol. Lett</em>. 19, 992–1006 (2016). 19.  B. H. Daru, J. Rodriguez. Mass production of unvouchered records fails to represent global biodiversity patterns. <em>Nat. Ecol. Evol</em>. 7, 816-831 (2023). 20.  T. Duong, M. L. Hazelton, Cross-validation bandwidth matrices for multivariate kernel density estimation. <em>Scand. J. Stat</em>. 32, 485-506 (2005). 21.  J. S. Evans, M. A. Murphy, <em>spatialEco</em> (<em>R package version 2.0-1</em>, &lt;https://github.com/jeffreyevans/spatialEco&gt;, 2023). 22.  L. R. Hertzog, A. Besnard, P. Jay-Robert, Field validation shows bias-corrected pseudo-absence selection is the best method for predictive species-distribution modelling. <em>Divers. Distrib</em>. 20, 1403–1413 (2014). 23.  G. F. Ficetola <em>et al.</em>, An evaluation of the robustness of global amphibian range maps. <em>J. Biogeogr.</em> 41, 211–221 (2014). 24.  Y. Fourcade, Comparing species distributions modelled from occurrence data and from expert-based range maps. Implication for predicting range shifts with climate change. <em>Ecol. Inform</em>. 36, 8–14 (2016). 25.  B. H. Alhajeri, Y. Fourcade, High correlation between species-level environmental data estimates extracted from IUCN expert range maps and from GBIF occurrence data. <em>J. Biogeogr</em>. 46, 1329–1341 (2019). 26.  J. A. Velasco <em>et al.,</em> Synergistic impacts of global warming and thermohaline circulation collapse on amphibians. <em>Commun. Biol.</em> 4, 141 (2021). 27.  C. Ureta <em>et al.</em>, Evaluation of animal and plant diversity suggests Greenland’s thaw hastens the biodiversity crisis. <em>Commun. Biol.</em> 5, 985 (2022). 28.  C. Meyer, H. Kreft, R. Guralnick, W. Jetz, Global priorities for an effective information basis of biodiversity distributions. <em>Nat. Commun</em>. 6, 1–8 (2015). 29.  S. E. Fick, R. J. Hijmans. WorldClim 2: new 1km spatial resolution climate surfaces for global land areas. <em>Int. J. Climatol</em>. 37, 4302–4315 (2017). 30.  B. Naimi <em>et al.</em> Where is positional uncertainty a problem for species distribution modelling. <em>Ecography</em> 37, 191–203 (2014). 31.  S. J. Phillips, R. P. Anderson, R. E. Schapire, Maximum entropy modeling of species geographic distributions. <em>Ecol. Model</em>. 190, 231–259 (2006). 32.  B. H. Daru, P. Karunarathne, K. Schliep, phyloregion: R package for biogeographic regionalization and macroecology. <em>Methods Ecol. Evol.</em> 11, 1483–1491 (2020). 33.  A. Radosavljevic, R. P. Anderson, Making better MAXENT models ofspecies distributions: Complexity, overfitting and evaluation. <em>J. Biogeogr</em>. 41, 629–643 (2014). 34.  C. Liu, P. M. Berry, T. P. Dawson, R. G. Pearson, Selecting thresholds of occurrence in the prediction of species distributions. <em>Ecography</em> 28, 385–393 (2005). 35.  M. H. Zweig, G. Campbell, Receiver-operating characteristic (Roc) plots—a fundamental evaluation tool in clinical medicine. <em>Clin. Chem</em>. 39, 561–577 (1993). 36.  A. H. Fielding, J. F. Bell, A review of methods for the assessment of prediction errors in conservation presence/absence models. <em>Environ. Conserv</em>. 24, 38–49 (1997). 37.  O. Allouche, A. Tsoar, R. Kadmon, Assessing the accuracy of species distribution models: prevalence, kappa and the true skill statistic (TSS). <em>J. Appl. Ecol</em>. 43, 1223–1232 (2006). 38.  M. S. Boyce, P. R. Vernier, S. E. Nielsen, F. K. A. Schmiegelow, Evaluating resource selection functions. <em>Ecol. Model</em>. 157, 281–300 (2002). 39.  A. H. Hirzel, G. Le Lay, V. Helfer, C. Randin, A. Guisan, Evaluating the ability of habitat suitability models to predict species presences. <em>Ecol. Model</em>. 199, 142–152 (2006). 40.  J. R. Landis, G. G. Koch, The measurement of observer agreement for categorical data. <em>Biometrics</em> 33, 159–175 (1977). 41.  D. Zurell <em>et al.</em>, A standard protocol for reporting species distribution models. <em>Ecography</em> 43, 1261–1277 (2020). 42.  R. Hijmans, <em>terra: Spatial Data Analysis</em> (R package version 1.7-3, &lt;https://CRAN.R-project.org/package=terra&gt;, 2023). 43.  M. Strimas-Mackey, <em>smoothr: Smooth and Tidy Spatial Features</em> (R package version 1.0.0, &lt;https://CRAN.R-project.org/package=smoothr&gt;, 2023). 44.  The IUCN Red List of Threatened Species.<em>Version 2022-1</em> (IUCN, accessed 19 October 2022); https://www.iucnredlist.org 45.  R Core Team, <em>R: A language and environment for statistical computing</em> (R Foundation for Statistical Computing, Vienna, Austria. URL https://www.R-project.org/, 2022). 46.  B. J. Seliger, B. J. McGill, J. C. Svenning, J. L. Gill, Widespread underfilling of the potential ranges of North American trees. <em>J. Biogeogr</em>. 48, 359–371 (2021). 47.  J. Nowosad, T. F. Stepinski, Spatial association between regionalizations using the information-theoretical V-measure. <em>Int. J. Geogr. Inform</em>. <em>Sci</em>. 32, 12, 2386–2401 (2018) 48.  B. H. Daru, M. Van der Bank, T. J. Davies, Spatial incongruence among hotspots and complementary areas of tree diversity in southern Africa. <em>Divers. Distrib</em>. 21, 769–780 (2015). 49.  M. Coates Palgrave, <em>Keith Coates Palgrave trees of southern Africa, 3rd edn</em> (Struik, Cape Town, South Africa, 2002). 50.  B. Van Wyk, E. van der Berg, M. Coates Palgrave, M. Jordaan, <em>Dictionary of names for southern African trees, 1st edn</em> (Briza, Pretoria, South Africa, 2011). 51.  E. N. Lughadha <em>et al.</em>, Counting counts: revised estimates of numbers of accepted species of flowering plants, seed plants, vascular plants and land plants with a review of other recent estimates. <em>Phytotaxa</em> 272, 82–88 (2016). 52.  K. J. Willis, <em>The state of the world’s plants 2017 Report</em> (R. B. Gardens Ed) (Kew, UK: Royal Botanic Gardens, 2017). 53.  Bird Species Distribution Maps of the World. <em>Version 2020.1</em> (BirdLife International, 2020); http://datazone.birdlife.org/species/requestdis 54.  W. Jetz, R. A. Pyron, The interplay of past diversification and evolutionary isolation with present imperilment across the amphibian tree of life. <em>Nat. Ecol. Evol</em>. 2, 850–858 (2018). 55.  Jetz, W., Thomas, G. H., Joy, J. B., Hartmann, K. & Mooers, A. O. The global diversity of birds in space and time. <em>Nature</em> 491, 444–448 (2012). 56.  J. R. R. Tonini, K. H. Beard, R. B. Ferreira, W. Jetz, R. A. Pyron, Fully-sampled phylogenies of squamates reveal evolutionary patterns in threat status. <em>Biol. Conserv</em>. 204, 23–31 (2016). 57.  O. R. Bininda‐Emonds <em>et al</em>., The delayed rise of present‐day mammals. <em>Nature</em> 446, 507–512 (2007). 58.  D. P. Faith, Conservation evaluation and phylogenetic diversity. <em>Biol. Conserv</em>. 61, 1–10 (1992). 59.  D. Rosauer, S. W. Laffan, M. D. Crisp, S. C. Donnellan L. G. Cook, Phylogenetic endemism: a new approach for identifying geographical concentrations of evolutionary history. <em>Mol. Ecol.</em> 18, 4061–4072 (2009). 60.  M. N. Wright, A. Ziegler, ranger: A Fast Implementation of Random Forests for High Dimensional Data in C++ and R. <em>J. Stat. Softw</em>. 77, 1–17 (2017). 61.  R. Valavi, G. Guillera-Arroita, J. J. Lahoz-Monfort, J. Elith, Predictive performance of presence-only species distribution models: a benchmark study with reproducible code. <em>Ecol. Monogr</em>. 92, e01486 (2022). 62.  R. Vallejos, F. Osorio, M. Bevilacqua, <em>Spatial Relationships Between Two Georeferenced Variables: with Applications in R</em> (Springer, New York. ISBN   978-3-030-56680-7, 2020). 63.  UNEP-WCMC, <em>User Manual for the World Database on Protected Areas and world database on other effective area-based conservation measures</em> (1.6. UNEP-WCMC: Cambridge, UK. Available at: http://wcmc.io/WDPA_Manual, 2019). 64.  T. A. Pelletier <em>et al</em>., Predicting plant conservation priorities on a global scale. <em>Proc. Natl Acad. Sci. USA</em> 115, 13027–13032 (2018). 65.  D. W. Redding, A. Ø. Mooers, Incorporating evolutionary measures into conservation prioritization. <em>Conserv. Biol</em>. 20, 1670–1678 (2006). 66.  M. W. Cadotte, T. J. Davies, Rarest of the rare: Advances in combining evolutionary distinctiveness and scarcity to inform conservation at biogeographical scales. <em>Divers. Distrib</em>. 16, 376–385 (2010).
description: # Data from: Predicting undetected native vascular plant diversity at a
 global scale **Data description within the folder:** This folder includes
 the R scripts for the figures and a subset of the data for Daru (2024)
 paper. 1. `plant_presab_wag4.csv` : Community dataframe of species
 distributions in long format 2. `RESULTS`: Folder with all results **a.
 correlations** \* /wag4: Folder with results and maps projected under the
 Wagner IV projection. /spatialreg: Folder with spatial correlation tests
 for various diversity metrics \+ /endemism: Contains correlation test
 results for different test sets focused on weighted endemism
 (corr_test_set1.csv to corr_test_set28.csv). \+ /PD: Contains correlation
 test results for different test sets focused on phylogenetic diversity
 (corr_test_set1.csv to corr_test_set28.csv). \+ /PE: Contains correlation
 test results for different test sets focused on phylogenetic endemism
 (corr_test_set1.csv to corr_test_set28.csv). \+ /richness: Contains
 correlation test results for different test sets focused on species
 richness (corr_test_set1.csv to corr_test_set28.csv). **b.
 future_diversity**: Supplementary analysis that includes geographic change
 predictions for biodiversity under different climate scenarios. \*
 /geo_change: Folder with testing geographic change of diversity metrics \+
 /2021-2040: Predictions for SSP scenarios 126, 245, 370, and 585. \+
 /2041-2060: Predictions for SSP scenarios 126, 245, 370, and 585. \+
 /2061-2080: Predictions for SSP scenarios 126, 245, 370, and 585. \+
 /2081-2100: Predictions for SSP scenarios 126, 245, 370, and 585. **c.
 geodata**: Spatial files and geographic layers \* /wag4: Includes various
 geographic datasets in both raster and vector formats, such as
 wc2.1_10m_bio (climate-related data) and newRealms_knn.tif (geospatial
 data). \* /wgs84: Contains geographic projections like newRealms_knn.tif
 and various vector datasets such as ne_10m_lakes (lakes shapefiles) and
 worldmap_BEHRMANN.gpkg. **d. hotspots**: Contains analyses related to
 biodiversity hotspots. /wag4: Folder with hotspots results and maps
 projected under the Wagner IV projection. /protected_hotspots: Bar plot
 data and t-tests for metrics including endemism, PD, PE, and richness.
 /CSVs_4_barplots \* /emp: Empirical analysis \+ /plant_endemism.csv \+
 /plant_PD.csv \+ /plant_PE.csv \+ /plant_richness.csv \* /ext_rf: Random
 Forest extrapolations \+ /plant_endemism.csv \+ /plant_PD.csv \+
 /plant_PE.csv \+ /plant_richness.csv \* /t.test_emp.csv \* /t.test_ext.csv
 /emp: Empirical analysis \+ /plant_endemism.csv \+ /plant_PD.csv \+
 /plant_PE.csv \+ /plant_richness.csv /ext: Random Forest extrapolations \+
 /plant_endemism.csv \+ /plant_PD.csv \+ /plant_PE.csv \+
 /plant_richness.csv \* /YES_VERT: TIF files representing geographic
 predictions for endemism, PD, PE, and richness with vertical data. **e.
 LatDiversityGradient.csv** **f. PLANTS**: Spatial layers for plant
 biodiversity metrics \* /wag4: TIF files representing plant endemism,
 phylogenetic diversity, phylogenetic endemism, and species richness. **g.
 randomforest**: Modeling results using Random Forest methods. \* /wag4 \*
 /NO_VERT: Results without vertebrate data, including tuned model files and
 output raster files for biodiversity metrics. \* /model_results \*
 /endemism \+ /endemism_log_emp_noVERT.tif \+ /endemism_log_ext_noVERT.tif
 \* /PD \+ /PD_ext_noVERT.tif \+ /PD_noVERT.tif \* /PE \+
 /PE_log_emp_noVERT.tif \+ /PE_log_ext_noVERT.tif \* /richness \+
 /richness_ext_noVERT.tif \+ /richness_noVERT.tif \* /model_tune \*
 /endemism: Weighted endemism random forest model tune files \+ /mt1_cv.rds
 \+ /mt2_cv.rds \+ /mt3_cv.rds \+ /mt4_cv.rds \+ /mt5_cv.rds \+ /mt6_cv.rds
 \+ /mt7_cv.rds \+ /mt8_cv.rds \+ /mt9_cv.rds \+ /mt10_cv.rds \* /PD:
 Phylogenetic diversity random forest model tune files \+ /mt1_cv.rds \+
 /mt2_cv.rds \+ /mt3_cv.rds \+ /mt4_cv.rds \+ /mt5_cv.rds \+ /mt6_cv.rds \+
 /mt7_cv.rds \+ /mt8_cv.rds \+ /mt9_cv.rds \+ /mt10_cv.rds \* /PE:
 Phylogenetic endemism random forest model tune files \+ /mt1_cv.rds \+
 /mt2_cv.rds \+ /mt3_cv.rds \+ /mt4_cv.rds \+ /mt5_cv.rds \+ /mt6_cv.rds \+
 /mt7_cv.rds \+ /mt8_cv.rds \+ /mt9_cv.rds \+ /mt10_cv.rds \* /richness:
 Species richness random forest model tune files \+ /mt1_cv.rds \+
 /mt2_cv.rds \+ /mt3_cv.rds \+ /mt4_cv.rds \+ /mt5_cv.rds \+ /mt6_cv.rds \+
 /mt7_cv.rds \+ /mt8_cv.rds \+ /mt9_cv.rds \+ /mt10_cv.rds \* /YES_VERT:
 Similar structure as NO_VERT but with vertebrate (tetrapod) data. \*
 /model_fits \+ /endemism_imp.csv \+ /endemism_imp.rds \+ /PD_imp.csv \+
 /PD_imp.rds \+ /PE_imp.csv \+ /PE_imp.rds \+ /richness_imp.csv \+
 /richness_imp.rds \* /model_results \* /endemism: Model results for
 weighted endemism \+ /endemism_log_emp.tif \+ /endemism_log_ext.tif \*
 /PD: Model results for phylogenetic diversity \+ /PD_ext.tif \+ /PD.tif \*
 /PE: Model results for phylogenetic endemism \+ /PE_log_emp.tif \+
 /PE_log_ext.tif \* /richness: Model results for species richness \+
 /richness_ext.tif \+ /richness.tif \* /model_tune \* /endemism: weighted
 endemism random forest model tune files \+ /mt1_cv.rds \+ /mt2_cv.rds \+
 /mt3_cv.rds \+ /mt4_cv.rds \+ /mt5_cv.rds \+ /mt6_cv.rds \+ /mt7_cv.rds \+
 /mt8_cv.rds \+ /mt9_cv.rds \+ /mt10_cv.rds \* /PD: Phylogenetic diversity
 random forest model tune files \+ /mt1_cv.rds \+ /mt2_cv.rds \+
 /mt3_cv.rds \+ /mt4_cv.rds \+ /mt5_cv.rds \+ /mt6_cv.rds \+ /mt7_cv.rds \+
 /mt8_cv.rds \+ /mt9_cv.rds \+ /mt10_cv.rds \* /PE: Phylogenetic endemism
 random forest model tune files \+ /mt1_cv.rds \+ /mt2_cv.rds \+
 /mt3_cv.rds \+ /mt4_cv.rds \+ /mt5_cv.rds \+ /mt6_cv.rds \+ /mt7_cv.rds \+
 /mt8_cv.rds \+ /mt9_cv.rds \+ /mt10_cv.rds \* /richness: : Species
 richness random forest model tune files \+ /mt1_cv.rds \+ /mt2_cv.rds \+
 /mt3_cv.rds \+ /mt4_cv.rds \+ /mt5_cv.rds \+ /mt6_cv.rds \+ /mt7_cv.rds \+
 /mt8_cv.rds \+ /mt9_cv.rds \+ /mt10_cv.rds \* /predictions \+
 /endemism.tif \+ /endemism.tif.aux.xml \+ /PD.tif \+ /PD.tif.aux.xml \+
 /PE.tif \+ /PE.tif.aux.xml \+ /richness.tif \+ /richness.tif.aux.xml \*
 /RMSE_barplots: Folder with dataframes for plotting barplots for root mean
 square error \+ /endemism.csv \+ /PD.csv \+ /PE.csv \+ /richness.csv \*
 /RMSE_spatialBlocks: Error metrics from the spatial block cross-validation
 process for endemism, PD, PE, and richness. \+ /rmse_sb_endemism.gpkg \+
 /rmse_sb_PD.gpkg \+ /rmse_sb_PE.gpkg \+ /rmse_sb_richness.gpkg \*
 /spatialBlocks: Spatial blocks for tuning random forest models \+
 /spatialBlocks_endemism.rds \+ /spatialBlocks_PD.rds \+
 /spatialBlocks_PE.rds \+ /spatialBlocks_richness.rds **i.
 supplementary_analyses**: Contains supplementary analysis files: \+
 /SBMvsFORBIDDING.tif: Spatial comparison file. \+
 /SBMvsFORBIDDING.tif.aux.xml **j. TETRAPODS**: Data for amphibians, birds,
 mammals, and reptiles biodiversity metrics \* /wag4: TIF files related to
 endemism, phylogenetic diversity, phylogenetic endemism, and species
 richness for each tetrapod group. \* /wgs84: Similar structure but
 organized with geographic projections and threat data. \+
 /amphibian_endemism.tif \+ /amphibian_richness.tif \+
 /amphibian_threat.tif \+ /bird_endemism.tif \+ /bird_richness.tif \+
 /bird_threat.tif \+ /mammal_endemism.tif \+ /mammal_richness.tif \+
 /mammal_threat.tif \+ /reptile_endemism.tif \+ /reptile_richness.tif \+
 /reptile_threat.tif 3. `SDM_set1`: Raw output file from species
 distribution modeling for a subset of the data **a. polygons**: a folder
 with a subset of individual species polygons in geopackage format. Each
 species polygon has two fields: binomial, indicating the species binomial
 name, and origin corresponding to the native status of the species (1 for
 native species). **b. rasters**: folder with individual modeled spatial
 rasters in GeoTiff format at a grain resolution of 5–arc min (~9 km at the
 equator). **c. raw_rasters**: folder with raster layers in GeoTiff format
 at a grain resolution of 5–arc min (~9 km at the equator) of the raw
 suitability scores to explore alternative approaches of thresholding the
 predictions. 4. `CODES` /CODES.zip: Zipped folder with the R scripts used
 to generate the figures. \+ /Fig2_multipanel.R \+ /Fig3_cryptic_hotspots.R
 \+ /Fig4_protected_hotspots.R \+ /Fig5_functional_traits_CohensD.R \+
 /functions.R `phyloregion_1.0.9.tar.gz` (R package used to run most of the
 analyses) ## Funding: This work was supported by the US National Science
 Foundation (awards 2345994, and 2416314).