Modeling by Maxent software. Randomly select 75% samples as training data set to model, and the remaining 25% distributed samples as test data set to verify the model. The maximum number of iterations is set to 10000, and the modeling is repeated for 10 times. Bootstrap is selected as the repetition type, and the distribution value is output in Logistic form. The accuracy of the model is evaluated by the value of AUC. The range of AUC value is 0. 5-1.0. The closer the AUC value is to 1, the more accurate the prediction is; 0. 5-0.7 means that the prediction effect is poor; 0. 8~0.9 means that the prediction effect is good; and 0.9-1.0 means that the prediction effect is very accurate (Wang et al., 2020).

The feature combination (FC) and regularization multiplier (RM) are adjusted by ENMeval data package in R software to optimize the model. There are five feature combinations of the Maxent model, namely Linear features, Quadratic features, Product features, Threshold features and Fragment features (Zhao et al., 2021a). In this study, the default parameters of Maxent software are RM = 1, FC =LQHPT; in order to optimize the Maxent model, RM is set to 0. 5 ~ 4, and every time it is increased by 0. 5, a total of 8 regulated frequency doubling is made. At the same time, six combinations of one or more characteristics are adopted: L, L and Q, H, L, Q, H and H, L, Q, H and P, L, Q, H, P and T. According to the permutation and combination, 48 parameter combinations are calculated. The combination of the above 48 parameters is tested, and the complexity of the model is tested according to the value of delta. AICc and the value of (auc.train-auc.diff.avg). The lower these two value are, the more accurate the prediction result of the model is.

Spatial units with a probability of species existence ≥ 0.374 are classified as suitable areas, while space units with a probability of species existence< 0.374 are classified as unsuitable areas. Using the reclassification tool in ArcGIS10.5 software binarize the potential geographical distribution data of A. rugosa under current and future climate change scenarios, and establish the presence/absence (0,1) matrix for the potential geographical distribution of A. rugosa. The suitable area is represented by a numeric value of “1” representing the existence of species, unsuitable area is represented by a numeric value of “0” representing the non-existence of species. Based on this matrix, further analysis was conducted on the spatial pattern changes of suitable areas of A. rugosa under current and future climate scenarios, and three types of changes in suitable areas were defined: newly added suitable areas, retreated suitable areas, and stable suitable areas. The spatial pattern changes of potential suitable areas under current and future climate change are defined as follows: matrix values 0 → 1 represents newly added suitable growth areas, 1 → 0 represents retreated suitable areas, and 1 → 1 represents stable suitable areas (Zhao et al., 2021b; Yao et al., 2023).

In this study, the niche of A. rugosa in different periods and different climate modes was quantitatively analyzed. Based on the environmental factor layer and the distribution layer of A. rugosa in different periods and different climate modes, the niche was analyzed and visualized by using the “ecospat” package in R language (Yan et al., 2021), and the niche overlap index (percentage overlap index) was calculated, which was represented by the letter D, and its range was 0-1, indicating that the niche was from no overlap to complete overlap (Dan et al., 2008). Using the Niche breadth module in ENMTool, the niche width was calculated based on the potential distribution data onto A. rugosa in the current and future periods, with B1 representing the minimum niche width and B2 representing the maximum niche width.

Marxan model can be used to select the minimum cost area in the protection system. With the continuous improvement in the model, it has been widely used in the planning of land protection system (Xie et al., 2022). Taking the potential distribution area of A. rugosa in the current climate as the research object, the Marxan model was used to calculate the priority protected areas of A. rugosa in the current climate. Square Pu was used in the analysis, height and width were set to 25000. 1km² as the research unit, using the Zonal Statistics as Table tool for ArcGIS, the distribution area of target species in each planning unit is counted, and the species distribution matrix is constructed. The protection target is set to 30% of the total habitat area and the SPF value is 100. (Li et al., 2019; Zhang et al., 2011), The boundary length modifier (BLM) of the model is the correction parameter of the boundary length of the protected area. By modifying the BLM, we can analyze the relationship between the cost and the total length and total area of the boundary, thus finding a balance point, and obtaining a more reasonable spatial distribution pattern of the protected area through calculation(Mou et al., 2021). The final model uses the model boundary modifier of 25,000. The model is iterated for 100 times to get the optimal solution to the planning unit.

