Exposure of mammal genetic diversity to mid-21st century global change Research

Accelerating climate and land-use change are rapidly transforming Earth’s biodiversity. While there is substantial evidence on the exposure and vulnerability of biodiversity to global change at the species level, the global exposure of intraspecific genetic diversity (GD) is still unknown. Here, we assess the exposure of mitochondrial GD to mid-21st century climate and land-use change in terrestrial mammal assemblages at grid-cell and bioclimatic region scales under alternative narratives of future societal development. We used global predictions of mammal GD distribution based on thousands of georeferenced mitochondrial genes for hundreds of mammal species, the latest generation of global climate models from the ongoing sixth phase of the Coupled Model Intercomparison Project (CMIP6), and global future projections of land-use prepared for CMIP6. We found that more than 50% of the genetically poorest geographic areas (grid-cells), primarily distributed in tundra, boreal forests/taiga and temperate bioclimatic regions, will be exposed to mean annual temperature rise that exceeds 2 ° C compared to the baseline period under all considered future scenarios. We also show that at least 30% of the most genetically rich areas in tropical, subtropical and montane regions will be exposed to an increase of mean annual temperature > 2 ° C under less optimal scenarios. Genetic diversity in these rich regions is also predicted to be exposed to severe reductions of primary vegetation area and increasing human activities (an average loss of 5–10% of their total area under the less sustainable land-use scenarios). Our findings reveal a substantial exposure of mammal GD to the combined effects of global climate and land-use change. Meanwhile the post-2020 conservation goals are overlooking genetic diversity, our study identifies both genetically poor and highly diverse areas severely exposed to global change, paving the road to better estimate the geography of biodiversity vulnerability to global change.


