Potentially suitable habitat prediction of Pinus massoniana Lamb. in China under climate change using Maxent model

Background Pinus massoniana is an important timber species with high ecological and economic value in southern China and a pioneer species for the reforestation of barren mountains. The structure and function of the existing Masson pine ecosystem have been seriously affected by worsening habitats under current climate change. Its suitable habitat is likely to change greatly in the near future. Methods To estimate the potential geographic distribution of P. massoniana and its response to climate change, the Maxent model was selected to simulate the potentially suitable habitat and corresponding changes in the distribution pattern of P. massoniana under current and future climate scenarios (SSP1-2.6, SSP2-4.5, and SSP5-8.5) with two periods (2050s and 2090s). Results Under current climate conditions, the total suitable habitat area of P. massoniana was 2.08 × 106 km2, including 0.76 × 106 km2 of highly suitable habitat, concentrated mainly in Jiangxi and Zhejiang, central and southeastern Hunan, northern Fujian, central and western Chongqing, southern Anhui, central and surrounding areas of Guangdong, central and eastern Guangxi, and northern Taiwan. The areas of moderately suitable habitat and poorly suitable habitat were 0.87 × 106 km2 and 0.45 × 106 km2, respectively. Temperature and precipitation appear to be the most important predictors: precipitation of the driest month (14.7–215.6 mm), minimum temperature of the coldest month (−3.5–13.8°C), annual temperature range (8.1–32.9°C), and mean temperature of the warmest quarter (23.6–34.7°C). P. massoniana are predicted to expand their potential distribution under future climate change: by the end of this century, their total suitable habitat area increased 0.26 × 106 km2 (10.61%) and 0.45 × 106 km2 (17.05%) under the most moderate (SSP1-2.6) and severe (SSP5-8.5) warming scenarios, respectively, by mainly extending northward. Conclusion Under the different future climate scenarios, the total suitable habitat area of P. massoniana increased by mainly extending northward. Overall, our study clarifies the potential habitat distribution of P. massoniana and provides a critical empirical reference for future P. massoniana conservation and planting practices.