Using ArcGIS V10.5 software, the prediction results of Maxent model are transformed into raster data, and the values of the raster represent the survival probability of A. rugosa in this area. The natural discontinuous breakpoint grading method was used to divide the suitable areas into four grades: high suitable area (P≥0.7), moderate suitable area (0.57≤P< 0.7), low suitable area (0.374≤P<0.57) and unsuitable area (P<0.374). The data onto priority protected areas calculated by Marxan model are imported from ArcGIS, and the visualization of priority protected areas is realized by overlaying with the current distribution layer.

Using SDMTool data packet of R language, the centroid position of the suitable area of A. rugosa in current and future climate scenarios is calculated, and the migration direction of the spatial distribution of the suitable area of A. rugosa is reflected by the change of centroid position. Geosphere data packet of R language is used to calculate the centroid migration distance of A. rugosa in different climate scenarios.

Based on 390 distribution points of A. rugosa and 12 environmental factors selected for modeling, the potential distribution area of A. rugosa was simulated by Maxent model. When the model is the default parameter, ΔAICc =11.88; when the model parameters are set to FC = L and RM = 0.5, ΔAICc = 0 ( Table 2 ), which is the optimal model. Therefore, FC = L and RM = 0.5 are set as modeling parameters. Using the optimized parameters, the model was reconstructed to simulate the suitable area of A. rugosa, and the AUC value of simulated training for this parameter was 0.890 ( Figure 2 ), which indicated that the prediction result was accurate.

Using ArcGIS software to classify the current distribution area of A. rugosa, and get the distribution map under the current climate conditions ( Figure 3 ). As shown in Figure 3 , the distribution area of A. rugosa is mainly concentrated on the eastern and central regions of China, covering Liaoning, Shanxi, Henan, Sichuan, Yunnan, Guizhou, Hunan, Hubei and other provinces. Among them, the highly suitable areas are mainly located in northeastern Sichuan, southwestern Shaanxi and western Hubei, while the moderately suitable areas are mainly located in eastern Sichuan and Guizhou. The simulation and prediction results are basically consistent with the geographical distribution data of A. rugosa, which further shows that the simulation and prediction results are accurate to some extent. According to the number of grids occupied by the suitable area, the area occupied by the suitable area can be calculated ( Table 3 ). According to Table 3 , the suitable area of A. rugosa in the current climate is 1,761,700 km², accounting for 18.4% of China’s total land area, of which the highly suitable area is 138,700 km², accounting for 4.7% of China’s total land area. It can be seen that the highly suitable area of A. rugosa in China is small and the growth range is concentrated.

The forecast result of Maxent model is converted into raster data, and the suitable area of A. rugosa under different climate scenarios is calculated. it can be seen that the suitable area of A. rugosa has increased compared with the current climate scenarios, but the growth rate is different, with the smallest growth rate under SSP1-2. 6, followed by SSP2-4.5 and SSP5-8.5 ( Table 3 ). Taking 2090s as an example, it increased by 6.6%, 9.1% and 15.9% respectively.

Using ArcGIS software to draw the suitable zoning map of A. rugosa in the future climate ( Figure 4 ), by comparing the predicted suitable zones in 1950s and 1990s, we can see that there are differences in the response to A. rugosa to climate change in different scenarios. The change amplitude is the smallest in the SSP1-2.6 scenario, with an area of 239.4km² in 2050s and 239.8km² in 2090s, which is basically unchanged. However, under the climate scenarios of SSP2-4.5 and SSP5-8.5, the suitable area of A. rugosa increased greatly, increasing by 20 km² and 39 km² respectively. It can be seen that under the SSP1-2.6 scenario, the response to the suitable area of A. rugosa to climate change is not obvious, while under the SSP2-4.5 and SSP5-8.5 climate scenarios, the response of the suitable area of A. rugosa to climate change is more sensitive.

