To gather a comprehensive dataset of the plant’s distribution within its range, we conducted exhaustive literature reviews utilizing online databases, incorporating references sourced from the Global Biodiversity Information Facility (GBIF, https://www.gbif.org/) and the Chinese Virtual Herbarium (CVH, http://v5.cvh.org.cn/). In cases where distribution records lacked complete geographic coordinates, we employed Google Earth¹ to ascertain precise latitudes and longitudes.

Fieldwork was conducted from July 2022 to September 2023. It covered potential distribution zones such as the Qinling-Bashan Mountains and Daba Mountains (Supplementary Figure S1). Survey methods included setting up systematic transects in different ecological habitats, making direct observations to record individual plants, using high-precision GPS for precise localization, and, in some cases, collecting specimens for later morphological identification. To ensure data quality, duplicate records were eliminated, and spatial filtering using ENM Tools was applied to guarantee one locality per 1 km² grid cell. This resulted in a final dataset comprising 77 unique distribution points. ArcGIS 10.3 software was then utilized to create a distribution map based on these data, providing valuable insights for subsequent research endeavors and confirming the accuracy and reliability of the distribution dataset (Figure 2).

To fully consider the relationship between the distribution of suitable area for tree species and environmental factors, data on climate, soil, topography and vegetation classification were selected for model calculations, with a spatial resolution of 30s. Current and projected (2021–2,100) climate data were obtained from WorldClim.² Future climate scenarios were chosen from the BCC-CSM2-MR climate model within CMIP6, suitable for China’s geographical environment, which includes the sustainable development pathway limiting temperature increase to 2°C (SSP1-2.6), the moderate development pathway limiting temperature increase to 3°C (SSP2-4.5), and the conventional development pathway limiting temperature increase to 5°C (SSP5-8.5). Topographical variables, including elevation, slope, and aspect, were acquired from the Geographical National Condition Monitoring Cloud Platform³ and further processed using ArcGIS software. Additionally, soil variables comprising 56 factors were sourced from the Harmonized World Soil Database. The 1:1,000,000 Vegetation Data Set, sourced from the National Cryosphere Desert Data Center,⁴ was utilized for acquiring vegetation classification data for China. The distribution data of nature reserves in China were sourced from the Resource and Environment Data Cloud Platform of the Chinese Academy of Sciences.⁵ The study employed a 1:4,000,000 scale map of China along with its corresponding administrative division map, both sourced from the National Fundamental Geographic Information System.⁶

To mitigate the risk of overfitting in model predictions due to multicollinearity among environmental variables, correlation analysis was conducted using ENMTools. A Pearson correlation coefficient of 0.8 was defined as the threshold. For any two environmental factors with an absolute correlation coefficient ≥ 0.8, the one with the higher contribution rate based on the initial MaxEnt model results was retained (Supplementary Figure S2). Ultimately, 32 variables were selected from the initial 68 for evaluation purposes, as detailed in Supplementary Table S1.

The species distribution spot data, along with the processed environmental variables, were incorporated into the MaxEnt model for comprehensive analysis (Phillips et al., 2006) (Maxent version 3.4.4; http://www.cs.princeton.edu/wschapire/Maxent/). During model operation, data were randomly sampled and divided into a training subset (75% of the total data) and a testing subset (25% of the total data). The bootstrap method was chosen for replicated runs, while the jackknife technique was employed to quantify the individual contribution of environmental factors to model development (Peterson and Cohoon, 1999). Response curves were generated to show the relationship between species suitable area and climate factors (Moreno et al., 2011). The model was run 10 times, and the average of these ten runs was used as the final output. The accuracy of the model was evaluated with the area under the Receiver Operating Characteristic curve (AUC) (Supplementary Figure S3).

The final output is the potential habitat suitability distribution data for the species in the form of a habitat suitability index (HSI) ranging from 0 to 1. Regions with HSI values greater than or equal to the maximum test sensitivity plus specificity logistic threshold are considered suitable area for the species. These regions were subsequently categorized into four potential habitat classes: “high suitable” (>0.6), “medium suitable” (0.4–0.6), “low suitable” (0.2–0.4), and “not suitable” (<0.2).

Utilizing current climate data, this study simulated the suitable area range of P. neoveitchii and then conducted predictive simulations for three prospective climate scenarios (SSP1-2.6, SSP2-4.5, SSP5-8.5) to anticipate future habitat suitability. The simulation outcomes categorized future habitat polygons into three distinct groups: Expansion, indicating areas where suitability is projected to increase; Shrinkage, denoting regions where suitability is anticipated to decrease; and Stability, reflecting zones where suitability is forecasted to remain unchanged. Following this classification, the spatial extents of these categories were quantified and comprehensively described.

Utilizing the SDM Toolbox, a Python-powered GIS application, regions exhibiting HSI values exceeding 0.2 were designated as potentially suitable habitats, whereas those with HSI values below 0.2 were classified as unsuitable for the species’ persistence. This categorical distinction facilitated the subsequent analysis of habitat distributional changes. Subsequently, the changes in potential suitable area for P. neoveitchii under different climate scenarios were analyzed and categorized into expansion areas, stable areas, and contraction areas. The areas of these regions were then calculated (Cong et al., 2020).

Furthermore, utilizing the SDM Toolbox, we conducted a comparative analysis of the central points (centroids) within the current and projected suitable habitat areas to discern potential shifts in spatial distribution (Brown, 2014). Specifically, the centroids of P. neoveitchii distribution points under different climate scenarios were calculated first. Then, centroid vector files containing distance and direction information were created. Finally, the movement trajectories of the centroids were analyzed to confirm the change trends of the distribution of P. neoveitchii (Kong et al., 2021).

In this study, the ArcGIS platform was employed to quantitatively assess and interpret the spatial intersection between the projected suitable habitats of P. neoveitchii and the boundaries of designated nature reserves across varying temporal horizons. By adopting this methodology, we gain insights into the conservation priorities for P. neoveitchii, encompassing both contemporary climatic regimes and projected future scenarios, thereby facilitating the assessment of its vulnerability and resilience. Additionally, based on the overlap analysis results, this research identifies key areas for future conservation efforts, providing a basis for developing scientific conservation strategies.

The species distribution model was evaluated on the designated training and testing datasets. It exhibited an AUC score of 0.985, and the standard deviation is 0.004. This score is comparable to or even surpasses the typical AUC values reported in similar studies. For example, in a vulnerability assessment of Taxus wallichiana in the Indian Himalayan Region by Tiwary et al. (2024), the average training AUC was 0.96. The high AUC score of our model indicates excellent performance and reliable predictive results (Supplementary Figure S3), supporting the robustness and accuracy of the findings. Based on the current assessment, the total area of “high suitable” and “medium suitable” habitats for P. neoveitchii is approximately 1,740,000 square kilometers, primarily distributed in southeastern Gansu, southern Shaanxi, western Henan, and western Hubei. Observing the topographical map of China, the distribution areas are mainly located in the sunny slopes, semi-shaded slopes, and cliff areas of the Qinling-Bashan Mountains, Wushan Mountains, Wuyi Mountains, Hengduan Mountains, and their surrounding areas (Figure 3).

The MaxEnt framework evaluated the significance of various factors. The minimum temperature during the coldest month (Bio 6) emerged as the most influential, contributing 29.1% to model variation. Annual precipitation (Bio 12) followed closely, with a contribution of 25.1%. Temperature seasonality (Bio 4) and altitude also played significant roles, accounting for 9 and 7.7%, respectively. Together, these four factors comprised 70.9% of the model’s predictive power (Supplementary Table S1).

Response curves were utilized to determine the optimal bioclimatic conditions for P. neoveitchii. An HSI threshold of 0.2 was established. Climates conducive to P. neoveitchii’s growth are characterized by relatively cold temperatures that are not excessively low, moderate and stable precipitation, and distinct seasonal temperature variations. Higher altitudes, approximately 2,100 meters, provide a more favorable growth environment for the species (Figure 4).

The MaxEnt model simulations projected an estimated potential suitable habitat of approximately 4,226,000 km² for P. neoveitchii within China’s borders, constituting 4.48% of the country’s total territorial area (Supplementary Table S2). Under various projected future climatic scenarios, the extent of viable habitats exhibits a declining trend punctuated by periods of variability. Notably, the SSP5-8.5 scenario manifests the most pronounced decline in the extent of suitable habitats, particularly those of high suitability, while the SSP2-4.5 scenario exhibits a lesser yet notable decrease. In contrast, the SSP1-2.6 scenario demonstrates minimal alteration, preserving the largest extent of viable habitats. Remarkably, within this climate trajectory, the suitable habitat area experiences a marginal increase of 0.05% between 2080 and 2,100, as compared to the current projected distribution (Figure 5), underscoring its resilience under mitigated climate change. Regarding geographical shifts in habitat suitability, the SSP5-8.5 scenario elicits the most pronounced alterations in the extent of suitable habitat for P. neoveitchii, particularly evident in the northwestward and southeastward directions. The SSP2-4.5 scenario follows with a lesser degree of variation, while the SSP1-2.6 scenario exhibits minimal geographical redistribution (Figure 6).

The observed changes are closely linked to the degree of temperature increase under the three different climate scenarios. From SSP1-2.6 to SSP5-8.5, global warming increases significantly. Under the SSP1-2.6 scenario, global temperature increase is strictly limited to 2°C, resulting in relatively minor impacts on P. neoveitchii. Therefore, although the area of suitable habitats may decrease in the mid-term (2040–2060), it can recover to slightly above the current level in the long term (2080–2,100). This finding suggests that P. neoveitchii, under the auspices of a sustainable development trajectory, possesses a degree of climatic adaptability, potentially maintaining or even slightly expanding its suitable habitat area. In contrast, under the SSP5-8.5 scenario, global warming reaches 5°C, a relatively high level of temperature increase. The area of suitable habitats for P. neoveitchii decreases significantly, suggesting that under the conventional development pathway, P. neoveitchii will face immense survival pressure, leading to a substantial reduction in its suitable habitat range. The spatial distribution changes are also the most significant under this scenario, probably due to its sensitivity to extreme climate conditions such as high temperatures and droughts. In comparison, under the SSP2-4.5 scenario, the global temperature increase was limited to 3°C, falling between the other two scenarios. Therefore, both the area and spatial distribution changes of suitable habitats are intermediate between those observed in the SSP1-2.6 and SSP5-8.5 scenarios.

Currently, the centroid of P. neoveitchii distribution is situated in Ziyang County, Ankang City, Shaanxi Province, China (108°18′E, 32°33′N). Predictions under three distinct climate scenarios indicate that by 2080–2,100, the centroid’s migration trend predominantly leans northwestward, albeit with scenario-specific variations in path and distance. Specifically, under SSP1-2.6, the centroid is projected to relocate to Huixian County, Longnan City, Gansu Province (106°24′E, 33°40′N), approximately 215 km from its current position. In SSP2-4.5, the migration path is basically the same as in SSP1-2.6, but extends further, concluding in Huixian County, Longnan City, Gansu Province (105°43′E, 34°9′N), about 300 km away. SSP5-8.5 forecasts a more significant northward shift, settling in Lingtai County, Pingliang City, Gansu Province (107°9′E, 34°59′N), approximately 290 km distant from the current position (Figure 7). Furthermore, the importance of geographical variability, particularly in relation to conservation efforts at higher altitudes, is underscored by the observed migration trends of P. neoveitchii. As temperatures rise, the species is projected to migrate further towards higher latitudes and altitudes, seeking cooler environments. This migration towards higher altitudes not only highlights the species’ adaptive capacity but also poses significant conservation challenges. Effective conservation strategies must consider the shifting habitat requirements of P. neoveitchii and prioritize the protection of habitats at higher altitudes, which are likely to become crucial refuge areas under future climate scenarios.

It is noteworthy that across the four time periods from present to 2,100 under the three scenarios, the centroid exhibits a significant fluctuation towards the northwest direction during the period of 2060–2080. However, the amplitude of these fluctuations varies among scenarios, with SSP5-8.5 showing the largest fluctuation. Under this scenario, the centroid even shifts to the northeast of Gansu Province during 2060–2080, before settling approximately 107 kilometers northeast of its position during 2040–2060 by 2080–2,100. SSP2-4.5 exhibits a lesser fluctuation, with the centroid ending up approximately 40 kilometers away from its position during 2040–2060. SSP1-2.6, on the other hand, displays the smallest fluctuation, with the centroid in 2080–2,100 closely approximating its position during 2040–2060, just 5 kilometers away. These fluctuations are speculated to correlate with future climate variability. P. neoveitchii shows stronger adaptability and resilience under scenarios with lower temperature increases. Conversely, under SSP5-8.5, characterized by a higher temperature increase within the conventional development pathway, the centroid’s fluctuation is the most significant. This reflects significant challenges P. neoveitchii faces in adapting to new environments under extreme climate conditions. While the centroid returns somewhat by 2080–2,100, its final position is considerably farther northeast, approximately 107 kilometers away from its original location. This suggests that high temperatures pose substantial pressure and uncertainty during the migration of P. neoveitchii, making it difficult to find stable habitats and moving away from its ecological resources.

The protection work of P. neoveitchii currently faces severe challenge, with protected areas covering only 6.48% of its distribution range. These protected areas are mainly concentrated in central Sichuan Province and the border areas of Chongqing, Hubei, and parts of Hubei. However, in critical regions of its native distribution, such as western Gansu, southern Shaanxi, western Hubei, and western Henan, the coverage of nature reserves is notably insufficient. Moreover, these protected areas are severely fragmented, posing serious challenges to the effectiveness in their conservation (Figure 8).

Facing the uncertainties brought by future climate change, the conservation status of P. neoveitchii appears to be undergoing a positive transformation. Across three distinct climate scenarios (SSP1-2.6, SSP2-4.5, SSP5-8.5), the projected gradual reduction in the suitable habitat area for P. neoveitchii is counterbalanced by an impending increase in the proportion of national nature reserves encompassed within these areas (Supplementary Figure S2). Of particular note is the projection that, within the timeframe spanning 2080 to 2,100, under the most severe climate scenario (SSP5-8.5), the proportion of national natural reserves residing within the suitable habitat area of P. neoveitchii is anticipated to reach 10.2%. For SSP2-4.5 and SSP1-2.6, these proportions are projected to be 8.66 and 7.52%, respectively. This positive trend may be associated with the potential northwestward migration of P. neoveitchii as predicted by the models. Natural reserves in northwest China are relatively more intact and more extensive. Thus, although the suitable habitat area for P. neoveitchii decreases, they may migrate to these more comprehensive natural reserves, which provides a relatively optimistic outlook for the conservation of P. neoveitchii.

The distribution patterns of species are shaped by a complex interplay of both abiotic and biotic factors. Abiotic factors include climate, terrain, and soil properties, while biotic factors encompass competition dynamics, growth traits, dispersal and migration capabilities, natural regeneration processes, herbivory, and various human-induced impacts such as deforestation, overgrazing, tourism pressures, and industrial and mining activities (Ballardin et al., 2024; Campedelli et al., 2015; Wesemeyer et al., 2023). Exemplified by the Qilian Mountains, where from 1989 to 1998, a surge in human activities resulted in a marked decline of 33.5% in forest vegetation, thereby contributing to the degradation of forest ecosystems (Guo et al., 2019). Furthermore, limitations in species occurrence data can introduce deviations between simulated results and actual species locations. Additionally, the modeling accuracy for estimating species’ geographical distribution ranges may be impacted by variations in the resolution of selected environmental variables. This study aimed to mitigate such errors by utilizing field-based distribution data collected since 2022, eliminating redundant distribution data, removing some environmental data with high correlation coefficients, and performing multiple iterations. The evaluation metrics, notably the Area Under the Receiver Operating Characteristic Curve (AUC) and the training-based omission rate, underscore the reasonability of the delineated distribution range for P. neoveitchii in this investigation. These findings carry substantial implications for the strategic cultivation, effective conservation, comprehensive management, and sustainable utilization strategies of this species.

The spatial distribution of plant species at a given scale is primarily governed by climatic determinants, with hydrothermal parameters (pertaining to water and heat) occupying a paramount position in influencing this distribution (Tiwary et al., 2024). Model simulations have revealed that Min Temperature of Coldest Month, Annual Precipitation and Temperature Seasonality are the primary climatic factors influencing the distribution of P. neoveitchii. Specifically, lower temperatures, moderate and stable precipitation and distinct seasonal temperature variations are particularly favorable for the growth of P. neoveitchii. Our discovery converges significantly with the findings presented by Dad et al. (2023), who identified Annual Precipitation and annual mean temperature as the three pivotal environmental factors that exert the most profound influence on the distribution of wild gymnosperms in the northwestern Himalayan region. The slight differences between the two conclusions can be attributed to variations in plant species and study areas.

In the modeling calculations, apart from climate, soil conditions also have a significant, albeit relatively lesser, impact on the growth of P. neoveitchii. The species prefers gravelly alluvial soil, which offers excellent drainage and aeration properties, essential for healthy root growth. From a topographical perspective, the formation of alluvial soil is closely associated with hydrodynamic processes such as rivers and floods, indicating that P. neoveitchii may be more adapted to environments like wide river valleys, riverbanks, and floodplains adjacent to plains. Simultaneously, the survival probability of P. neoveitchii is closely associated with the vegetation type in which it resides. According to model analyses, the most suitable vegetation type for P. neoveitchii is subtropical coniferous and broad-leaved mixed forests. However, this vegetation type has a relatively small distribution area in China, primarily located in the south-central and western regions of Yunnan, the western part of Sichuan, the southeastern section of Tibet, and the central mountainous areas of Taiwan. These vegetation types are mainly distributed in the mountainous environments of the mid-subtropical and temperate zones, which aligns with their actual distribution. This discovery not only reveals the ecological habits and distribution patterns of P. neoveitchii, but also provides important insights for the development of relevant protection measures. In future conservation efforts, we should prioritize areas with suitable vegetation types and employ measures such as afforestation and ecological restoration to create a better habitat for rare plants like P. neoveitchii.

The future climate scenarios are projected to have significant impacts on both the size and distribution range of viable habitats for P. neoveitchii. Specifically, adhering to the Sustainable Development Scenario (SSP1-2.6), the predicted trend indicates that the suitable habitat area for P. neoveitchii will maintain stability or undergo a modest expansion. However, under the Intermediate Development Pathway (SSP2-4.5) and the Conventional Development Pathway (SSP5-8.5), its suitable habitat area will be reduced and disrupted to varying degrees. Our forecast concurs with the conclusions outlined in the Sixth Assessment Report of the Intergovernmental Panel on Climate Change (IPCC) published in 2023 (Intergovernmental Panel on Climate, 2023), reinforcing a shared understanding. Additionally, it harmonizes with the research findings reported by Gao et al. (2021), who explored the evolving patterns of suitable habitat for European spruce under future climatic scenarios. The modeling results indicate that under current climatic conditions, the distribution center of P. neoveitchii is located in Ziyang County, Ankang City, Shaanxi Province (108°18′E, 32°33′N). Under three different future climate scenarios, the migration direction of P. neoveitchii is generally northwest, though the specific paths and distances vary depending on the scenario. Our findings are congruent with the conclusions presented in the IPCC’s AR6 (Intergovernmental Panel on Climate, 2023), indicating that climatic alterations are anticipated to exert pronounced effects on ecological systems located within specific low-latitude areas of the Northern Hemisphere. Many plant species may adapt to climate change by migrating to higher latitudes or altitudes. Furthermore, our findings resonate with the projections made by He et al. (2023) and Xiao et al. (2024), as well as those reported in another study (Duan et al., 2022), which showed that the distribution center of P. densiflora shifts to a small northward migration from current to future climate scenarios. While all these studies anticipate a northwestward migration of alpine coniferous forests in response to ongoing climate warming trends, the specific migration distances and angles may vary slightly among different plant species.

Facing the uncertainties brought by future climate change, the conservation status of P. neoveitchii shows a positive shift. Despite the projected gradual reduction in the extent of suitable habitat for P. neoveitchii, the fraction of this habitat encompassed within national nature reserves demonstrates an upward trajectory, affording a moderately optimistic perspective on its conservation prospects. Nonetheless, it is crucial to acknowledge that the model’s projections solely rely on climatic simulations and may therefore possess inherent limitations. Given that P. neoveitchii is currently endangered and the populations are highly isolated from each other, it remains uncertain whether P. neoveitchii can migrate northwestward through its own dispersal mechanisms and the extent of such migration.

To address this scenario, strategic interventions encompassing in-situ and ex-situ conservation efforts, safeguarding germplasm repositories, implementing stringent enclosure measures, and enhancing reproductive programs on a larger scale ought to be prioritized. (1) Firstly, efforts should be intensified to protect the native distribution areas of P. neoveitchii, especially the regions where protected areas are severely fragmented. By restoring and connecting fragmented reserves, a more stable and continuous habitat for P. neoveitchii can be provided. (2) Secondly, a germplasm bank for P. neoveitchii should be established. The favorable conditions of mature individuals, such as good seed production and high germination rates, should be utilized optimally by timely and appropriate seed collection. Nurseries should be established locally to expand the artificial populations and study rapid propagation methods for P. neoveitchii. Additionally, transplantation and seedling cultivation among natural populations should be conducted to ensure genetic exchange and maximize the genetic diversity of P. neoveitchii. This would provide germplasm and materials for further breeding and in-depth conservation efforts. (3) Furthermore, the migration process of P. neoveitchii needs to be monitored and studied to understand its migration routes, speed, and quantity. This information would help in accurately predicting its future distribution and guide the planning and construction of protected areas to ensure comprehensive protection in new distribution areas. (4) Finally, promoting the implementation of sustainable development strategies to limit global warming levels is also crucial for preserving the habitat of P. neoveitchii and other organisms.

Currently, research on P. neoveitchii primarily focuses on population size, community characteristics, and endangerment mechanisms, with relatively little attention given to habitat simulation, cultivation methods, and seed germination. Research on P. neoveitchii should be deepened further into ecological physiology and molecular levels, thereby gaining a detailed understanding of their growth dynamics and other related issues. An in-depth study of P. neoveitchii would allow for a comprehensive analysis of the causes of its endangerment, intra-population relationships, seed germination patterns, and genetic differentiation. These insights could then help to develop more effective conservation measures, ultimately aiming to rescue this “living fossil” from endangerment.