Introduction
Pinus massoniana Lamb. (Masson pine), a common and dominant tree species characterized by strong tolerance and adaptability to depleted soil, with a northern boundary of the natural distribution of the southern Henan and Shandong Province in China. It has been used as an important afforestation tree species in southern China considering its high economic (e.g., timber resource) and ecological (e.g., carbon sequestration, water, and soil conservation, forest development) values Shao et al., 2022), which accounting for 20% of the total afforestation area (National forestry and Grassland administration, 2019). However, P. massoniana is sensitive to climate change , for example, absolute minimum temperature, precipitation, and soil acidity-alkalinity (Jiang et al., 2019). Due to human or environmental factors, P. massoniana has presented problems such as low productivity and poor ecological function in the current situation in China.
Climate, especially temperature and precipitation, are confirmed the main environmental variables affecting the geographical distribution and growth of species worldwide (Jochum et al., 2007;Dyderski et al., 2018). The global mean temperature is predicted to increase by at least 0.3 • C and up to 4.8 • C by the end of this century as global warming continues (Stocker et al., 2013). There was enough evidence that current climate warming has caused many plants to migrate to higher altitudes and latitudes (Mckenney et al., 2007;Bertin, 2008;Thurm et al., 2018). Meanwhile, seasonal drought and heavy rain, hot summer and severe winter are occurring more and more frequently in subtropical areas, which may expose species to extinction, regional biodiversity and ecological security patterns change, and ecosystem functions reduction, forcing them to adapt to new condition or shift of their geographic distributions (Alexander et al., 2018;Feeley et al., 2020;He et al., 2021;Jian et al., 2022;Zhang et al., 2022). Hence, under the context of global change, increasing global warming and frequency of extreme events (extreme heat, cold, drought, etc.), how P. massoniana changes its distribution area to adapt to future climate change is a problem worth discussing. These changes will challenge forest management, utilization and preservation in the future (Dyderski et al., 2018).
Species distribution models (SDMs) or Ecological Niche Models (ENM), are mathematical models that reflect the distribution of suitable habitats of species on a large scale based on species presence or richness data and environmental factor data, which include but are not limited to genetic algorithms for rule-set production (GARP) (Qin et al., 2015), generalized linear models (GLM) (Lopatin et al., 2016), classification and regression tree (CART) (Cao et al., 2005), Bioclim (Zhu et al., 2020b), random forest (RF) (Kolli et al., 2022) and maximum entropy (Maxent) (Yan et al., 2020). The Maxent model is one of the most widely used models because of its good performance in accuracy and friendly interface, even if there are a few occurrence points or training data (Elith et al., 2006;Hernandez et al., 2006;Phillips et al., 2006;Kaky et al., 2020). It has been widely used in potential habitat predicting or geographic distribution migration under global change of invasive species, Chinese Herbs, commercial crops, desert plants, and endangered species (Bose et al., 2016;Guo et al., 2017;Liu et al., 2019;Guan et al., 2022), covering plants, animals, bacteria, and fungi Prevéy et al., 2020;Ye et al., 2022), and achieved satisfactory results.
In recent years, studies on P. massoniana have been conducted in physiological and ecological processes (Lei et al., 2022;Wang et al., 2022;Yuan et al., 2022), nutrient cycle (Dong et al., 2021), and forest management , little attention has been paid to the range shift of this widespread species in China. While this task will be a key step in the development of climate change adaptation strategies and forest management for P. massoniana. For this purpose, we used the Maxent model and ArcGIS to predict the potential suitable habitats for P. massoniana under current and various possible future climatic conditions according to P. massoniana distribution data and associated environmental variables. The main objectives of this research included the following: (1) predicting the potential suitable habitat for P. massoniana under the current climate condition; (2) screening the main environmental variables limiting the geographic distribution of P. massoniana; (3) predicting the potential suitable habitats for P. massoniana under future climate conditions and exploring changes in geographic distribution patterns of P. massoniana in the present and future. Our results provide practical guidance for the future management conservation and planting practices of P. massoniana in China.

Species occurrence data collection
The distribution records of P. massoniana we used were obtained from the Global Biodiversity Information Facility (GBIF) 1 and the National Plant Specimen Resource Center (CVH) 2 . These data were integrated into the present dataset, and only the records within the territory of China were kept. We collected 371 P. massoniana distribution records after deleting invalid and duplicate records. Only one distribution record was chosen from each 2.5' × 2.5' grid to minimize sampling bias. Ultimately, 264 useful records were selected to create the Maxent model (Figure 1).