In terms of spatial pattern, there are some differences in the migration positions of the suitable areas of A. rugosa under different climate scenarios, but the overall migration trend is consistent, and it generally migrates to the northeast ( Figure 5 ). At present, the center of mass of the suitable area of A. rugosa is located in Zigui County, Yichang City, Hubei Province (110.75°E, 30.70°N). When the climate scenario is SSP1-2.6-2090s, the center of mass of the suitable area of A. rugosa moves to the northeast, while the center of mass of the suitable area of A. rugosa is located in Lushan County, Pingdingshan City, Henan Province (112.8°E, 33.60°N) When the climate scenario is ssp 2-4.5-2090s, the center of mass of the suitable area moves to the northeast. At this time, the center of mass of the suitable area of A. rugosa is located in xinmi city, Zhengzhou City, Henan Province (113.50°E, 34.42°N), with a migration distance of 485,880 m; When the climate scenario is SSP5-8.5-2090s, the center of mass of the suitable area moves to the northeast. At this time, the center of mass of the suitable area of A. rugosa is located in Anyang County, Anyang City, Henan Province (114.20°E, 36.25°N), and the migration distance is 693,417 m. It can be seen that under the future climate change scenario, global warming and humidification will make the center of mass of the suitable area of A. rugosa in China move to the northeast as a whole, and the migration position will further expand upon the north.

The distribution patterns of A. rugosa in different periods under different climate scenarios were compared and analyzed with the current period ( Table 4 and Figure 6 ), and the results showed that the suitable areas of A. rugosa in the future period and climate scenario were mainly increased. Among them, the expansion rate, retreat rate and stability rate of suitable area remained the highest in the scenario of 2081-2100 –SSP5-8.5: the area from unsuitable area to suitable area was 175.14km², and the expansion rate was 99.42%. The area from suitable area to unsuitable area is 22.77km², and the lost rate is 12.93%. The unchanged area is 153.4km², and the stability rate is 86.49%. As can be seen from Figure 6 , the expanded suitable areas are mainly concentrated in the northeast of China, such as Jilin, Hebei, Shanxi, Shaanxi and other provinces. The lost areas are mainly concentrated in the southeast of China, such as Fujian, Jiangxi, Hunan and other provinces. Generally speaking, the suitable area of A. rugosa changes with climate change, and the distribution pattern of A. rugosa in different periods under different climate scenarios has basically the same response to climate change. Among them, the change of the suitable area of A. rugosa is the most significant under the climate scenario of SSP5-8.5, and the change of the suitable area of A. rugosa is the least obvious under the climate scenario of SSP1-2.6.

The calculated habitat suitability and environmental factors are analyzed by knife cutting method. From the table of contribution rates of various environmental factors of the distribution of A. rugosa ( Table 5 ), the top five environmental factors are the wettest monthly precipitation (59.1%), the lowest temperature in the coldest month (11.9%), altitude (9%), slope variability (6.6%) and the highest temperature in the hottest month (4.7). The top five environmental factors are the wettest monthly precipitation (42.5%), the lowest temperature in the coldest month (14.3%), the slope variability (9.8%), the highest temperature in the hottest month (9.2%) and the altitude (7.1%), accounting for 82.9% in total. Considering the contribution rate and important value, the dominant environmental factors restricting the distribution of A. rugosa are precipitation, temperature, altitude and slope.

The suitable area of species mainly depends on the dominant environmental factors. We can clearly understand the relationship between the distribution of A. rugosa and environmental factors by modeling 12 environmental factors and drawing the single factor response curve ( Figure 7 ). It is generally believed that when the survival probability is greater than 0.5, the corresponding environmental factors are suitable for plant growth. It can be seen from Figure 4 that when the wettest monthly precipitation (bio13) is less than 100ml, the survival probability of A. rugosa is less than 10%. After that, with the increase of precipitation, the survival probability of A. rugosa also increases. When the precipitation reaches 210ml, the survival probability reaches the peak; when the precipitation exceeds 210ml, the survival probability of A. rugosa decreases with the increase of precipitation. When the precipitation is between 160 ml and 300 ml, the survival probability of A. rugosa is greater than 0.5. Similarly, similar to the wettest monthly precipitation, there is an optimum range for the highest temperature in the hottest month, the lowest temperature in the coldest month, altitude and slope variability, and the survival probability of A. rugosa is less than 0.5 if it is too large or too small. Therefore, it can be inferred that under the current climate conditions, the wettest monthly precipitation in the suitable area of A. rugosa is 160ml-300ml, the highest temperature in the hottest month is 24°C-32°C, the lowest temperature in The coldest month is -5°C-7°C, the altitude is 100m-1200m, and the slope variability is 0.2×10⁶-1.3×10⁶. This is basically consistent with the environmental characteristics of the current potential distribution area of A. rugosa calculated in the result analysis of main environmental factors ( Table 6 ), indicating that the prediction results are accurate.

