Introduction
Threats to biodiversity can be classified as ultimate or proximate. Ultimate threats include factors such as population growth and the expansion of global markets, the social, economic and political drivers of which are often beyond the direct control of conservation planners (Lambin et al. Reference Lambin, Turner, Geist, Agbola, Angelsen and Bruce2001). Proximate threats are defined as localized expressions of ultimate threats that directly affect biodiversity at regional or local scales (Pressey et al. Reference Pressey, Cabeza, Watts, Cowling and Wilson2007). Proximate threats are also dynamic, as they can increase or decrease in extent and intensity, operating as space–time mosaics of impact followed by recovery (Lambin et al. Reference Lambin, Turner, Geist, Agbola, Angelsen and Bruce2001, Pressey et al. Reference Pressey, Cabeza, Watts, Cowling and Wilson2007).
One dynamic proximate threat is land-use and land-cover change (LULCC), in which natural vegetation is transformed for alternative uses commonly related to the agricultural sector, introduced vegetation, forest plantations, deforestation and urban development (Velázquez et al. Reference Velázquez, Durán, Ramírez, Mas, Bocco, Ramírez and Palacio2003, Sullivan et al. Reference Sullivan, Binder and Burch2018). Moreover, conservation areas are open systems that could be affected by LULCC both within and outside their boundaries (Janzen Reference Janzen1983, Bengtsson et al. Reference Bengtsson, Angelstam, Elmqvist, Emanuelsson, Folke, Ihse and Nyström2003, Figueroa & Sánchez-Cordero Reference Figueroa and Sánchez-Cordero2008, Figueroa et al. Reference Figueroa, Illoldi-Rangel, Sánchez-Cordero and Linaje2011, Aguilar-Tomasini et al. Reference Aguilar-Tomasini, Escalante and Farfán2020). Currently in Mexico, there is no effective strategy to address dynamic threats such as LULCC. In particular, protected areas are not effective at stopping LULCC (Aguilar-Tomasini et al. Reference Aguilar-Tomasini, Escalante and Farfán2020) in the Trans-Mexican Volcanic Belt (TMVB) biogeographical province (Morrone Reference Morrone2014, Morrone et al. Reference Morrone, Escalante and Rodriguez-Tapia2017) in Central Mexico. As a consequence, LULCC models have been built to address this issue (Mendoza-Ponce et al. Reference Mendoza-Ponce, Corona-Núñez, Kraxner, Leduc and Patrizio2018, Reference Mendoza-Ponce, Corona-Núñez, Kraxner and Estrada2020).
A dynamic conceptual framework may be incorporated into systematic conservation planning (SCP) through LULCC models. SCP is a structured, stepwise protocol for prioritizing conservation area networks that also use surrogates as a representation of the whole biodiversity of an area (Margules & Sarkar Reference Margules and Sarkar2007). The selection of such surrogates is a critical component in conservation prioritization (Margules & Sarkar Reference Margules and Sarkar2007), and the process includes feedback, revision and reiteration of every step at all stages of the process (Margules & Pressey Reference Margules and Pressey2000, Margules & Sarkar Reference Margules and Sarkar2007). SCP also recognizes that tangible results and costs can be associated with particular actions, such as prioritizing certain land-use practices outside protected areas, which could be an efficient complementary strategy to reduce biodiversity loss (Tallis et al. Reference Tallis, Fargione, Game, McDonald, Baumgarten, Bhagabati and Possingham2021, Staden et al. Reference Staden, Lötter, Holness and Lombard2022, Giakoumi et al. Reference Giakoumi, Richardson, Doxa, Moro, Andrello, Hanson and Katsanevakis2025).
Most of the strategies involved in the design of protected areas are treated primarily as static entities, which limits their capacity to respond to rapid environmental changes (Harris et al. Reference Harris, Hobbs, Higgs and Aronson2006). Area-based dynamic management has been proposed to address this challenge (Bengtsson et al. Reference Bengtsson, Angelstam, Elmqvist, Emanuelsson, Folke, Ihse and Nyström2003). Dynamic conservation areas are temporarily protected places that are specifically designed to protect the increasingly dynamic nature of landscapes, as well as population, community and ecosystem-level processes (D’Aloia et al. Reference D’Aloia, Naujokaitis-Lewis, Blackford, Chu, Curtis, Darling and Fortin2019). This approach involves temporarily protecting certain areas from anthropogenic pressures, followed by the withdrawal of formal protection when it is no longer needed (Wen et al. Reference Wen, Bennett, Rytwinski, Karimi, Spetka, Taylor and Smith2024). Although dynamic conservation areas have gained greater acceptance in marine environments, they are also applicable to terrestrial and freshwater systems (D’Aloia et al. Reference D’Aloia, Naujokaitis-Lewis, Blackford, Chu, Curtis, Darling and Fortin2019).
Dynamic areas can be especially useful tools in fragmented landscapes with high human populations, where high-conservation-value areas are limited and change over time (Reynolds et al. Reference Reynolds, Sullivan, Hallstein, Matsumoto, Kelling, Merrifield and Morrison2017). They are also particularly well suited for species with migratory behaviours or other specific ecological traits (Reynolds et al. Reference Reynolds, Sullivan, Hallstein, Matsumoto, Kelling, Merrifield and Morrison2017). The inclusion of dynamic conservation areas within permanent networks could provide temporary protection to mitigate LULCC and biodiversity redistributions (D’Aloia et al. Reference D’Aloia, Naujokaitis-Lewis, Blackford, Chu, Curtis, Darling and Fortin2019). The eastern sub-province of the TMVB in Central Mexico is a crucial location for a dynamic strategy because it harbours endemic and rare surrogate species (Escalante et al. Reference Escalante, Varela-Anaya, Noguera-Urbano, Elguea-Manrique, Ochoa-Ochoa and Gutiérrez-Velázquez2020) while containing some of the most populous cities in the country, including Mexico City, Cuernavaca, Puebla and Toluca, with high LULCC (Mendoza-Ponce et al. Reference Mendoza-Ponce, Estrada, Campo, Calderón, Del Carmen and Ramírez2025).
In this context, we aim to identify the potential areas that may facilitate the conservation of rare and endemic biodiversity in the face of ongoing LULCC in the eastern sub-province of the TMVB under a dynamic conceptual framework and using mammals and amphibians as surrogates.
Materials and methods
Study area
The eastern sub-province of the TMVB covers c. 62 848 km2 and harbours some of the highest mountains in Mexico (Fig. 1; Morrone Reference Morrone2020). For operational reasons, we delimited the eastern sub-province of the TMVB (Fig. 1) based on the polygon of Morrone et al. (Reference Morrone, Escalante and Rodriguez-Tapia2017), refined according to the biogeographical regionalization of Gámez et al. (Reference Gámez, Escalante, Rodríguez, Linaje and Morrone2012) and the updated descriptions in Morrone (Reference Morrone2020).
Study area showing the eastern sub-province of the Trans-Mexican Volcanic Belt.
Geographical distribution areas and species distribution models
We selected 13 endemic species (seven amphibians and six mammals) from the eastern sub-province of the TMVB that were previously identified as good biodiversity surrogates by Escalante et al. (Reference Escalante, Varela-Anaya, Noguera-Urbano, Elguea-Manrique, Ochoa-Ochoa and Gutiérrez-Velázquez2020) using cross-taxon congruency based on their complementary richness and rarity patterns. Scientific names were verified for accuracy and taxonomic validity using the Amphibian Species of the World database (Frost Reference Frost2024) and the Mammal Diversity Database (2024).
Occurrence data for the 13 species were obtained from the database of Flores-Villela and Ochoa-Ochoa (Reference Flores-Villela and Ochoa-Ochoa2024) and the Global Biodiversity Information Facility (GBIF; www.gbif.org). The occurrences were reviewed for errors in R (R Core Team 2023) using the ‘Coordinate Cleaner’ R package (Zizka et al. Reference Zizka, Silvestro, Andermann, Azevedo, Duarte Ritter and Edler2019), removing empty fields, duplicates and invalid coordinates. For data downloaded from GBIF, we removed occurrences with spatial uncertainty greater than 1000 m. We also identified and excluded outliers and filtered the occurrences based on expert range polygons from the International Union for Conservation of Nature (IUCN 2024).
Because there were species with <10 records, which is fewer than the minimum number of occurrence records required for reliable species distribution models (SDMs), and in order to obtain a comparable inference of a geographical distribution area (GDA), we used those expert-derived distribution maps from the IUCN (2024), which were rasterized at c. 1-km2 resolution. Because SDMs built with very small sample sizes can produce unstable predictions, high uncertainty and overfitted responses that poorly represent the environmental niche of the species (Pearson et al. Reference Pearson, Raxworthy, Nakamura and Townsend Peterson2007, Shcheglovitova & Anderson Reference Shcheglovitova and Anderson2013), we used expert-based range maps as a more robust alternative for poorly sampled or geographically rare taxa. This allowed the inclusion of highly vulnerable species despite their limited occurrence data in conservation prioritization analyses.
We used the ‘Wallace2’ package (Kass et al. Reference Kass, Pinilla-Buitrago, Paz, Johnson, Grisales-Betancur, Meenan and Anderson2023) to apply the maximum entropy algorithm to build the SDM. The 19 bioclimatic variables from CHELSA (Climatologies at High Resolution for the Earth’s Land Surface Areas) and the three topographical variables of altitude, aspect and slope (United States Geological Survey 2021) were used at c. 1-km2 resolution. To produce the SDMs, we defined the calibration areas for each species using QGIS (http://qgis.org) through a minimum convex hull and applied a buffer based on the average distance between occurrences. We then filtered the data to retain only one occurrence within a 1-km radius to reduce sampling bias and over-representation of environmental conditions (Aiello-Lammens et al. Reference Aiello-Lammens, Boria, Radosavljevic, Vilela and Anderson2015). This ensured spatial independence among the records (Pearson et al. Reference Pearson, Raxworthy, Nakamura and Townsend Peterson2007).
To reduce potential sampling bias and the overrepresentation of environmental conditions in some areas, the cleaned occurrence records were spatially filtered using a minimum inter-record distance of 1 km (Aiello-Lammens et al. Reference Aiello-Lammens, Boria, Radosavljevic, Vilela and Anderson2015). This threshold ensured sufficient environmental heterogeneity among the records, allowing them to be treated as spatially independent observations (Pearson et al. Reference Pearson, Raxworthy, Nakamura and Townsend Peterson2007). To reduce multicollinearity among predictors, we calculated the variance inflation factor (VIF) and removed variables with VIF > 10 (Montgomery & Peck Reference Montgomery and Peck1992, Mandeville Reference Mandeville2008) using the ‘usdm’ package in RStudio (Naimi et al. Reference Naimi, NA, Groen, Skidmore and Toxopeus2014). In the ‘Wallace2’ package (Kass et al. Reference Kass, Pinilla-Buitrago, Paz, Johnson, Grisales-Betancur, Meenan and Anderson2023), we tested the regularization multipliers 0.5, 1.0, 1.5, 2.0, 2.5, 3.0, 3.5 and 4.0 and the feature classes linear (L), linear-quadratic (LQ) and hinge (H).
To partition the training and validation data, we used the jackknife method for species with 10–23 occurrence records. This method is a special case of k-fold cross-validation and is recommended for datasets with small sample sizes. For species with ≥24 occurrence records, we used the random k-fold method, which creates four random groups of equal size and trains the model using k – 1 of those subsets, validating it with the remaining subset not used in training (Muscarella et al. Reference Muscarella, Galante, Soley-Guardia, Boria, Kass, Uriarte and Anderson2014).
For species with 10–39 occurrences, the best model was selected based on the highest area under the curve (AUC) value and the lowest omission rate (OR). For species with ≥40 occurrences, the four models with the best AUC and OR values were evaluated using the partial receiver operating characteristic (ROCp) test implemented in the NicheToolbox platform (Osorio-Olvera et al. Reference Osorio-Olvera, Lira-Noriega, Soberón, Townsend-Peterson, Falconi and Contreras-Díaz2020), and the model with the highest ROCp value was selected as the best. The parameters used for the ROCp test were an omission proportion of 0.1, 50% random points and 500 bootstrap iterations.
Finally, binary models were generated (0 = unsuitable habitat, 1 = suitable habitat) using the quantile training presence thresholding method, setting the first quartile as the cutoff point.
Land-use and land-cover change models
To identify priority areas for conservation, we used land use and land cover (LULC) maps from 1985 to 2018 developed by the Instituto Nacional de Estadística y Geografía (INEGI) and the LULCC map models of Mendoza-Ponce et al. (Reference Mendoza-Ponce, Corona-Núñez, Kraxner, Leduc and Patrizio2018, Reference Mendoza-Ponce, Corona-Núñez, Kraxner and Estrada2020). These models were constructed under three scenarios: a ‘Business-As-Usual’ (BAU) scenario combining SSP2 (Shared Socioeconomic Pathways) and RCP4.5 (Representative Concentration Pathways); a ‘Green’ scenario combining SSP1 and the RCP2.6; and a ‘Worst’ scenario, combining SSP3 and the RCP8.5. Each scenario was modelled using four general circulation models (GCMs) for the intermediate time horizon 2045–2069: CNRMC-M5, GFDL-CM3, HADGEM2-ES and MPI-ESM-LR. This resulted in a total of 12 future (2050) LULCC projections at c. 1-km2 resolution.
The selected GCMs were originally chosen by Mendoza-Ponce et al. (Reference Mendoza-Ponce, Corona-Núñez, Kraxner, Leduc and Patrizio2018, Reference Mendoza-Ponce, Corona-Núñez, Kraxner and Estrada2020) to incorporate variability among contrasting climate projections for Mexico and to maintain consistency with national climate change assessments developed by the Instituto Nacional de Ecología y Cambio Climático (INECC). The use of multiple GCMs allowed the incorporation of inter-model climatic variability into the LULCC projections and reduced dependence on a single climate trajectory.
To evaluate uncertainty among GCMs, we performed spatial agreement analyses among the four GCM-derived LULCC projections within each SSP–RCP scenario (BAU, Green and Worst). Agreement was assessed pixel-wise by calculating the proportion of GCMs sharing the dominant projected land-use transition for each raster cell, resulting in agreement values of 25%, 50%, 75% or 100%, corresponding to agreement among one, two, three or four models, respectively. Areas with agreement values ≥75% were interpreted as spatially consistent projected transitions among at least three of the four GCMs. Spatial agreement analyses revealed high consistency among GCMs, with areas showing agreement of ≥75% representing 87.66%, 88.08% and 88.59% of the study area under the BAU, Green and Worst scenarios, respectively, indicating that the main projected patterns of land-use change were robust across climate models. Although multiple GCM-derived projections were available, we selected the HADGEM2-ES-derived BAU, Green and Worst scenarios for subsequent analyses to maintain consistency and comparability with previous studies evaluating climate change impacts on biodiversity and vegetation dynamics in Mexico (e.g., Manzanilla-Quiñones et al. Reference Manzanilla-Quiñones, Aguirre-Calderón, Jiménez-Pérez, Treviño-Garza and Yerena-Yamallel2019, Prieto-Torres et al. Reference Prieto-Torres, Nori, Rojas-Soto and Navarro-Sigüenza2021, Armenta-Montero et al. Reference Armenta-Montero, Menchaca-García, Pérez-Silva and Velázquez-Rosas2022, Estrada et al. Reference Estrada, Mendoza, Murray, Calderón-Bustamante, Botzen, De León Escobedo and Velasco2023).
Evaluation of the effect of LULCC on richness and rarity patterns
We measured the endemic species richness and endemic rarity (also known as micro-areality, referring to species with the most restricted distributional ranges) using the ‘letsR’ package (Vilela & Villalobos Reference Vilela and Villalobos2015). We calculated the total endemic species richness by counting the number of species in each pixel within the eastern sub-province of the TMVB.
Rare species were identified using the ‘Rarity’ package (Leroy Reference Leroy2023), which is based on a scoring procedure developed by Leroy et al. (Reference Leroy, Petillon, Gallon, Canard and Ysnel2012, Reference Leroy, Canard and Ysnel2013). This method includes the selection of a rarity cutoff, expressed as a percentage of the maximum occurrence (Leroy et al. Reference Leroy, Petillon, Gallon, Canard and Ysnel2012). Here, we selected 25% of species with the lowest occurrences as the rarity cutoff point (Gaston Reference Gaston1994). Then, we quantified the richness of rare species in each pixel and made maps of them as raster files.
Finally, we overlaid the geographical distribution patterns of both metrics with LULCC maps, including the historical scenario (1985–2018) and projections for 2050 using QGIS. We also calculated the percentage of loss of natural vegetation, translated here as loss of habitat quality for each species due to anthropogenic factors.
Identification of priority conservation areas
To identify priority conservation areas, we used Zonation 4.0.0 (Moilanen et al. Reference Moilanen, Pouzols, Meller, Veach, Arponen, Leppänen and Kujala2014), which is commonly used to prioritize geographical areas based on their biodiversity importance and other ecological and socioeconomic values. Zonation performs a hierarchical spatial prioritization by progressively removing landscape cells that cause the smallest loss in conservation value (Moilanen et al. Reference Moilanen, Pouzols, Meller, Veach, Arponen, Leppänen and Kujala2014). The result is a continuous ranking of conservation priority across the entire landscape.
We used the following settings to generate the conservation model: first, we selected the core-area zonation (CAZ) removal rule, which prioritizes even generally feature-poor cells that can be identified as priorities (Di Minin et al. Reference Di Minin, Veach, Lehtomäki, Pouzols and Moilanen2014). Thus, even small contributions to species with highly restricted ranges strongly influence the ranking; this approach emphasizes the protection of geographically rare species, the distributions of which are more vulnerable to habitat loss due to their specialized ecological and habitat requirements.
To account for spatial connectivity within the prioritization process, we implemented Zonation’s built-in spatial aggregation settings, which promote structurally connected conservation networks during cell removal and ranking. Specifically, we activated the edge removal rule to maintain spatial cohesion and avoid premature fragmentation of priority areas. In addition, we applied a boundary length penalty (BLP = 0.5), which discourages dispersed or highly irregular spatial configurations by favouring compact and contiguous priority areas. We also set the boundary quality penalty (BQP = 1) to enhance the selection of neighbouring cells with a similar conservation value. Together, these parameters incorporated structural connectivity directly into the prioritization procedure, promoting coherent conservation areas (Di Minin et al. Reference Di Minin, Veach, Lehtomäki, Pouzols and Moilanen2014). Existing Protected Natural Areas (PNAs; CONANP 2024) within the study region were not masked in the analysis. This decision was based on the findings of Aguilar-Tomasini et al. (Reference Aguilar-Tomasini, Escalante and Farfán2020), who assessed the effectiveness of PNAs in the TMVB with respect to deforestation and LULCC. Their results indicated that only three PNAs are currently effective, and only two of them (Ciénegas de Lerma and Iztaccíhuatl-Popocatépetl) are located within the eastern sub-province of the TMVB. As such, masking these areas was not expected to influence the prioritization outcomes.
For the prioritization, we included the binary species maps, which were weighted according to their conservation status on the IUCN Red List (IUCN 2024): Critically Endangered (CR) = 1; Endangered (EN) = 0.75; Vulnerable (VU) = 0.5; Least Concern (LC) = 0.25. We also incorporated LULCC data from 2018 developed by INEGI (2021), as well as a future LULCC scenario (BAU, because it represents the most plausible future trajectory), a layer of communal lands obtained from the Registro Agrario Nacional (RAN 2023), a layer of ecosystem service valuation derived from the Co$ting Nature platform (www.policysupport.org/costingnature) and a layer of human population density from Schiavina et al. (Reference Schiavina, Freire, Carioli and MacManus2023; see Appendices S1 & S2 for explanations of these layers).
Finally, we established a conservation target of 30% of the study area to align with international biodiversity commitments (Convenio sobre la Diversidad Biológica 2022). Based on this threshold, we retained only the upper 30% of the landscape identified as having the highest conservation priority.
Results
Geographical distribution areas and species distribution models
We obtained GDA rastering from the expert range maps from the IUCN for four species: Chiropterotriton orculus, Neotoma nelsoni, Thorius munificus and Xerospermophilus perotensis. SDMs were obtained for nine other species: Chiropterotriton lavae, Cratogeomys merriami, Dryophytes plicatus, Neotomodon alstoni, Peromyscus bullatus, Pseudoeurycea altamontana, Pseudoeurycea gadovii, Pseudoeurycea robertsi and Romerolagus diazi (Appendix S1 Table S1).
The total endemic richness pattern revealed higher species concentrations on some volcanic mountain systems such as the Sierra del Ajusco-Chichinautzin (located in the states of Morelos and Mexico), Tepozteco (Morelos), Iztaccíhuatl-Popocatépetl (Mexico, Puebla and Morelos) and Cofre de Perote (Veracruz; Fig. 2a). The rare endemic species identified were C. lavae, N nelsoni, P. robertsi and T. munificus, and their richness patterns also overlapped with volcanic mountain systems, including Nevado de Toluca (Mexico State) and Cofre de Perote (Veracruz; (Fig. 2b).
Richness patterns in the eastern sub-province of the Trans-Mexican Volcanic Belt in Central Mexico: (a) total richness (number of endemic species per pixel); (b) richness of rare endemic species (number of rare endemic species per pixel).
Land-cover and land-use change models
Rainfed agriculture has been the dominant land-use class in the study area since 1985, covering c. 40% of the territory. Under the Worst scenario, this cover type could increase to nearly 70% of the study area. Together with pastureland for livestock, which could represent an additional 15% of the area, anthropogenic land covers would account for c. 85% of the projected landscape. In contrast, temperate forests are projected to represent only 1% of the study area under this scenario (Fig. 3).
Percentage of land-cover types in the eastern sub-province of the Trans-Mexican Volcanic Belt for different years and future land-use change scenarios: 2050 Green = optimistic scenario with reduced greenhouse gas emissions and deforestation; 2050 BAU = a scenario that follows the current trends of greenhouse gas emissions and deforestation (business as usual); 2050 Worst = pessimistic scenario with accelerated increase of greenhouse gas emissions, deforestation and environmental degradation.
Figure 3. Long description
A stacked bar graph comparing the percentage of land-cover types over various years and future scenarios. The horizontal axis represents the years and scenarios: 1985, 1993, 2002, 2007, 2011, 2015, 2018, 2050 Green, 2050 BAU, and 2050 Worst. The vertical axis represents the percentage of land-cover types, ranging from 0 percent to 100 percent. Each bar is divided into segments representing different land-cover types, with each segment color-coded. The colors represent various land-cover types: Temperate forests, Tropical dry forests, Rainfed agriculture, Irrigated agriculture, Cloud forests, Hydrophilic vegetation, Natural grasslands, Other kind of natural vegetation, Shrubland, Pastures for livestock, Human settlements (urban), Water bodies, and Barren land. The graph shows changes in land-cover types over time and under different future scenarios, with notable increases or decreases in specific land-cover types in the future scenarios.
LULCC impacts on the species richness and rarity patterns
In the 2018 scenario, habitat degradation affected c. 65% of the pattern of total richness of endemic species and c. 43% of the richness of rare endemic species. The degradation decreased under the Green scenario to 34% for total richness and to 18% for rare species. Degradation increased considerably under the BAU scenario, reaching up to 75% for total richness and up to 64% for rare species. Maximum degradation was projected under the Worst scenario, reaching 98% for total richness and 97% for rare species (see Appendix S3 for a graphical representation).
Priority conservation areas
In the first-priority rank map, we used the 2018 LULCC map (Fig. 4a) to identify priority areas; those were concentrated within the remaining native vegetation patches of the study area, with several of them associated with volcanic mountain systems such as the Sierra del Ajusco-Chichinautzin (Morelos, Mexico City and Mexico State), El Tepozteco (Morelos), Iztaccíhuatl-Popocatépetl (Mexico State, Puebla and Morelos), La Malinche (Tlaxcala and Puebla) and Cofre de Perote (Veracruz; see Fig. 1). The second-priority rank map was produced using the 2050 LULCC map under the BAU scenario (Fig. 4b). This map highlighted the aforementioned mountain systems, along with smaller natural patches surrounding them and areas projected to undergo agricultural expansion. Compared to the 2018 scenario, priority areas under this projection appeared more spatially dispersed and less continuous across the landscape. Some priority areas persisted within heterogeneous anthropogenic landscapes, which may be consistent with reports of certain taxa, such as C. merriami and D. plicatus, occurring in disturbed environments (Marquez-Villalba & Martínez-Coronel Reference Marquez-Villalba and Martínez-Coronel2024, Ramírez et al. Reference Ramírez, Marín, Mendoza, Bautista, Zamora and Jiménez2025).
Priority rank map in the eastern sub-province of the Trans-Mexican Volcanic Belt. (a) Prioritization using the 2018 land-use and land-cover change (LULCC) scenario. (c) Zonation performance curve of (a), showing the relationship between the loss of priority rank areas across the landscape and the corresponding decrease in species distribution coverage; the dotted grey line represents the 30% conservation target and the corresponding average proportion of species’ ranges protected. (b) Prioritization using the 2050 Business-As-Usual (BAU) LULCC scenario. (d) Zonation performance curve of (b), showing the relationship between the loss of priority rank areas across the landscape and the corresponding decrease in species distribution coverage; the dotted grey line represents the 30% conservation target and the corresponding average proportion of species’ ranges protected for the future.
Performance curves summarized the conservation coverage achieved within each top-priority fraction derived from the priority rank maps. In particular, mean performance curves provided an overall assessment of the prioritization outcomes. When considering average performance, the top 30% of prioritized areas captured c. 65% of the mean surrogate species distributions in the 2018 scenario and c. 70% for the 2050 BAU scenario, representing the smallest area required for protection (Fig. 4c,d).
Discussion
Land-use transformation has been widely documented as one of the main drivers of biodiversity loss in Mexico, particularly in the TMVB (Arriaga et al. Reference Arriaga, Espinoza, Aguilar, Martínez, Gómez and Loa2000, Suárez-Mota & Téllez-Valdés Reference Suárez-Mota and Téllez-Valdés2014, Mendoza-Ponce et al. Reference Mendoza-Ponce, Corona-Núñez, Galicia and Kraxner2019, Bonilla-Moheno & Aide Reference Bonilla-Moheno and Aide2020). Our results are consistent with this pattern and indicate that anthropogenic land covers already dominate much of the eastern sub-province of the TMVB. By 2018, c. 43% of the habitat associated with geographically rare species and c. 65% of the habitat associated with species richness had already been degraded. These values suggest that land-use change is not only reducing the extent of natural habitats but is disproportionately affecting areas that concentrate biodiversity value.
Future projections indicate that the magnitude of this impact depends strongly on the development pathway considered. Under the Green scenario, degraded habitat decreases substantially (18% for rarity and 34% for richness), suggesting that more sustainable land-use trajectories could partially alleviate biodiversity loss. In contrast, the BAU and Worst scenarios project extensive habitat degradation, reaching 64% for rarity and 75% for richness in the BAU scenario and 97% for rarity and 98% for richness in the Worst scenario. These trends coincide with the projected expansion of rainfed agriculture, which could increase from c. 40% of the landscape in 2018 to nearly 70% under the Worst scenario, while temperate forests decline to c. 1% of the study area. Together, these results suggest that agricultural expansion is likely to become the dominant driver of habitat loss across the region.
Additionally, previous studies have shown that endemic species are particularly vulnerable to land-use change in Mexico (Escalante et al. Reference Escalante, Sánchez-Cordero, Morrone and Linaje2007, Sánchez-Cordero et al. Reference Sánchez-Cordero, Illoldi-Rangel, Escalante, Figueroa, Rodríguez, Linaje, Columbus and Kuznetsov2009). For example, Aguilar-Ortega and Escalante (Reference Aguilar-Ortega and Escalante2021) reported that nearly 10% of the endemic species evaluated in the TMVB had their entire distribution affected by land-use transformation. Our results support and extend these findings by showing that both endemic and geographically rare species are expected to experience substantial reductions in suitable habitat under future land-use scenarios. This pattern reinforces the importance of implementing conservation actions before projected landscape changes become irreversible.
On the other hand, Aguilar-Tomasini et al. (Reference Aguilar-Tomasini, Escalante and Farfán2020) demonstrated that the current protected area network in the TMVB has been insufficient to halt land-use change across the region. Several studies have identified priority areas for conservation using different approaches, including panbiogeographical analyses (Torres-Miranda & Luna-Vega Reference Torres-Miranda and Luna-Vega2006) and SCP (Fuller et al. Reference Fuller, Munguía, Mayfield, Sánchez-Cordero and Sarkar2006, Morales et al. Reference Morales, Escalante, Noguera-Urbano, Gámez and Rodríguez-Tapia2016, Juárez-Ramírez et al. Reference Juárez-Ramírez, Rodriguez-Soto, Estrada-Torres, Vazquez, de la Vega and Arias-Del Razo2025). Despite differences in datasets, conservation targets and analytical methods, these studies consistently highlight a common set of mountain systems as conservation priorities, including Iztaccíhuatl-Popocatépetl, Sierra de Tlaxco, Nevado de Toluca, Pinal Mountain, Cofre de Perote and Pico de Orizaba. The recurrence of these areas across independent analyses reinforces their ecological importance and suggests that they represent persistent biodiversity strongholds within the TMVB.
Our prioritization results are broadly consistent with these previous studies, but they also identify additional priority areas. Unlike most previous studies, our analysis integrated not only species distributions and land-use patterns but also socioeconomic variables, including land tenure, population density and ecosystem service valuation. The inclusion of these variables shifted conservation priorities towards areas where biological importance coincides with greater implementation feasibility. As a result, several priority areas emerged outside the main volcanic systems, potentially increasing opportunities for landscape-scale conservation planning.
Comparison between the current and future prioritizations revealed notable changes in the spatial distribution of conservation priorities. The principal mountain systems remained among the highest-ranked areas under both scenarios, indicating that their conservation importance is robust to projected land-use changes. However, under the 2050 BAU scenario, priority areas became more spatially dispersed across the landscape and increasingly associated with smaller remnants of natural vegetation embedded within anthropogenic matrices. This pattern is consistent with the projected expansion of agricultural land and the reduction of natural habitats described earlier. Although we did not explicitly quantify fragmentation metrics, the spatial redistribution of high-priority cells suggests that future conservation actions may need to operate across a more heterogeneous landscape than at present.
The prioritization curves provide additional insights into these patterns. Under the 2018 scenario, c. 60% of the average species representation was achieved by protecting 30% of the landscape, whereas under the BAU scenario the same proportion of protected area represented c. 70% of species distributions. This result suggests that, in future, suitable habitats may become increasingly concentrated within a smaller set of locations, thereby increasing the conservation value of the remaining natural areas. At the same time, this concentration may also increase the vulnerability of biodiversity to additional habitat loss if these areas are not effectively conserved.
Previous studies have suggested that conservation strategies capable of responding to changing environmental conditions can be valuable in landscapes experiencing rapid ecological transformation. For example, McComb et al. (Reference McComb, Powers and Larsen2022) proposed a temporally dynamic conservation framework in California that improved habitat connectivity for Vulpes macrotis through the opportunistic use of agricultural lands. Similarly, Edwards et al. (Reference Edwards, Lisgo, Leroux, Krawchuk, Cumming and Schmiegelow2022) developed dynamic reserve designs to address changing disturbance regimes in boreal forests. Our results support the exploration of dynamic conservation approaches within the TMVB, providing a spatially explicit basis for identifying where such approaches could be explored in future conservation planning efforts.
Proposal of a conceptual framework for biodiversity conservation
We propose a conservation network with the intention to integrate it with the current network of PNAs, including a new set of nodes and corridors under a dynamic conceptual framework (Fig. 5). In that sense, the current PNAs that overlap with our prioritization could act as permanent nodes because they provide long-term legal stability, and, as we can see in our results, they contain persistent natural cover, harbour multiple species of conservation interest and act as large ecological cores with established governance and monitoring systems (Fig. 5a). The main PNAs include Cerro de Garnica (Michoacan), Mariposa Monarca (Michoacan and Mexico State), Nevado de Toluca (Mexico State), Corredor Biológico Chichinautzin (Mexico City, Morelos and Mexico State), Cumbres del Ajusco (Mexico City), Iztaccíhuatl-Popocatépetl (Mexico State, Puebla and Morelos), La Malinche (Tlaxcala and Puebla), Pico de Orizaba (Veracruz and Puebla) and Cofre de Perote (Veracruz).
Proposed dynamic conservation network. (a) Permanent nodes (black polygons: 1 = Cofre de Perote, 2 = Corredor Biológico Chichinautzin, 3 = Iztaccíhuatl-Popocatépetl, 4 = Nevado de Toluca, 5 = Cerro de Garnica, 6 = Mariposa Monarca, 7 = Cumbres del Ajusco, 8 = La Malinche, 9 = Pico de Orizaba) and type I dynamic nodes (purple polygons). (b) Type II (green polygons) and type III (blue polygons) dynamic nodes and corridors (pink polygons).
Figure 5. Long description
Panel A: A map highlighting permanent conservation nodes in black polygons and type I dynamic nodes in purple polygons across central Mexico. Key locations include Cofre de Perote, Corredor Biológico Chichinautzin, Iztaccíhuatl-Popocatépetl, Nevado de Toluca, Cerro de Garnica, Mariposa Monarca, Cumbres del Ajusco, La Malinche, and Pico de Orizaba. Panel B: A map showing type II dynamic nodes in green polygons, type III dynamic nodes in blue polygons, and corridors in pink polygons. The maps cover states such as Querétaro, Guanajuato, Hidalgo, México, Morelos, Guerrero, and Puebla. Priority areas for conservation are indicated with a color gradient from green to red, representing increasing priority.
Regarding the dynamic nodes, we propose three types (Fig. 5a). The first type is represented by priority areas that expand the functional range of the permanent nodes (surrounded by purple polygons in Fig. 5a). Their location around the current PNAs is important because they overlap with extensive natural cover in the foothills of the mountains, potentially acting as broader refuges for multiple species. These nodes could be withdrawn from formal protection once the species currently listed as CE by the IUCN (e.g., R. diazi, X. perotensis) move to a safer conservation category (e.g., LC). After the species’ status improves, these dynamic conservation areas could be reactivated only during key ecological periods, such as the rainy season for amphibians, which are critical for reproduction.
The second type of dynamic node corresponds to priority areas located within highly fragmented landscapes dominated by agricultural activities (surrounded by green polygons in Fig. 5b). One of these areas (the northern node) overlaps with the distribution of geographically rare species such as C. lavae, N. nelsoni and T. munificus. Because these species are categorized as CE, protecting this node is crucial, and habitat restoration should be a key component of its management. We therefore propose a strict conservation scheme for this northern node until the species move to a safer conservation category; after that, the scheme could become more flexible. In contrast, the second node within this group may follow a more flexible management framework, where ecological restoration would remain a priority but could be combined with sustainable agricultural practices.
The third of the dynamic node types corresponds to areas that are not highly anthropized and that can be temporarily conserved during ecologically important periods for the species (blue polygons in Fig. 5b). Those nodes also offer opportunities to develop sustainable management plans for the resources they provide.
Finally, patches of natural land cover could act as ecological corridors connecting the network (pink polygons in Fig. 5b). These corridors could operate dynamically, being activated during periods when species movement is most intense and relaxed during less critical periods. This temporal activation would maintain connectivity while enabling sustainable resource-use practices when ecological sensitivity is lower.
We propose a spatial–temporal and politically dynamic conceptual framework (Fig. 5a,b). Our approach is spatially dynamic, as it allows the expansion of some NPAs that still retain natural vegetation or high-quality habitat in their surroundings, thereby strengthening ecological integrity. It is also temporally dynamic, incorporating the protection of specific sites during ecologically important seasons for the species, when their habitat requirements are at their peak. Finally, it is politically dynamic, enabling greater flexibility in management once species show improvements in their conservation status, thus allowing the sustainable use of resources without compromising long-term biodiversity outcomes.
Limitations and future directions
Despite the insights provided by our analyses, several limitations should be acknowledged. First, future land-use projections are inherently uncertain because they depend on assumptions regarding climate, socioeconomic trajectories and policy decisions. Although the agreement analysis among GCMs revealed high spatial consistency, future landscape trajectories may still differ from those projected here. Second, the prioritization analyses were based on modelled species distributions and landscape conditions and therefore do not explicitly incorporate demographic processes, dispersal dynamics, species interactions or evolutionary responses. We acknowledge that this approach does not explicitly represent functional connectivity or species-specific movement processes, which could be explored in future studies using resistance surfaces or graph-theoretical approaches. Third, the dynamic conservation framework proposed in this study should be interpreted as a conceptual framework rather than a validated management strategy. Its effectiveness remains to be tested through empirical applications, stakeholder participation processes and long-term monitoring programmes. Implementation may also face practical challenges associated with governance structures, funding availability and social acceptance, all of which can strongly influence conservation outcomes in human-dominated landscapes.
Future research should evaluate the performance of dynamic conservation strategies under alternative socioeconomic and climate scenarios, incorporate connectivity and explore participatory approaches involving local communities and decision-makers (e.g., through a multi-criteria analysis). Such efforts would help determine the feasibility and effectiveness of translating conceptual dynamic conservation frameworks into operational conservation actions.
Conclusions
Endemic and rare species are particularly affected by LULCC in the eastern TMVB; nearly half of the potential habitat for those species has already been degraded. The BAU projections for the year 2050 indicated an even more fragmented landscape, emphasizing the urgency of implementing more effective and forward-looking conservation measures.
Although various studies have used diverse datasets and analytical approaches to find priority areas for conservation, our findings converge with previous research in identifying a recurrent group of mountains and surrounding areas as essential conservation nodes. This consistency strongly suggests that these areas hold high ecological value. By incorporating evolutionary (through the patterns of endemism and rarity), ecological (distributional patterns), land-use and socioeconomic variables (land tenure, human population density and economic valuation of ecosystem services), our framework also identifies additional priority areas and potential corridors that may improve connectivity and social feasibility.
Based on our analyses, we propose a dynamic conceptual framework with spatial, temporal and political flexibility. This approach allows NPAs to be adjusted according to habitat quality, protects key sites during ecologically important periods and enables adaptable management when species or local conditions change.
Supplementary material
The supplementary material for this article can be found at https://doi.org/10.1017/S0376892926100587.
Acknowledgements
We thank Leopoldo Galicia Sarmiento from IG-UNAM and Julián A Velasco Vinasco from ICACC-UNAM for their valuable feedback. We are also grateful to Lynna M Kiere for her assistance in improving the English language, style and grammar of this manuscript.
Financial support
Oscar Campos is particularly grateful to the Secretaría de Ciencias, Humanidades, Tecnología e Innovación (SECIHTI) for funding and support of this research through a graduate scholarship for his MSc (CVU: 1281678). We acknowledge the National Autonomous University of Mexico (UNAM) for the economic support provided through the Research and Technological Innovation Projects Support Program (PAPIIT) for the project Planeación Sistemática de la Conservación en Querétaro: incendios forestales y Biodiversidad (IN221225).
Competing interests
The authors declare that there are no known competing financial interests or personal relationships that could have appeared to influence the work reported.
Ethical standards
Not applicable.
Open access