Environmental variables and processing
Besides climatic factors, there are many other factors that influence the geographic distribution of species, such as topography and soil factors (Bradie and Leung, 2017). Numerous studies have used only climatic variables for simulation, yet these variables cannot fully describe all of the complex processes that affect species distribution (Thuiller et al., 2019). Although these models perform well, their predicted species distributions may have large biases from the true distributions and carry some prediction risk (Beale et al., 2008). Thus, a total of 37 environmental variables, including bioclimatic, topographic, and soil variables, were considered in the prediction model to better identify which variables have a greater impact on the potential distribution of P. massoniana ( Table 1). The 19 bioclimatic variables were obtained with a resolution of 2.5' from the Worldclim database 3 for both current  and future (2050s, 2090s) scenarios. For future scenarios, 2050s bioclimatic data represented the average of 2041 to 2060, and 2090s data represented the average of 2081 to 2100. To estimate future climate change, we used general circulation model (GCM) forecasts under Shared Socio-economic Pathways (SSPs) scenarios, which were started by the IPCC's Coupled Model Intercomparison Project Phase 6 (CMIP6) (Popp et al., 2017). The CMIP6 models have improved the simulation capability of the climatological precipitation distribution and East Asian summer monsoon (EASM) climate patterns, and also have strong simulation capability in Chinese climate indices (Xin et al., 2020;Zhu et al., 2020a). The GCM we selected for the present study was the BCC-CSM2-MR climate system model created by the National Climate Center, which performs best in reproducing the anomalous rainfall pattern in eastern China with the highest spatial correlation and can realistically simulate the rainfall anomalies in eastern China (Xin et al., 2020). It includes the following three emission scenarios: SSP1-2.6, SSP2-4.5, and SSP5-8.5. (Fick and Hijmans, 2017). SSP1-2.6 is a sustainability scenario with low levels of greenhouse gas emissions. SSP2-4.5 indicates that greenhouse gas emissions are at a moderate level, which means that the future socio-economic development pattern will continue to follow the current pattern. In contrast, SSP5-8.5 is a rapid socio-economic development scenario with high levels of greenhouse gas emissions (Eyring et al., 2016). Digital elevation model data were obtained from Geospatial Data 3 http://www.worldclim.org/ Cloud 4 in order to extract elevation, slope, and aspect data. Soil variables were derived from the Harmonized World Soil Database, 5 including 15 basic soil indicators.
Potential correlations between environmental variables may affect the accuracy of prediction (Pearson et al., 2006). Studies that ignore the redundant information that highly correlated variables introduced during species distribution modeling typically have low accuracy (Duflot et al., 2018). Eliminating redundant variables can reduce multicollinearity among variables and thus give the model stronger predictive power (Ashcroft et al., 2011;Yi et al., 2016). To improve the accuracy of the model, we selected environmental variables according to the following principles: (1) the variables imported into Maxent with contribution rates < 1% were removed; (2) after calculation of Pearson correlation coefficients between the remaining variables, those with correlation coefficients (r) less than 0.8 were retained, while variable with less contribution was eliminated when the r ≥ 0.8 (Lu et al., 2020). This approach effectively eliminated the environmental variables with low contribution and high correlation. Ultimately, nine main environmental variables were selected for modeling (Table 1).

Parameter setting and accuracy verification of Maxent model
The Maxent model was created and replicated 10 times based on the environmental variables and distribution data; 25% of the distribution data were randomly selected for model testing, while the rest 75% were utilized for model training (Phillips, 2008). The replicated run type was set as "Cross-validate, " the maximum iteration was set to 500, and the output format was logistic. To Reliability test of the distribution model created for Pinus massoniana.
identify the key environmental variables affecting P. massoniana distribution, the jackknife test and response curves were used to assess the importance of environmental variables and their effects on P. massoniana distribution (Narouei-Khandan et al., 2016). We evaluated the model performance using the area under the receiver operating curve (AUC) (Qin et al., 2017). In this evaluation, poor performance falls in 0.5 < AUC ≤ 0.7, moderate performance falls in 0.7 < AUC ≤ 0.9, and excellent performance falls in 0.9 < AUC ≤ 1 (Peterson et al., 2011).

Division of suitable habitat grades
We used ArcGIS 10.4 to process the suitability classification of P. massoniana. We quantified the potential distribution probabilities (P) of P. massoniana from 0 to 1 and classified them according to the method used by Zhao et al. (2021), and the habitat suitability of P. massoniana was classified into four categories: unsuitable habitat (P < 0.1), poorly suitable habitat (0.1 ≤ P < 0.3), moderately suitable habitat (0.3 ≤ P < 0.5), and highly suitable habitat (0.5 ≤ P < 1).
We utilized the SDMToolbox in ArcGIS to visualize habitat changes of P. massoniana in different climate conditions, and we compared current and future potential habitat distributions (Brown, 2014;Hu et al., 2017). The results were divided into three components: decreased area, increased area, and invariant area . The centroids of suitable habitats and their migration distances were also calculated for current and future climate scenarios using the SDM Toolbox in ArcGIS.