From the results analysis of main environmental factors ( Table 6 ), it can be seen that the maximum temperature in the hottest month decreases from the passage of time for the SSP1-2.6 scenario, while the maximum temperature in the hottest month increases in the passage of time for the other two climate scenarios. The lowest temperature in the coldest month, the wettest monthly precipitation and the habitat suitability of species all show an increase in time for three climate scenarios. However, the altitude and slope variability remain basically unchanged. It can be inferred that in the future climate scenario, the temperature in the suitable area of A. rugosa will increase compared with the current overall, and the precipitation will increase.

The niche space of A. rugosa under different climate models was visualized. As shown in Figure 8 , the niche overlap between A. rugosa in different climatic backgrounds is relatively large (the maximum niche overlap D_50126 = 0.837; The minimum niche overlap D_90585 = 0.666), so there is no significant niche differentiation. In addition, in the same period, with the increase of radiation intensity, the niche overlap of A. rugosa tended to decrease gradually, among which the niche overlap in the climate scenario of SSP5-8.5 was significantly lower than that in the other two climate scenarios. Compared with 2050s, the niche overlap of A. rugosa in 2090s decreased to different degrees under each radiation intensity. This shows that the resources that can be used together in the future and the current period are reduced.

ENMTools software package was used to calculate the niche width of A. rugosa in different climate backgrounds. As shown in Table 1 supplementary materials, the maximum value of B1 is 0.557 and the minimum value is 0.436; the maximum value of B2 is 0.967, and the minimum value is 0.951. Therefore, there is no obvious difference between B1 and B2 in each period, which indicates that A. rugosa is more inclined to be a generalized species. In addition, compared with the current period, B1 and B2 have increased in other climate scenarios, which indicates that all kinds of resources that A. rugosa can use in the future climate scenarios have increased and are widely distributed, and have strong adaptability to the future environment.

The priority protected areas of A. rugosa were calculated by marxan model, and the results were imported from ArcGIS software to generate a system protection plan of A. rugosa as the main protection target. As shown in Figure 9 , the priority protected areas of this species are concentrated on Chongqing, eastern Sichuan, southern Guizhou, western Hunan and Hubei, and southwestern Shaanxi, which is basically consistent with the highly suitable areas of A. rugosa predicted by Maxent model, indicating that the prediction results are accurate. In addition, the priority protected area of A. rugosa occupies a small proportion of the land area, and the distribution is concentrated, which is conducive to the formulation of targeted protection and management.

Based on environmental variables such as climate, topography and soil factors, this study applies ENMeval data package to optimize the model. This method limits the background data onto the area corresponding to the calibration position, so that the potential geographical distribution area simulated by Maxent covers the current distribution point. This method allows the model parameters to be adjusted to improve the performance of the Maxent model, and its accuracy can be measured by improving the fitting degree between the prediction results and the actual distribution area and by visual inspection of the geographical prediction map (Guo et al., 2018). The Maxent model with optimized parameters can effectively reduce the complexity of the model, improve the fitting degree between the predicted results and the actual situation, predict the species distribution effectively, and the response curve obviously becomes smooth, and it is close to the normal distribution curve, which conforms to Shelford’s tolerance law (Phillips et al., 2017; Li et al., 2018; Ouyang et al., 2019). Many related studies have confirmed that the concentration of sample size and distribution points of species distribution data will significantly affect the accuracy of species model simulation results. In general, with the increase of sample size, the simulation accuracy of species distribution model increases, and the increasing range gradually decreases until it no longer increases, and finally tends to reach the maximum accuracy of the model (Chen et al., 2012). Too concentrated distribution points will increase the over-fitting of environmental deviation caused by spatial autocorrelation of distribution points, which will have a certain impact on the simulation results of the model (Yao et al., 2023). In this study, a total of 505 distribution points of A. rugosa were collected, and the correlation between the distribution point data and environmental variables was analyzed and screened, and finally 390 distribution point data were obtained, which solved the inaccurate modeling results caused by the small sample size and strong multicollinearity among environmental factors. Some studies show that the change of research scale will lead to the change of background data, and then affect the model construction (Ji et al., 2019; Ma and Li, 2023). In this study, the potential distribution area of A. rugosa is predicted at the national scale. If the study is conducted at a larger or smaller scale, the selection of background data points should be effectively compared with the occurrence points of species, rather than simply randomly extracting background points from the research scope. However, how to choose a reasonable background data point is still worth further study. At present, there is no unified standard for the classification of the grade in the suitable area, and the selection of the threshold directly affects the area of the suitable area, which may lead to a big difference between the divided range of the highly suitable area and the actual distribution range, so it is necessary to select an appropriate threshold for the classification and compare it with the actual distribution range (Wang et al., 2020).