Introduction
The exposure of terrestrial biodiversity to major anthropogenic pressures, such as climate and land-use change, is steadily increasing (Newbold et al. 2015, Foden et al. 2019, and its adverse effects on species and ecosystems are predicted to rapidly escalate in the absence of drastic mitigation efforts (Tilman et al. 2017, Powers and Jetz 2019, Trisos et al. 2020). Intensified land-use change, including deforestation, agricultural expansion and urbanization, in the last decades has been the main driver of terrestrial biodiversity decline, with its impacts possibly surpassing the proposed safety boundaries for the long-term maintenance of ecosystem functions at global scale (Steffen et al. 2015, Newbold et al. 2016. Additionally, climate change has already impacted biodiversity (Scheffers et al. 2016), while its effects will likely match or exceed those of land-use in the coming decades (Newbold 2018). Global annual temperature is currently rising by 0.2°C (± 0.1°C) per decade, with many land regions already experiencing annual temperatures that exceed 1.5°C above pre-industrial levels (Allen et al. 2018). An increase of more than 2°C in annual temperature would significantly elevate species abrupt exposure (Trisos et al. 2020) and extinction risk (Urban 2015, Warren et al. 2018) and could trigger continued warming of the Earth system, causing serious disruptions to ecosystems globally (Steffen et al. 2018).
In the past decades, increased data availability on species distributions and their evolutionary relationships has enabled thorough assessments of the exposure and vulnerability of terrestrial species and their diversity to future environmental change, both at regional (Pio et al. 2014, Dullinger et al. 2020) and global scales (Hof et al. 2011, Pacifici et al. 2018, Carvalho et al. 2019b, Trisos et al. 2020). Yet, a significant gap of knowledge exists regarding the exposure of biodiversity below the species level, and in particular for intraspecific genetic diversity (GD). Genetic variation can facilitate adaptation of species under rapidly changing environments (Zheng et al. 2019, Bitter et al. 2019) and promote ecosystem resilience (Oliver et al. 2015). Moreover, the rapidly accumulating variation in animal mitochondrial genes has been traditionally considered 'neutral' and is extensively used to estimate evolutionary histories of individuals and populations (Avise 2009). However, this assumption of neutrality is increasingly being challenged. Specifically, with regard to global warming, a growing body of research highlights the dual role of variation in the mitochondrial genome both in reflecting local adaptation to climate and in determining the thermal range that an organism can tolerate (Camus et al. 2017, Lasne et al. 2019, Li et al. 2019. High levels of mitochondrial GD, particularly in warmer climates, are indicative of potentially beneficial standing variation that could encode for increased heat tolerance and buffer against rapid warming (Camus et al. 2017, Li et al. 2019). On the contrary, a homogeneous pool of genetic variation adapted to colder climates suggests reduced evolutionary potential of species to respond to rising temperatures (Hoffmann andSgrò 2011, Camus et al. 2017). Since GD is a vital component of biodiversity, there is a pressing need to assess its exposure under ongoing and future climate and landuse change. Recent studies have assessed the vulnerability of GD to climate change focusing either on single species (Theodoridis et al. 2018), or on small species assemblages at regional scale (Carvalho et al. 2019a), yet global assessments within large taxonomic groups are lacking.
For terrestrial mammals, an iconic group in conservation (McGowan et al. 2020), large efforts in data collection and compilation have identified conservation priority areas and hotspots of vulnerability to future environmental change for three key dimensions of mammal biodiversity, i.e. taxonomic, phylogenetic and phenotypic diversity (Brum et al. 2017, Pacifici et al. 2018). Yet, only recently a description of the global distribution of genetic variation in terrestrial mammals became available (Miraldo et al. 2016). The generation of large-scale geo-referenced sequences for two widely-used mitochondrial genes, i.e. cytochrome b (cytb) and cytochrome oxidase 1 (co1), revealed a high degree of covariation between GD and interspecific diversity and the major role of rapid climate change during the last millennia in reducing current levels of diversity in co1 (Theodoridis et al. 2020). These results coupled with the suggested relevance of standing mitochondrial variation to climatic adaptation supports the interest for a first evaluation of the global exposure of mammal GD to future climate and land-use change, potentially providing crucial information to anticipate the fate of the most basal dimension of biodiversity.
To date, most of the assessments of biodiversity exposure and vulnerability to climate change have relied on climate simulations from either the most recently completed phase (fifth) of the Coupled Model Intercomparison Project (CMIP5, IPCC 2014), or previous phases. Yet, the latest generation of global climate models (GCMs) from the ongoing sixth phase (CMIP6) shows a substantially increased sensitivity to atmospheric CO 2 concentration than did the previous generation of models, resulting both in wider inter-model variation and overall greater projected warming through the twenty-first century (Forster et al. 2020, Tokarska et al. 2020, Zelinka et al. 2020. These new climate projections are driven by a new set of energy, land-use and emissions scenarios produced with integrated assessment models that combine narratives of future societal development, namely shared socioeconomic pathways (SSPs), and climate mitigation policies, namely representative concentration pathways (RCPs, O'Neill et al. 2016). Multiple combinations of SSPs and RCPs were developed in preparation for CMIP6 spanning a wide range of scenarios, from sustainable growth (SSP1) coupled with low warming (RCP2.6), to fossil-fueled development (SSP5) coupled with high warming by the end of the century (RCP8.5, O'Neill et al. 2016). Importantly, the updated SSP-driven climate projections are in line with the most recent future projections of land-use change, the land use harmonization 2 (LUH2) dataset, prepared for CMIP6 (Lawrence et al. 2016), thus enabling an integrated assessment of biodiversity exposure to alternative scenarios of future environmental change.
Here, we assess the exposure of GD in terrestrial mammal assemblages to mid-21st century climate and land-use change under the four basic Tier 1 scenarios, namely SSP1-2.6, SSP2-4.5, SSP3-7.0 and SSP5-8.5, developed for CMIP6. Following Foden et al. (2019), we define exposure as the magnitude of environmental change experienced by living organisms within a geographic area. We take advantage of the recently-published modeled distribution of mitochondrial GD at global scale (Theodoridis et al. 2020) to estimate exposure in regions of the world that lack adequate genetic data. To further account for the uncertainty in the predicted distribution of GD, we also estimate the exposure of measured GD in geographic areas (i.e. grid-cells) based on thousands of mitochondrial sequences and hundreds of mammal species across the globe. We estimated exposure across grid-cell and bioclimatic regions and show that regions of the world that harbor high GD, such as the tropical bioclimatic regions, are largely exposed both to climate and land-use change, while genetically poor regions, such as the boreal forests, are mainly exposed to climate change.

Mammal genetic diversity
Due to the advantageous properties (e.g. rapid evolution and low recombination) of mitochondrial DNA in mammals, we utilized two mitochondrial genes, i.e. cytb and co1, that have been extensively used in taxonomic, phylogenetic and phylogeographic studies (Avise 2009), and therefore constitute the richest resource of genetic data with available spatial information. Additionally, standing genetic variation in these genes has been recently suggested to affect molecular mechanisms that could provide evolutionary adaptation under climatic stress (Camus et al. 2017, Li et al. 2019. Mitochondrial sequences were retrieved from GeneBank and BOLD and geographic coordinates to the sequences that were not already annotated with latitude and longitude were assigned using the API tool provided by GeoNames. org (<http://api.geonames.org>). A detailed description of the retrieval and georeferencing of genetic sequences is provided in Miraldo et al. (2016) and Theodoridis et al. (2020). Following Miraldo et al. (2016), GD was defined as the average number of nucleotide differences per geographic area (grid-cells) across all pairwise sequence comparisons within a species. For a particular assemblage of species, GD was then defined as the average nucleotide diversity per grid-cell (385.9 × 385.9 km spatial resolution) across all species present in that assemblage. The choice of the spatial resolution of grid-cells was based on the following criteria: 1) maximization of the number of sequences from each species that are included in each cell, 2) minimization of the difference in the number of sequences between species and 3) minimization of the georeferencing errors deriving from the broad locality descriptions attached to each retrieved sequence (Miraldo et al. 2016).
To estimate the exposure of mammal GD in future climate and land-use change, we utilized the recently-published predicted distribution of mitochondrial GD in terrestrial mammal assemblages at global scale. These global maps (1375 cells, Fig. 1a) are based on the predictions made by long-standing biodiversity theories that invoke past population divergence and extirpations as the intermediate drivers of the distribution of GD at global scale. Multiple predictors, including interspecific diversity and past climate change, were evaluated in a multi-model inference framework and models were fitted using selected grid-cells based on data availability (Theodoridis et al. 2020 for details). The most parsimonious model for cytb was fitted to a total of 185 grid-cells (15 706 sequences, 581 species, Supporting information) and included phylogenetic diversity as the only predictor, while the model for co1 was fitted to a total of 76 grid-cells (16 842 sequences, 492 species, Supporting information) and included both phylogenetic diversity and trends in late quaternary temperature change as predictors. Utilizing the predicted distribution of GD allowed us to include data-poor regions of the world, such as the tropical and subtropical regions (Theodoridis et al. 2020), thus enabling assessments of the exposure of mammal GD at global scale. To account for uncertainty (i.e. the unexplained variance) in the predicted global distribution of GD and validate our assessments, we also assessed exposure based only on the aforementioned observed distribution of GD for both genetic markers. The maps of predicted and observed GD for cytb and co1 were downloaded from <https://github. com/spyrostheodoridis/Genetic-geography-of-terrestrial-mammals>. In order to assess the exposure to climate and land-use change of grid-cells with higher or lower genetic diversity than expected given their species diversity, we further regressed GD against mammal species richness (as estimated in Theodoridis et al. 2020) and mapped the residuals for both mitochondrial genes.
Since bioclimatic regions are primarily defined by climate and land-cover (i.e. vegetation type, Olson et al. 2001), mammal assemblages within a region may be subject to similar exposure (i.e. dependent with respect to exposure). Therefore, we additionally assessed the exposure of GD at the regional scale using the classification of the WWF Terrestrial Ecoregions (Olson et al. 2001, 14 terrestrial bioclimatic regions, Fig. 1b). The predicted GD per region was calculated by intersecting the grid-cells from the predicted maps described above with the bioclimatic region map. In cases where several regions intersected with a grid-cell, the cell was assigned to the region with the greatest area overlap with the respective cell.

Future climate change scenarios
We estimated mid-21st century climate change as the magnitude of change in annual temperature and precipitation between the baseline period 1961-1990 and the future period 2040-2069. These reference periods are adopted by the Intergovernmental Panel on Climate Change (IPCC, <www.ipcc-data.org/ddc/ddc_faq.html>) as representative of 20th and mid-21st century climate.
Increasing annual temperature has been recently suggested to be the primary driver of exposure of species assemblages to climate change (Trisos et al. 2020), while long-term historical variation in both annual temperature and precipitation has been shown to play a major role in shaping current biodiversity patterns at multiple dimensions (Brown et al. 2020, Theodoridis et al. 2020. To estimate annual climate for the baseline period, we first sourced observed land-surface monthly temperature and precipitation from the latest version of the Climatic Research Unit Time-Series (CRU TS Ver. 4.03, <http://data.ceda. ac.uk/badc/cru/data/cru_ts/cru_ts_4.03>) at 0.5° × 0.5° (latitude/longitude) spatial resolution (Harris et al. 2014) and then calculated the mean annual temperature and total annual precipitation for the years 1900-2018. To match the spatial resolution of future climate (below), we re-gridded the CRU annual values to a 1.4° × 1.4° using distance-weighted average remapping. Baseline climate was defined as the average values of annual temperature and precipitation across the years .
Future projections of monthly surface temperature and precipitation were sourced from the latest generation of GCMs that are being currently released as part of the ongoing sixth phase of the Coupled Model Intercomparison Project (CMIP6, available at <https://esgf-node.llnl.gov/ search/cmip6>, accessed on 19 February 2020). We selected a total of 19 GCM's with at least one simulation/variant available for each of the four scenarios (i.e. SSP1-2.6, SSP2-4.5, SSP3-7.0 and SSP5-8.5), aiming at maximizing the number of models and thus the overall multi-model performance under each scenario (Flato et al. 2013). The 19 selected GCMs represent different modelling institutions so as to minimize model interdependence in the final model ensembles (Abramowitz et al. 2019). For models with multiple climate simulations per scenario (Table 1 for details), we averaged across all simulations to minimize internal model variability and uncertainty (Eyring et al. 2019). Monthly temperature and precipitation data were obtained for each of the 19 GCMs and the four basic scenarios and re-gridded to a 1.4° × 1.4° (latitude/longitude, the median longitude resolution across all models, Table 1) global grid using bilinear interpolation. Monthly climate variables were then converted to °C (from Kelvin) and mm/month (from kg m 2 s 1 ), respectively, and summarized to mean annual 821 temperature and total annual precipitation. To preserve the internal climate variability and uncertainty across GCMs, we first calculated the magnitude of climate change (i.e. annual anomalies) for each individual model and scenario by subtracting the baseline climate  from future climate and then averaged the anomalies across all 19 GCMs (i.e. mutli-model ensemble of annual anomalies, Flato et al. 2013). Finally, we averaged the ensemble anomaly values across the years 2040-2069 to obtain the predicted mid-21st century magnitude of climate change.
We evaluated the overall performance of the multi-model ensemble in predicting observed climate by comparing historical runs  for all 19 GCMs with observed climate. We calculated annual temperatures predicted by each model, as described above, and then averaged the annual values across all 19 models (multi-model ensemble of annual temperatures). Following the IPCC procedures (Flato et al. 2013), we estimated ensemble bias as the difference between the multi-model ensemble and the observed temperature from the CRU dataset for the baseline years . Additionally, we visually compared the global weighted mean of the mutli-model ensemble of annual anomalies (see above for calculation), with the observed anomalies for the years 1900-2015 (Fig. 2a).

Future land-use change scenarios
We estimated GD exposure to mid-21st century landuse change as the difference in major land-use/land-cover categories (below) between current conditions (i.e. the year 2015 sourced from the common historic dataset v2h) and mid-21st century (i.e. 2050, v2f ) predictions under different scenarios. Global land-use data for current conditions and 2050 projections were sourced from the landuse harmonization 2 (LUH2) dataset prepared for CMIP6 (Lawrence et al. 2016, <https://luh.umd.edu/data.shtml>). The LUH2 project models the states (fractions ranging from 0 to 1) of 12 land-use categories for the period spanning 2015-2100 at 0.25° × 0.25° spatial resolution. The four land-use projections used in this study were generated under the same future scenarios that were used to generate the CMIP6 climate projections For each scenario, we summed the states of the 12 available land-use categories in three general categories (Supporting information) based on their similarities with regard to changes in vertebrate biodiversity (Newbold et al. 2015, Newbold 2018: primary vegetation (pristine habitat with no record of destruction), secondary vegetation (natural habitat recovering after some recorded historical destruction) and human activities (plantations, croplands, pastures and areas of human settlement). Here, we only consider primary vegetation and human activities because changes in these two major classes capture the primary land-use threats to vertebrate biodiversity globally (Newbold et al. 2015, Carvalho et al. 2019b) and are correlated with changes in secondary land. We then calculated the difference in the frequency of primary vegetation and human activities between 2015 and 2050.

Exposure of mammal GD to climate and land-use change
The exposure of GD to future climate change was assessed by intersecting the predicted mid-21st century anomalies in annual temperature and precipitation (1.4° × 1.4°) and landuse (0.25° × 0.25°) with the maps of GD. Exposure within each GD grid-cell was calculated as the area-weighted mean of all anomaly values intersecting the respective grid-cell. The weight for each anomaly grid-cell was calculated using the following formula: where y is the absolute latitude of the cell centroid.
Because of reduced sample size (i.e. number of grid-cells) for the observed GD, we additionally conducted all statistical analyses classifying the observed GD in two percentiles (0-50% and 50-100%). Due to prominent unequal variance (i.e. heteroskedasticity) for all exposure variables across the GD quartiles (and percentiles for the observed GD), we first applied one-way Welch's ANOVA across classes, and upon rejection of the null hypothesis of equal means we further applied the post hoc Welch's t-test with Bonferroni p value corrections for significant differences between the lowest and highest class (hereafter low and high GD respectively). To further test for significant differences in the medians, we applied the non-parametric ANOVA (Kruskal-Wallis test) across GD classes followed by the post hoc Dunn-Bonferroni pairwise test.
We further assessed the exposure of GD to climate and land-use change at the bioclimatic region scale by applying mixed linear models treating exposure as response variable, and GD and region as fixed and random effects respectively. To this end, each grid-cell was assigned a bioclimatic  . Numbers in parentheses indicate the total number of model variants used for each future scenario and the common historical ensemble. Grey zones indicate the baseline and the future reference periods . Shades around the lines represent one standard deviation (across all 19 GCMs) from the multi-model ensemble mean. (b) Multi-model ensemble bias in annual temperature (the difference between the multimodel ensemble mean and the observed temperature from the CRU TS dataset for the baseline period. (c) Global predicted distribution of mean annual temperature and annual precipitation change between the baseline and future periods. region as described above ('Mammal genetic diversity' section). We fitted the models using the raw exposure variables, with the exception of primary vegetation change, where we used the square root of its absolute value across grid-cells to better meet the normality assumptions of the mixed linear model. The significance of the variation in GD exposure across regions was assessed using likelihood ratio test between the mixed model and the simple model that included only GD as independent variable (Harrison et al. 2018). We additionally assessed the mixed model fit for all exposure variables using the conditional R 2 (R 2 c, Nakagawa and Schielzeth 2013).

Global patterns of mammal genetic diversity
For both mitochondrial genes, the grid-cells in the lowest GD quartile (low GD) are primarily distributed across high latitudes of the Northern hemisphere and in the deserts and xeric shrublands of North Africa and Arabian Peninsula, while grid-cells in the highest GD quartiles (high GD) are distributed in tropical and subtropical latitudes (Fig. 1a, 3, 4, Supporting information). Additionally, for cytb, regions that encompass low GD include the deserts of Australia and southern Patagonia (Fig. 1a, 3, 4a). For co1, low GD is also distributed in Northern Europe, including Scandinavia and the British Isles, and in North America (Fig. 1a, 4b, Supporting information).
In agreement with the reported patterns at the grid-cell scale, tropical and subtropical bioclimatic regions showed a substantially higher number of predicted average mutations per base pair compared to temperate and polar regions (Fig. 1b, Supporting information). For both mitochondrial genes, mangroves harbored the highest median and mean genetic diversity across regions, followed by the tropical and subtropical coniferous forests and the tropical and subtropical grasslands, savannas and shrublands. Additionally, for both genes the montane grasslands and shrublands exhibited high GD, with maximum values exceeding those of tropical and subtropical forests (Fig. 1, Supporting information). The bioclimatic regions harboring the lowest GD for both genes were boreal forests/taiga and tundra. Genetic diversity for co1 was noticeably lower compared to cytb for all temperate and polar regions. Figure 3. Exposure of mammal cytb genetic diversity to climate and land-use change across three future scenarios. In the color scales, the x-axis indicates genetic diversity split in quartiles. The y-axis indicates the change (i.e. exposure) between future and baseline values for each climate and land-use variable.

Future climate change
The multi-model ensemble (19 GCMs, Table 1) of the new CMIP6 climate projections predict an increase in global terrestrial annual temperature of > 2°C by mid-21st century under all four future scenarios (Fig. 2a). Overall, the multimodel ensemble of historical runs performed well in predicting the observed baseline annual temperature (Fig. 2b). For the majority of the land grid-cells, the multi-model ensemble bias in annual temperature falls within −1.79 and 1.43°C (mean = −0.18, SD = 1.61, Fig. 2b, Supporting information).
Deviations from observed values were particularly pronounced in high-latitude regions, such as Greenland (negative bias), and central East Siberia (positive bias, Fig. 2b). High ensemble biases were also present in some mountainous regions, such as the central and southern Andes (positive bias) and the Tibetan plateau (negative bias, Fig. 2b).
Mid-21st century change (i.e. average values of anomalies between 2040 and 2069) in annual temperature showed a consistent spatial pattern under all four future scenarios (Fig. 2c). Overall, the model ensembles predict higher temperature change at higher latitudes, with the exception of southern Greenland where ensembles predict lower temperatures compared to the baseline period (Fig. 2c). The magnitude of change increased gradually from SSP1-2.6 (mean anomaly across grid-cells = 2.32°C, SD = 1.88), to SSP5-8.5 (mean = 3.51°C, SD = 2.08, Supporting information). Mid-21st century annual precipitation change also showed a consistent spatial pattern under the four scenarios (Fig. 2c). Overall, model ensembles predicted a global increase in annual precipitation, with the magnitude of change being largely similar under the four scenarios Supporting information), although areas like the north-eastern Amazon are predicted to experience pronounced reductions in precipitation.

Exposure of genetic diversity to global change across grid-cells
Highly diverse regions in parts of tropical Africa, South America and Southeast Asia were predicted to lose > 5% of their primary vegetation area with many grid-cells, particularly in tropical Africa, being also exposed to severe warming of > 2°C under most considered scenarios and for both mitochondrial genes (Fig. 3, Supporting information). The majority of the less genetically diverse regions in temperate and high latitudes are exposed to climate warming of > 2°C, with many boreal regions also predicted to lose > 5% of their primary vegetation area by the mid-21st century (Fig. 3, Supporting information).
We identified significant differences in mid-21st century annual temperature change between low and high GD regions under all four future scenarios and for both mitochondrial genes (post hoc Welch's t-test p values < 0.05, Kruskal-Wallis test p values < 0.05, Fig. 4). For both genes, the vast majority of the low GD grid-cells were found to be exposed to temperatures > 2°C compared to the baseline period under the four scenarios (Fig. 3, 4, Supporting information). For cytb, the magnitude of exposure for low GD regions increased from SSP1-2.6 (median = 2.38), to SSP2-4.5 (median = 2.77), to SSP3-7.0 (median = 3.03), to SSP5-8.5 (median = 3.37), while the majority of grid-cells with high GD for cytb were found to be exposed to annual temperature changes of > 1.5°C under SSP2-4.5 (median = 1.75), SSP3-7.0 (median = 1.99) and SSP5-8.5 (median = 2.25, Fig. 3, 4a). These patterns were similar for co1 (Fig. 4b, Supporting information).
Significant differences in the means in annual precipitation change between low and high GD for cytb (post hoc Welch's t-test p value < 0.05) were observed only under SSP1-2.6 (mean exposure of low GD = −21.69 mm, mean exposure of high GD = 58.07 mm) and SSP2-4.5 (mean exposure of low GD = −16.13 mm, mean exposure of high GD = 51.89 mm, Fig. 4a), and only under SSP1-2.6 for co1 (mean exposure of low GD = −2.51 mm, mean exposure of high GD = 53.99 mm, Fig. 4b). Additionally, under all four future scenarios the majority of both low and high GD grid-cells are exposed to higher precipitation levels compared to the baseline period (Fig. 3, 4, Supporting information), with median values in precipitation change ranging from 92.13 mm (low GD under SSP1-2.6) to 120.56 mm (high GD under SSP5-8.5) for cytb ( Fig. 4a), and from 88.68 mm (high GD under SSP1-2.6) to 152.51 mm (low GD under SSP5-8.5) for co1 (Fig. 4b).
We also found that high GD regions are significantly more exposed to increased human activity by mid-21st century compared to low GD regions across most pairwise comparisons, i.e. future scenarios and genes (pos hoc Welch's t-test p values < 0.05, Kruskal-Wallis test p values < 0.05, Fig. 3, 4, Supporting information). The exceptions were the pairwise comparisons under the scenario SSP1-2.6, where the majority of high GD grid-cells for both genes were predicted to exhibit a small reduction in human activities (cytb: high GD median = −0.006, co1: high GD median = −0.005, Fig. 4). The highest exposure in human activities within high GD grid-cells was observed under SSP3-7.0, where most of these regions were predicted to experience an increase in human activities of > 0.035 (3.5%) of their total area (cytb high GD median = 0.036, co1 high GD median = 0.039, Fig. 4). Under this scenario, the vast majority of regions in tropical Africa were predicted to experience an increase in human activities of > 0.5 (i.e. 5%) of their total area (Fig. 3, Supporting information).
The exposure patterns to climate and land-use change identified using the predicted GD distribution matched, to a large extent, the exposure patterns identified using the observed GD distribution (Supporting information). Moreover, when using the observed GD, some exposure patterns were even more pronounced, such as pairwise differences in annual precipitation change for cytb and in annual temperature change for co1, most probably reflecting the low sample size (i.e. number of grid-cells) for observed GD (Supporting information).
Parts of the tropical and subtropical regions of the planet, mostly distributed in central, eastern and southern Africa, and south Asia harbor disproportionally higher levels of GD than expected by species richness (Supporting information). Those regions were found to be highly exposed to temperature rise of > 2°C for both mitochondrial genes. A large proportion of grid-cells in these regions was also found to be exposed to loss of primary vegetation and increased human activities, particularly under the SSP3-7.0 scenario (Supporting information). In addition, eastern North America is also found to harbor higher than expected GD and exposed both to severe warming and reduced precipitation by mid-21st century. Highly exposed regions to temperature rise that harbor lower levels of GD than predicted by species richness are primarily distributed in East Asia and western North America (Supporting information).

Exposure of genetic diversity to global change across bioclimatic regions
We identified significant variation in GD exposure to annual temperature change (p lrt < 0.01) across bioclimatic regions under all scenarios and for both genes (Fig. 5). Regions harboring the lowest levels of GD, including tundra, boreal forests/taiga, temperate broadleaf and mixed forests, temperate conifer forests, Mediterranean forests, woodlands and scrub, deserts and xeric shrublands, temperate grasslands, savannas and shrublands and temperate conifer forests, exhibited higher mean exposure (> 2°C under SSP2-4.5, SSP3-7.0 and SSP5-8.5 for both genes) compared to tropical and subtropical regions (Fig. 5). The fit of the mixed models increased gradually from SSP1-2.6 (cytb R 2 c = 0.169, co1 R 2 c = 0.163) to SSP5-8.5 (cytb R 2 c = 0.26, co1 R 2 c = 0.257). Moreover, GD at the bioclimatic region scale showed significant variation in exposure to annual precipitation change (p lrt < 0.01) across scenarios and genes, with tundra, boreal forests/taiga, temperate conifer forests and montane grasslands and shrublands showing on average an increase of > 150 mm by mid-21st century and under all considered scenarios (Fig. 3, 5, Supporting information). Several regions harboring high GD, including mangroves, tropical and subtropical grasslands, savannas and shrublands, tropical and subtropical moist broadleaf forests and flooded grasslands and savannas exhibited on average a reduction in precipitation (0-150 mm, Fig. 3, 5, Supporting information). The fit of the mixed models for annual precipitation exposure was similar under all four future scenarios ranging for cytb from R 2 c = 0.178 (SSP5-8.5) to R 2 c = 0.188 (SSP1-2.6, Fig. 5a), and from for co1 R 2 c = 0.205 (SSP5-8.5) to R 2 c = 0.216 (SSP1-2.6, Fig. 5b).
Genetic diversity at the bioclimatic region scale showed significant variation in exposure to primary vegetation change (p lrt < 0.01) under all scenarios and for both genes, with a mean reduction in primary vegetation predicted across all regions (Fig. 5). The fit of mixed models was higher under SSP2-4.5 (cytb R 2 c = 0.148, co1 R 2 c = 0.67), SSP3-7.0 (cytb R 2 c = 0.175, co1 R 2 c = 0.191) and SSP5-8.5 (cytb R 2 c = 0.167, co1 R 2 c = 0.189), where the majority of regions harboring high GD, including mangroves, tropical and subtropical grasslands, savannas and shrublands, tropical and subtropical coniferous forests, tropical and subtropical moist broadleaf forests and flooded grasslands and savannas, exhibited a mean loss of primary vegetation between 0.05 (i.e. 5%) and 0.1 (i.e. 10%) of their total area (Fig. 5). Additionally, we identified significant variation across regions (p lrt < 0.01) in the exposure of GD to increased human activities under SSP2-4.5 (cytb R 2 c = 0.079, co1 R 2 c = 0.08), SSP3-7.0 (cytb R 2 c = 0.154, co1 R 2 c = 0171) and SSP5-8.5 (cytb R 2 c = 0.104, co1 R 2 c = 0.116). Under these three less sustainable scenarios, tropical and subtropical bioclimatic regions rich in GD showed an average increased exposure to human activities (with the exception of tropical and subtropical coniferous forests under SSP5-8.5), while under SSP1-2.6 all regions exhibited a reduction in areas impacted by human activities (Fig. 5).

Discussion
Much of the global pool of mitochondrial GD is found to be exposed to annual temperature change above the critical threshold of 2°C degrees by mid-21st century, and this exposure is significantly higher for genetically poor areas mostly located at high latitude regions. Additionally, tropical and subtropical bioclimatic regions rich in GD are predicted to suffer severe losses of primary vegetation area due to deforestation and cropland expansion, while many of these regions are also exposed to severe future warming. Evidence buried in paleo-archives, including ancient DNA (Lorenzen et al. 2011, Botta et al. 2019, Fordham et al. 2020), but also in historical records over the last century (Rubidge et al. 2012) shows that mammal species have experienced significant bottlenecks and losses of genetic diversity, under periods of rapid warming. If similar declines are to be expected in the future, global conservation efforts should also attempt to maximize the protection of mammal intraspecific diversity under global environmental change (Laikre et al. 2020).

Exposure of genetic diversity to climate change
The exposure of GD to future climate change decreases from the poles to the tropics. High latitude regions of the world were exposed to rapid past climatic changes (Brown et al. 2020), eroding their biodiversity at all levels, including genetic diversity (Theodoridis et al. 2020). Tundra and boreal forests/taiga will experience a significant exposure to climate change, with mean annual temperature exceeding the 2°C threshold above the baseline period under all considered future scenarios. Warming in these low GD regions is twice the global average, and the overall low intraspecific genetic variation in these mammal communities poses great challenges for their adaptation to climate change, likely leading to severe local extinctions (Mills et al. 2018). Climate warming, and the resulting declines in suitable habitat (e.g. sea-ice extent, seasonal snow duration, declining tundra flora) have already affected the survival of terrestrial mammals inhabiting these regions (Pagano et al. 2018, Post et al. 2019, Richter-Menge et al. 2019. The predicted warming in these regions, coupled with the predicted changes in annual precipitation (increase of > 150 mm on average under all scenarios) and the subsequent rise in humidity and frequency of rain-on-snow episodes, will further reduce food availability of mammal species at all trophic levels, from herbivores to carnivores (Gilg et al. 2009, Domine et al. 2018. Additionally, mammal assemblages in the Sahara Desert, Arabic Peninsula and the deserts of Australia were also found to harbor low levels of GD. Despite being arid-adapted and exposed to less severe warming compared to higher latitude assemblages, desert mammals may already be living close to their upper thermal limits and the lack of GD may well hinder further adaptation to rising temperatures, increasing vulnerability and extinction risk (Vale and Brito 2015). Within each plot, the R 2 c value indicates the conditional variance explained by a linear mixed model including exposure (y-axis) as independent variable, and GD (x-axis) and bioclimatic region (color) as fixed and random effects respectively. p-values (P lrt ) indicate significance of the likelihood ratio test between the mixed and the simple (i.e. without the random effect/bioclimatic region) model. The error bars around the points indicate the 25 and 75 percentiles.
The temperate zones of North America, Europe and central Asia, regions of intermediate levels of GD, have been characterized by the lowest species extinction risks globally (Urban 2015). Yet, owing to ongoing climate change, an accountable proportion of mammals in these regions has already experienced range contractions and local extirpations (Pacifici et al. 2017), and these negative effects are predicted to increase under current warming trends (Maiorano et al. 2011). These temperate bioclimatic regions are predicted to experience an overall exposure in temperature change above 2°C under the three less sustainable, but more realistic, scenarios SSP2-4.5, SSP3-7.0 and SSP5-8.5, indicative of the increased vulnerability of mammal GD in these regions.
Highly genetically diverse tropical and subtropical mammal assemblages are less exposed to rising annual temperatures compared to temperate and polar assemblages, yet many of these regions, particularly in Africa and South America, are exposed to annual temperature change above 2°C. Additionally, many tropical mammal assemblages are predicted to experience an overall reduction in annual precipitation, particularly those occurring in the Amazon basin. Mammal assemblages in these regions have been reported to face the highest extinction risk globally (Pacifici et al. 2018), primarily due to simultaneous exposure of these assemblages to mid-21st century temperature and precipitation conditions beyond their realized niche limits (Trisos et al. 2020). The high intraspecific mitochondrial variation identified in these regions may serve as evolutionary reservoir for future adaptation to climate change, potentially reducing extinction risk, yet phylogenetic constraints and niche conservatism may counteract the beneficial effects of high genetic variability (Hoffmann and Sgrò 2011). The potential consequences of climate change impacts in GD may be even more relevant in those regions with disproportionally higher levels of mitochondrial variation, including central, eastern and southern Africa, and south Asia. These regions are highly exposed to global warming and may warrant further research and particular conservation attention.
We herein have explored the exposure of genetic diversity in mammals to climate change at global scale. The geography of species exposure to climate change is however detached from their sensitivity as the result of the spatial variation on species traits or physiological limits (Dickinson et al. 2014). Tropical mammal species are for example closer to experience climatic conditions beyond their thermal limits than high temperate species (Khaliq et al. 2014, Trisos et al. 2020, thus being potentially more vulnerable to extinction risk even under relatively low exposure to climate change. The ever-increasing availability of population level genetic data in conjunction with species ecological, life history and physiological data will pave the road to better estimations of species' sensitivity and vulnerability to climate change.

Exposure of genetic diversity to land-use change
Human impacts, especially conversion and degradation of natural habitats, have already caused severe biodiversity declines across all terrestrial regions (Newbold et al. 2015, Tilman et al. 2017, Pacifici et al. 2020. Primary vegetation in highly genetically diverse regions will decrease in the coming decades (average loss of 5-10% of their total area in tropical, subtropical and montane regions), although reduction of primary vegetation is also projected in low genetic diversity regions under all considered future scenarios. The loss of primary vegetation in genetically rich regions is coupled with increased human activities, primarily under the less sustainable future scenarios (SSP2-4.5, SSP3-7.0 and SSP5-8.5). These results suggest that the most evolutionary diverse terrestrial bioclimatic regions will be exposed to significant habitat loss, increasing species extinction risk. Tropical and subtropical regions have been indeed identified as the most exposed and vulnerable to 21st century human activities, such as cropland expansion and intensification (Kehoe et al. 2017), particularly for mammal diversity (Carvalho et al. 2019b, Powers andJetz 2019). Although high levels of standing GD in these regions, and to a great extent a higher level of GD than expected given their diversity of species, may provide adaptive potential to mammal assemblages in the face of future climate change, it is highly unlikely that they can mitigate future habitat loss and human exploitation, especially among large mammals. These direct human pressures will probably lead to severe reductions in intraspecific genetic variation in the tropical and montane regions through the extirpation of populations with unique haplotypes, significantly reducing the ability of mammal species to adapt to accelerating global change. Northern latitudes, and specifically boreal forests and taigas, which encompass ~30% of the global forest area, were also found to be significantly exposed to primary vegetation loss, without an identified increase in croplands. This pattern likely reflects a predicted increase in boreal deforestation due to logging (Gauthier et al. 2015), putting additional pressure, apart from climate change, to the already low mammal GD in this bioclimatic region.

Biases in CMIP6 models and spatial resolution limitations
Overall, the multi-model ensemble of the new CMIP6 GCMs showed similar spatial biases in annual land-surface temperature to the previous CMIP5 simulations (Flato et al. 2013). So far, both the previous and new generations of climate models show consistent biases in regions with steep topography (e.g. Himalayas and Andes), and in high latitude areas (e.g. Greenland and central East Siberia) where snow or ice cover influence the regional climate (Knutti et al. 2010, Flato et al. 2013. These biases may result in underestimation (in the case of negative bias) or overestimation (in the case of positive bias) of the exposure of GD to temperature change. Specifically, for southern Greenland and parts of the Himalayas and the Tibetan plateau, the multi-model ensembles predict a lower annual temperature by mid-21st century compared to the baseline period under all future scenarios. Yet, recent studies in these regions reported ongoing warming and extensive glacier retreats (Maurer et al. 2019, Richter-Menge et al. 2019, suggesting that our assessments may underestimate the exposure of GD to climate change in these regions. Nevertheless, the new CMIP6 GCMs show an overall increased sensitivity to emissions, simulating greater global warming over the 21st century (Forster et al. 2020, Zelinka et al. 2020, and future work in evaluating and constraining identified biases (Tokarska et al. 2020) will further improve our ability to assess biodiversity exposure to future climate change.
Given the relatively limited availability of genetic sequences with spatial information at global scale, and in order to minimize the uncertainty in the calculations of GD stemming from this limitation, a relatively coarse spatial resolution of grid-cells was favoured over finer resolutions. Additionally, the theories used in Theodoridis et al. (2020) to model the global distribution of GD, including evolutionary speed, the Red Queen and Late Quaternary climate stability hypotheses, invoke long-term metapopulation processes of persistence and extinction within relatively broad geographic areas. On the contrary, at finer spatial resolution, e.g. intrapopulation levels across elevational gradients, one would need to invoke a different set of theories and drivers more relevant to population or site-specific genetic diversity, including gene flow/inbreeding, linked selection and recent demographic history. Thus, denser sampling of wild populations globally is needed to model GD at finer scales and better inform conservation planning and management of the species genetic substrate. Population-specific genomic information is starting to unlock our potential to understand evolutionary responses of species and ecological assemblages to environmental change, including the ability to estimate relevant conservation metrics, such as the magnitude of mutation load in wild populations or the adaptive genetic variation revealed using landscape genomic analyses (Fitzpatrick andKeller 2015, Capblancq et al. 2020). While the coarse spatial grain used in our assessment cannot directly inform conservation strategies for GD at local population levels, it contributes towards understanding of global adaptive capacity and the potential evolutionary responses of mammal assemblages to climate and land-use change.

Conclusion
The key role of genetic diversity to safeguard biodiversity under global change was early recognized by the Convention on Biological Diversity (CBD 1992), and it is explicitly acknowledged in the Aichi 2020 strategic goals, Goal C and targets, Target 13 (SCBD 2010). Our assessment shows that a large proportion of the earth's mammal genetic diversity pool is predicted to experience a significant exposure to both climate and land-use changes. These results, coupled with the strong covariation between GD and phylogenetic diversity (Theodoridis et al. 2020), unveil the geographical variation of the exposure of mammal evolutionary history to global change, both at the intra-and interspecific levels. The inception of macrogenetics resulting from the availability of large-scale population-level genomic information (Blanchet et al. 2017), plus the accumulation of species trait information will soon make possible to link haplotypes and populations to more biologically meaningful climate and land-use change metrics (Garcia et al. 2014) for anticipating biodiversity sensitivity and vulnerability to global change.

Data availability statement
This study did not generate new data. The sources of the data used in this study are cited in the paper.