Model performance and major environmental variables
Our validation results showed a high accuracy of the P. massoniana distribution prediction in China (Figure 2), with an AUC value of 0.914 for habitat suitability maps. Thus, the model  Table 2. The four environmental variables with the highest contribution rate were precipitation of the driest month (Bio14, 70.5%), minimum temperature of the coldest month (Bio6, 15%), isothermality (Bio3, 4.9%), and annual temperature range (Bio7, 4.2%). When only a single variable was used for modeling, the four variables with the highest regularized training gain were minimum temperature of the coldest month (Bio6), precipitation of the driest month (Bio14), annual temperature range (Bio7), and mean temperature of the warmest quarter (Bio10). The regularized training gains for the four variables mentioned above were higher than others significantly, indicating that these variables contain information absent from other variables (Figure 3).
We used Maxent to plot response curves with single environmental variables (Figure 4) to elucidate the relationship between the probability of P. massoniana existence and key environmental variables, as well as the range of suitability of environmental variables determined by applying a threshold based on the moderately suitable habitat. The suitable habitat conditions relating to the key environmental variables for P. massoniana are: precipitation of the driest month (Bio14), 14.7-215.6 mm ( Figure 4A); minimum temperature of the coldest month (Bio6), -3.5 to 13.8 • C ( Figure 4B); annual temperature range (Bio7), 8.1-32.9 • C ( Figure 4C); mean temperature of the warmest quarter (Bio10), 23.6-34.7 • C ( Figure 4D).

Potential distribution under current climate conditions
The potentially suitable habitat for P. massoniana distribution under the current climatic scenarios is shown in Figure 5. The total suitable habitat area was estimated to be 2.08 × 10 6 km 2 , accounting for 21.67% of the whole study area. The highly suitable habitat was mainly located in Jiangxi and Zhejiang, central and southeastern Hunan, northern Fujian, central and Jackknife test of variable importance for Pinus massoniana ("T_USDA_TEX" is topsoil USDA texture classification; "T_pH_H 2 O" is topsoil pH (H 2 O); "Slope" is slope; "Bio7" is annual temperature range; "Bio6" is minimum temperature of the coldest month; "Bio3" is isothermality; "Bio14" is precipitation of the driest month; "Bio10" is mean temperature of the warmest quarter; "Aspect" is aspect). western Chongqing, southern Anhui, central and surrounding areas of Guangdong, central and eastern Guangxi, and northern Taiwan, with a combined area of 0.76 × 10 6 km 2 , accounting for 7.92% of the study area. The moderately suitable habitat was mainly located in Hubei and Guizhou, southern Fujian, eastern Guangdong, western Guangxi, eastern Sichuan, southern Shaanxi, central Anhui, and central Taiwan, with a combined area of 0.87 × 10 6 km 2 , accounting for 9.06% of the study area. The poorly suitable habitat was mainly located in Jiangsu, western Guizhou, southeastern Tibet, western Guangxi, southern Gansu, central Shaanxi, central and southwestern Henan, eastern Shandong, southern Taiwan, and northern Hainan, with a combined area of 0.45 × 10 6 km 2 , accounting for 4.69% of the study area.