Under the current climate conditions, the main environmental factors affecting the distribution of A. rugosa are precipitation and temperature, followed by topographic factors, and soil factors have the lowest influence on A. rugosa. Based on the response curves of various factors, the wettest monthly precipitation in the suitable area of A. rugosa can reach 300ml, the highest temperature can reach 32°C and the lowest temperature can reach -5°C, which is consistent with the biological characteristics of A. rugosa that it likes high temperature, humidity and cold resistance of roots. Although hydrothermal conditions play a major role in the potential geographical distribution pattern of A. rugosa in China, the constraints of topographic factors and soil factors cannot be ignored. The contribution rate of environmental factors shows that many factors jointly affect the potential geographical distribution of A. rugosa, such as temperature, moisture, altitude, slope and so on.

Under different climate scenarios in the future, the potential distribution of A. rugosa in China is quite different from that at present, indicating that future climate warming will have a certain impact on the geographical distribution of A. rugosa. Compared with the current climate scenario, the area and spatial pattern of the highly suitable area of A. rugosa have changed significantly. In 2050s, the area of highly suitable areas under the three climate scenarios did not change much, but in 2090s, the area of highly suitable areas under the three climate scenarios showed an increasing trend in turn, indicating that the growth rate of the distribution area of A. rugosa accelerated with the passage of time, so it can be inferred that the global warming situation is also gradually intensifying. The response of the spatial pattern of the suitable area to climate change is consistent on the whole, that is, with the intensification of climate warming, the overall migration range of the spatial position of the suitable area of A. rugosa becomes larger. The migration trend of A. rugosa is consistent with the geographical distribution of temperate tree species migrating to high latitudes under future global warming (Zhang et al., 2019). The size of niche overlap reflects the similarity of plant utilization of resources. The large niche overlap shows that they have similar ecological requirements, resource utilization and biological characteristics under certain circumstances (Yuan et al., 2021). The decrease of resources that can be used together in the future and the current period indicates that the geographical distribution area of A. rugosa has changed under the influence of global climate change, which is consistent with the overall migration of the center of mass of the suitable area of A. rugosa predicted by Maxent model to the northeast. Compared with the current period, the niche width of A. rugosa in the future climate scenario has increased, indicating that the area of its suitable area has increased under the influence of global change, which is consistent with the prediction results of Maxent model and proves the accuracy of Maxent model. The prediction of the future suitable area of A. rugosa is based on the assumption that the soil and topographic factors have not changed for 70 years(Zhang et al., 2018). The research results are only valid at the national and provincial levels, and cannot be applied to local microclimate. Agastache species all have similar morphological and physiological characteristics, and its geographical distribution may also have some similarities in response to climate change, so this model has certain guiding significance for the prediction of distribution areas of other species of Agastache in future climate scenarios.

From the perspective of nature reserves, the systematic protection zoning for A. rugosa is explored. The results show that the priority protected areas are concentrated on Chongqing, eastern Sichuan, southern Guizhou, western Hunan and Hubei, and southwestern Shaanxi. These areas are mostly warm and humid mountainous areas with high altitude and terrain, good natural conditions such as light, moisture and temperature, which are consistent with the growth habits of most medicinal plants and are very suitable for the growth of Chinese herbal medicines (Liu et al., 2021). The establishment of protected areas in mountainous areas has less human interference and is more conducive to plant protection. At the same time, the priority protected areas suitable for the growth of A. rugosa have certain guiding significance for the industrialization development of A. rugosa in the future, and the development of A. rugosa planting base can also be taken as the primary task to promote the local economic development.