Potential distribution under future climate conditions
The suitability of P. massoniana was classified to obtain predictions of suitable habitats for P. massoniana in the 2050s and 2090s under future climate scenarios (Figure 6 and Supplementary Figure 1). Overall, the total area of suitable habitats under future climate scenarios was predicted to increase compared to the suitable habitat under current climate scenarios, with a significant increase in the highly and poorly suitable habitats and a decrease in the moderately suitable habitat.
Under the SSP1-2.6 climate scenario, the total suitable habitat for P. massoniana in the 2050s was predicted to be 2.22 × 10 6 km 2 , accounting for 23.13% of the whole study area. The areas with highly, moderately, and poorly suitable habitats were predicted to be 1.09 × 10 6 km 2 , 0.62 × 10 6 km 2 , and 0.51 × 10 6 km 2 , respectively, accounting for 11.35, 6.46, and 5.31% of the study area. The highly and poorly suitable habitat areas were predicted to increase by 0.32 × 10 6 km 2 and 0.061 × 10 6 km 2 , respectively, in comparison to the corresponding suitable habitat areas in the current climate, whereas the moderately suitable habitat area was predicted to decrease by 0.25 × 10 6 km 2 . By the 2090s, under the SSP1-2.6 climate scenario, the total suitable habitat for P. massoniana was predicted to be 2.33 × 10 6 km 2 , accounting for 24.27% of the whole study area. The areas with highly, moderately, and poorly suitable habitats were predicted to be 1.08 × 10 6 km 2 , 0.72 × 10 6 km 2 , and 0.53 × 10 6 km 2 , respectively, accounting for 11.25, 7.5, and 5.52% of the study area. The highly and poorly suitable habitat areas were predicted to increase by 0.32 × 10 6 km 2 and 0.083 × 10 6 km 2 , respectively, relative to the suitable habitat areas in the current climate, whereas the moderately suitable habitat area was predicted to decrease by 0.15 × 10 6 km 2 (Figure 6). Under the SSP5-8.5 climate scenario, the total suitable habitat for P. massoniana in the 2050s was predicted to be 2.40 × 10 6 km 2 , accounting for 25% of the whole study area. The areas with highly, moderately, and poorly suitable habitats were predicted to be 1.31 × 10 6 km 2 , 0.58 × 10 6 km 2 , and 0.51 × 10 6 km 2 , respectively, accounting for 13.65, 6.04, and 5.31% of the study area. The highly and poorly suitable habitat areas were predicted to increase by 0.55 × 10 6 km 2 and 0.064 × 10 6 km 2 , respectively, in comparison to the corresponding suitable habitat areas in the current climate, whereas the moderately suitable habitat area was predicted to decrease by 0.29 × 10 6 km 2 . By the 2090s, under the SSP5-8.5 climate scenario, the total suitable habitat for P. massoniana was predicted to be 2.53 × 10 6 km 2 , accounting for 26.35% of the whole study area. The areas with highly, moderately, and poorly suitable habitats were predicted to be 1.66 × 10 6 km 2 , 0.44 × 10 6 km 2 , and 0.43 × 10 6 km 2 , respectively, accounting for 17.29, 4.58, and 4.48% of the study area. The highly suitable habitat area was predicted to increase by 0.90 × 10 6 km 2 , relative to the suitable habitat areas in the current climate, whereas the moderately and poorly suitable habitat areas were predicted to decrease by 0.43 × 10 6 km 2 and 0.016 × 10 6 km 2 , respectively (Figure 6).

Potential habitat changes under future climate conditions
The projections indicated that the total suitable habitats for P. massoniana in the future periods showed an increasing trend compared to the current potential distribution, with a general expansion to the north (Figure 7 and Supplementary  Figure 2). The newly suitable habitat areas were mainly located in Shandong and Hebei, in addition, northern Jiangsu, northeastern Henan, northern Anhui, southwestern Shanxi, central Shaanxi, southeastern Gansu, especially Xinjiang and Tibet gradually suitable for P. massoniana distribution ( Figure 7D). The decreased area was small and mainly concentrated in a small part of Hainan. The predicted suitable habitat for P. massoniana differed less between the two different future periods under the same SSP climate scenarios, while it differed more between the different SSP climate scenarios within the same future periods. The suitable habitat changed more substantially under the SSP5-8.5 climate scenario.
To further evaluate the effects of different climatic scenarios on the potential geographic distribution of P. massoniana in the 2050s and 2090s, the area of suitable habitat change was predicted for each of the periods and each of the scenarios ( Table 3). Under the SSP1-2.6 climate scenario, the increased area, decreased area, invariant area, and total change in the 2050s were predicted to be about 0.15 × 10 6 , 0.015 × 10 6 , 2.18 × 10 6 , and 0.14 × 10 6 km 2 , accounting for 6.44, 0.64, 93.56, and 6.01% of the total suitable area, respectively. By the 2090s, under the same climate scenario, the increased area, decreased area, invariant area, and total change were predicted to be about 0.26 × 10 6 , 0.0076 × 10 6 , 2.19 × 10 6 , and 0.25 × 10 6 km 2 , accounting for 10.61, 0.31, 89.39, and 10.20% of the total suitable area, respectively.

Suitable habitat centroid trends under future climate conditions
The suitable habitat of P. massoniana in different periods and different climate scenarios was extracted by ArcGIS, and the location and migration distance of the total suitable habitat centroids were calculated accordingly (Figure 8)

Main environmental variables affecting the occurrence of P. massoniana
According to the Maxent model predictions, climatic variables are the main environmental variables affecting the potential geographic distribution of P. massoniana. The response curves of environmental variables were used to illustrate the impacts of the main environmental variables on the potential geographic distribution of P. massoniana, and it was concluded that the minimum temperature of the coldest month and the precipitation of the driest month (i.e., temperature and precipitation in winter) were the key environmental variables influencing the distribution of P. massoniana. Physiological ecology studies of P. massoniana have come to similar conclusions that higher winter temperatures and adequate winter precipitation play a critical role in tree growth. Low winter temperatures can lead to wilting of needle tips and may also reduce the activity of tree roots and buds, or even cause more damage to tree roots, buds or branches, to the detriment of plant growth (Pederson et al., 2004;Zhu et al., 2009;Shi et al., 2010). On the contrary, higher winter temperatures can avoid frost damage to foliar tissues and ensure normal metabolic activities of trees, and also enable trees to enter the growing season earlier, thus extending the growing season and improving photosynthetic efficiency, which is beneficial to tree growth (Northcote and Hartman, 2008;Chen et al., 2016). On the other hand, sufficient winter precipitation increases the soil water content, which facilitates the accumulation of carbohydrates, thus providing more nutritional support for the cambium activities of trees in the following year and for tree growth in the early part of the next growing season (Jiao et al., 2021), which is also consistent with its growth in warm and moist environments (Yan et al., 2019). We assessed the contributions of the main environmental variables, finding that the contributions related to temperature and precipitation were 20.2 and 70.5%, respectively, indicating that the effect of precipitation was greater than that of temperature. This supports previous findings that precipitation plays a large impact on the potential distribution of P. massoniana than temperature (Zhang et al., 2018;Jiang et al., 2019). Water availability is essential for plant seed germination, and both long-term and short-term water deficits can constrain the growth and development of trees and seed germination (Han and Ding, 2012). Under drought stress, P. massoniana suffers from different degrees of injury, resulting in corresponding physiological and biochemical changes. Water deficit significantly affects root vigor and photosynthetic rate, which eventually causes changes in external morphology and inhibits the growth of P. massoniana (Tan et al., 2017). Therefore, precipitation should be considered a critical environmental variable informing the introduction and cultivation of P. massoniana.

Analysis of current potential suitable habitat
The main potential habitat areas of P. massoniana in our study are south of the Qinling Mountains-Huaihe River Line, and spread outward in progressively decreasing habitat classes toward northern and northwestern China. The Qinling Mountains-Huaihe River Line is the dividing line in China between the warm temperate and subtropical zones as well as the humid and semi-humid zones, with distinctly different climate characteristics (such as temperature and precipitation) between the areas north and south. The suitable areas of P. massoniana are dominated by subtropical monsoon or subtropical monsoon humid climates, which are characterized by hot summers, warm winters, and four distinct seasons . The warm and humid climatic factors may be the main factors affecting the distribution of P. massoniana. In addition, the jackknife test showed that soil type and soil pH were second only to climatic variables in their importance to the P. massoniana distribution (Figure 3). P. massoniana is a species that prefers acidic soil and is intolerant of saline soils . Most of the soil types in the suitable areas mentioned above are dominated by red, yellow, or brown loams, with varying levels of aluminum in the soil resulting in neutral to acidic soil pH values, conditions that are favorable for the growth of P. massoniana (Meng and Shi, 2017). In addition, soil salinization mainly occurs in arid areas lacking moisture or in low-lying coastal areas subject to seawater intrusion (Lin and Tang, 2005;Guo, 2018), and salinized lands in China are mainly distributed in the northwest, northeast, north, and coastal areas of China, with relatively little land salinization occurring in inland areas south of the Qinling Mountains-Huaihe River Line, and less salinized soils are also favorable for the growth of P. massoniana. This suggests that soil variables also influence the potential distribution of P. massoniana to some extent, which is consistent with our findings. In summary, it may be the suitable climate and soil environment that have ultimately influenced the formation of the modern geographical distribution pattern of P. massoniana. Meanwhile, we suggest the development of P. massoniana forests in the inland areas south of the Qinling Mountains-Huaihe River Line. However, P. massoniana afforestation should not be attempted in coastal areas, despite their suitable climatic conditions, because soil salinization may not be conducive to the growth of P. massoniana.

Analysis of future potential suitable habitat
Under future climate scenarios, the geographic distribution pattern of P. massoniana habitat generally showed a northward migration trend. The area of suitable habitat tended to increase compared to the current period, particularly under the SSP5-8.5 climate scenario in the 2090s, in which the predicted area of new suitable habitat was the largest. Numerous studies have shown that species distribution changes are closely related to climate warming and that climate change may lead to significant changes in species distribution patterns (Parmesan and Yohe, 2003;Broennimann et al., 2006;Chen et al., 2011;, resulting in their migration to the north (Hickling et al., 2006;Bertrand et al., 2011).
Long-term climate observations suggest that under future climate conditions, the annual mean precipitation in China will increase (0-20%), especially in northern and northwestern China; moreover, the annual mean surface temperature in China is expected to increase by 2.7-2.9 • C (Jiang and Fu, 2012). Additionally, the temperature increase would be more pronounced under the SSP5-8.5 climate scenario with the highest greenhouse gas emissions, which would further promote the northward migration of tree species Xie et al., 2021). Under the three future climate scenarios considered (i.e., , the new suitable habitats in the 2050s and 2090s would be located mainly in northern and northwestern China. These results indicate that these regions will be affected by climate change with increased temperature and precipitation, resulting in an increase in the total suitable area for P. massoniana and a shift toward higher latitudes, consistent with its preference for warm and humid environments.
Therefore, in the context of climate change, it is important to achieve the reasonable utilization of P. massoniana resources in different suitable habitats. We should adopt a more active management mode to promote the development and utilization of P. massoniana resources in highly suitable habitats. Whereas in low suitable areas, we should control the exploitation of P. massoniana to ensure ecological benefits, considering masson pines of poor growth are more susceptible to diseases, such as Bursaphelenchus xylophilus. While expanding the forest resources of P. massoniana, the ecological service function should be improved to finally achieve a win-win situation for both economic and ecological benefits . Total suitable habitat centroid distribution shifts under climate change for Pinus massoniana.

Frontiers in Forests and Global Change
The predicted spatial patterns of change in this study were based only on the response of P. massoniana to future climate change without taking into account other factors that limit its migration. Migration capacity has important effects on the adaptation of tree species to future climate change (Corlett and Westcott, 2013), which is limited by plant community competition, geographic barriers, human activities, and land use (Zhang et al., , 2021. However, it is currently difficult to integrate all influencing factors into one model to simulate the potential distribution of species . On the other hand, introducing too many variables into the model may lead to more collinearity-related problems, the goodness-of-fit of the model may not be better than the current model, and the effects of key variables may be weakened . In addition, wind dispersal and animal dispersal are the main modes of seed dispersal in P. massoniana, and their spread may be limited (Yang et al., 2013). Seeds are affected by wind speed, wind direction and topography during wind dispersal (Pan et al., 2014). Seed dispersal is also related to the behavior of dispersal animals, and the movement patterns of dispersal animals can vary in different habitats, thus affecting the distribution pattern of seeds (Barraquand and Benhamou, 2008;Van Loon et al., 2011). Because we did not take these constraints into account, the increase in suitable habitat predicted by the model and the distance that suitable habitat centroids migrate may have been overestimated. Nevertheless, the migration changes in suitable habitats predicted by this study in future periods are consistent with its growth habit and point the way to future migration trends of P. massoniana. Consequently, this study is still an important reference value for the potential suitable habitat migration of P. massoniana in the context of climate change. However, these limiting factors should be accounted for in future studies to more accurately predict the distribution of P. massoniana.

Conservation and cultivation recommendations
Climate change has already led to reduced resistance and resilience of boreal and tropical forests (Forzieri et al., 2022), and increased tree mortality in Australia (Bauman et al., 2022). The two most dominant factors of climate change, uneven precipitation distribution and increased temperature, are simultaneously the key factors affecting the distribution of P. massoniana in our predicted results. Under the background of global change, extreme events become frequently in the northern hemisphere, such as extreme low-temperature in the northern hemisphere in winter (especially the coldest month in the middle and high latitudes) and high temperature (mean temperature of warmest quarter increase) and drought (precipitation of driest month decrease) in summer (especially in subtropical China) (Voosen, 2021). Therefore, we have reasonable grounds to suspect that in the near future, the potentially suitable habitats for P. massoniana in China will be threatened (northward expansion limited by low temperature in the coldest month and southward expansion limited by high temperature and drought in summer). The potential distribution maps we simulated can serve as an important but relatively conservative reference for the protection strategy of P. massoniana, and can also protect its potential habitat more effectively (Carroll et al., 2017).
The present study supports three main recommendations for the protection and cultivation of P. massoniana. (1) Based on the potentially suitable habitat under the current climatic scenarios, P. massoniana should be protected by preserving its local hydrological and geological conditions. (2) Based on the predicted future distribution pattern of P. massoniana, we should strengthen monitoring and reduce human interference in the potential new habitats to protect the gradual migration of P. massoniana.
(3) Planting activities should be tailored to moderately suitable habitats and highly suitable habitats; specifically, considering the influence of different habitats on the growth potential of plantations, we recommended developing short-term and small-diameter wood plantings in moderately suitable habitats and cultivating highquality and large-diameter wood in highly suitable habitats.

Conclusion
In this study, the Maxent model was used to predict the potentially suitable habitat of P. massoniana in China based on 264 distribution records and nine environmental variables. Under the current climatic conditions, the potentially suitable habitat of P. massoniana was mainly distributed in the subtropical region south of the Qinling Mountains-Huaihe River Line. The main environmental variables that affect the potential distribution of P. massoniana were temperature and precipitation, including precipitation of the driest month, minimum temperature of the coldest month, annual temperature range, and mean temperature of the warmest quarter. Under future climate conditions, the area of potentially suitable habitat for P. massoniana was projected to increase, and the centroid of suitable habitat would thus move to higher latitudes. Under the SSP5-8.5 climate scenario, climate warming would increase due to higher greenhouse gas emissions, further promoting the migration of P. massoniana to even higher latitudes. Our study on the potential distribution of P. massoniana can serve as an important empirical reference for future conservation and planting of P. massoniana. Moreover, the methodology of this study can be applied to predict the potentially suitable habitats for other species worldwide.

Data availability statement
The original contributions presented in this study are included in the article/Supplementary material, further inquiries can be directed to the corresponding authors.

Author contributions
YC: writing-original draft, validation, investigation, and data curation. YHan and SJ: conceptualization, methodology, writingreview and editing, and funding acquisition. MZ: methodology. GW: writing-review and editing. PJ, YHu, PS, ZW, AF, and PQ: investigation and data curation. All authors contributed to the article and approved the submitted version.