Modeling the potential distribution of different types of Dendrocalamus sinicus, the strongest woody bamboo in the world, with MaxEnt model

Climate change will significantly affect the distribution area of species. Through establishing distribution model, we can simulate the current and future potential distribution range and provide reference for the introduction and cultivation planning of rare or economic plants. Dendrocalamus sinicus, endemic to Yunnan Province of China, is the strongest woody bamboo in the world. In the present study, the MaxEnt model was performed to simulate the distribution of different types of D. sinicus in China and neighboring countries or regions. The results suggested that the suitable distribution range of “straight type”, the main type for cultivation and utilization, was 8°–30°N and 73°–122°E under the current climate conditions, while the potential distribution range of “bending type” was 6°–31°N and 79°–109°E. The two most key climate variables associated with distribution of “straight type” were “Temperature Annual Range” with 36.6% contribution rate and “Temperature Seasonality” (32.4%), while “Isothermality” (47.8%) and “Precipitation of Driest Month” (24.8%) for “bending type”. Under different climate change scenarios (SSP1-2.6, SSP5-8.5) and periods (2050, 2090), the potential distribution area of the “straight type” were apparently different, indicating that the distribution area of D. sinicus will be affected significantly by climate changes in the future. Our findings would be not only beneficial to understanding limiting factors for natural distribution of D. sinicus, but also helpful for further germplasm conservation, introduction and cultivation planning of this rare woody bamboo.


INTRODUCTION
Nowadays, global warming and biodiversity conservation are some of the most crucial challenges to human sustainable development, and the impacts of climate change on species survival have become hot-spot issue accordingly (Harrington, Fleming & Woiwod, 2001;Kriegler et al., 2017;Wu et al., 2018). In general, the global climate change significantly affects the distribution range of species, especially for species sensitive to climate factors culms at the lower half (Fig. 1B), and its internodes are swollen, shortened or deformed, which is a little value in garden landscaping. The significant differences in ecological conditions were detected between the habits of two types (Hui, Yang & Du, 2006). The ''straight type'' usually distributes at subtropical mountains with 1,000-1,500 m above sea level, while the ''bending type'' occurs in habitats of marginal tropical montane regions at elevations of ca. 500-1,000 m (Gu et al., 2012;Yang et al., 2018). So far, due to the integrative influence of climate change and biological characteristics, e.g., sporadic flowering and low seed setting rate, the habitat of D. sinicus is scattered and high-quality germplasm resources are scarce (Gu et al., 2012;Xie et al., 2019). Therefore, it is necessary to choose suitable areas to protect and propagate this rare bamboo species. In the present study, we use the MaxEnt model to screen the key environmental factors affecting the distribution of D. sinicus, and to predict potential distribution areas. Our aims are to distinguish out the key natural factors limiting the distribution of D. sinicus, and to provide a scientific basis for its germplasm conservation and further introduction and cultivation planning.

Species data sources
By searching the Chinese Virtual Herbarium (https://www.cvh.ac.cn/), National Science and Technology Infrastructure (http://www.nsii.org.cn/2017/home.php), Flora Reipublicae Popularis Sinicae (http://www.iplant.cn/frps), Global Biodiversity Information Facility (http://www.gbif.org) as well as our investigation and sample collection in the field, a total of 25 distribution locations of D. sinicus were obtained (Table 1) after checking the integrity of coordinate information and eliminating the duplicate coordinates. In the previous molecular genetic study (Yang et al., 2018), all existing populations of D. sinicus were divided into two genotypes, i.e., ''straight type'' and ''bending type'', which was consistent with the morphological characters of culms (Gu et al., 2012). Therefore, according to the previous results, 25 coordinate points were divided into 13 ''straight type'', eight ''bending type'' and four mixed coordinate points. All coordinates were saved in the CSV (comma-separated values) format for MaxEnt model analysis according to culm type. Software ArcGIS 10.2 was used to reproduce the current distribution points of D. sinicus (Fig. 2).

Climate variables
In this study, bioclimatic variables from the WorldClim website (https://www.worldclim. org/data/worldclim21.html#google_vignette) (Table S1) was used in further analysis. We removed variables 8, 9, 18, and 19 because of spatial artifacts (Ashraf et al., 2017). The Yu n n a n , C h i n a L a n c a n g -Me k o n g R i v e r My a n ma r L a o s environmental variables, recorded from 1970 to 2000, were saved in file format of ASCII (American Standard Code for Information Interchange) converted by ArcGIS 10.2, and were used to establish the initial model (Fick & Hijmans, 2017). In order to eliminate the auto-correlation and colinearity between variables, the Pearson correlation coefficient (r) and principal component analysis (PCA) of 15 bioclimate variables from occurrence records were tested, and only one variable from each set of highly cross-correlated variables (|r| > 0.8) was kept for further analysis according to Yi et al. (2016). For instance, the variables bio6 was correlated with both bio5 (r = 0.8) and bio7 (r = − 0.9), then bio6 was dropped and both bio5 and bio7 were reserved based on the PCA result. Future climate projections were extracted from the Shared Socioeconomic Pathways (SSPs) of the BCC-CSM2-MR global climate (2021-2040, 2041-2060, 2061-2080, and 2081-2100) database in the Coupled Model Intercomparison Projects 6 (CMIP6) (https://www.worldclim.org/data/cmip6/cmip6_clim2.5m.html). The associated scenarios included the SSP1-2.6 (ssp126), SSP2-4.5 (ssp245), SSP3-7.0 (ssp370), and SSP5-8.5 (ssp585), following Liu et al. (2021a). The SSPs reflect four different developments of the world that were characterized by varying levels of global challenges. We selected SSP1-2.6 (low forcing scenario, radiation intensity reaches 2.6w/m 2 in 2100) and SSP5-8.5 (high forcing scenario, radiation intensity reaches 8.5 w/m 2 in 2100) as future test scenarios according to Zhang, Chen & Xin (2019), and simulated the suitable distribution of ''straight type'' in the 2050s (2041-2060) and 2090s (2081-2100). All environmental data used in this model were 2.5-arc minute spatial resolution (also referred to as 4.5 km spatial resolution).

Construction and test of the MaxEnt model
For each type, we created 493 candidate models by combining 17 values of regularization multiplier (0.1-1.0 at intervals of 0.1, 2-6 at intervals of 1, as well as 8 and 10), and all 29 possible combinations of five feature classes (linear = L, quadratic = Q, product = P, threshold = T, and hinge = H). The performance of the candidate model was evaluated on the basis of significance (partial ROC, with 500 iterations and 50 percent of data for bootstrapping), omission rates (E = 5%), and model complexity (Akaike Information Criterion corrected for small sample sizes, AICc). The best models were selected according to the significant models and the omission rate ≤ 5%. The models of Delta AICc ≤ 2 were selected as final models in this model set (Cobos et al., 2019;Alkishe & Peterson, 2022). The complete set of occurrences and the selected parameterizations were used to create the final models for the two types. We produced 10 replicates by bootstrap, with logistic outputs, and transferred these models to the 12 countries or regions for current and future scenarios. Partial ROC and omission rates (E = 5%) were evaluated for the final models (Cobos et al., 2019). The above analysis was performed through the kuenm R package (https://github.com/marlonecobos/kuenm). The receiver operating characteristic curve (ROC) of all environmental variables was calculated, and the importance of different environmental variables was measured by Jackknife test, following Li et al. (2019a). The areas under the receiver operating characteristic curve (AUC) were used to evaluate the accuracy of the model (Goldenberg, Reginato & Michelangeli, 2020). The evaluation criteria are as follows: 0.5 < AUC ≤ 0.6 (the model is failed), 0.6 < AUC ≤ 0.7 (''poor''), 0.7 < AUC ≤ 0.8 (''general''), 0.8 < AUC ≤ 0.9 (''better'') and 0.9 < AUC ≤ 1 (''excellent''). The closer the AUC value was to 1, the higher was the prediction accuracy of the model (Wang et al., 2018;Phillips et al., 2009). The suitability maps were calculated using the logistic output of MaxEnt, which ranges from 0 (lowest suitability) to 1 (highest suitability). For visualization and further analysis, the prediction results were imported into ArcGIS 10.2 and divided the habitat suitability maps into four levels according to expert experience and relevant literature (Hui, Yang & Du, 2006;Zhang et al., 2019): unsuitable habitat (0-0.35), poorly suitable habitat (0.35-0.55), moderately suitable habitat (0.55-0.75), and highly suitable habitat (0.75-1).

Model selection and accuracy evaluation
After kuenm R package screening, the feature class was LP and regularization multiplier was 2 in the final model parameter combination of the ''straight type''. As for the ''bending type'', the feature class was LQTH and regularization multiplier was 3. We rebuild the model using the optimized parameter combination. The AUC values of the ''straight type'' and ''bending type'' were 0.978 and 0.996 ( Fig. S1) respectively, indicating both models were ''excellent''. The above results suggested that the distribution area simulation using the MaxEnt model was reliable, and we further analyzed the impact of climate change on the distribution area of D. sinicus.  (Table 2), with a cumulative contribution rate of 97.8%. The results of Pearson correlation coefficient (r) analysis indicated that none of correlation coefficient between the environmental variables exceeded 0.8 (Table 3). Because altitude and solar radiation were highly correlated with many climate variables (Tables S2, S3), climate variables were preferentially selected as variables for the classification of suitable areas, rather than solar radiation and altitude. Finally, all four climate variables above were deemed as the main factors affecting the distribution of two types. On this basis, the MaxEnt models of the ''straight type'' and ''bending type'' distributions were established and were further evaluated.

Relationship between geographical distribution and environmental variables
The importance values of key climate variables to the distribution of D. sinicus were analyzed using the Jackknife test in the MaxEnt model. Within the ''straight type'', the order of the four key climate factors was Temperature Seasonality (bio4), Temperature Annual Range (bio7), Annual Precipitation (bio12) and Precipitation of Wettest Month (bio13) (Fig. S2), based on importance. And their appropriate variation ranges were <390, <23.9 • C, >1,418 mm and > 282 mm, respectively (Table 4). As for the ''bending type'', the order of four key climate factors was Temperature Seasonality (bio4), Isothermality (bio3), Mean Temperature of Coldest Quarter (bio11) and Precipitation of Driest Month (bio14) (Fig. S2), and their appropriate variation ranges (optimal values) were 210-410   (Table 4). The response curves of environmental variables to the distribution probability were shown in Figs. S3 and S4. When the environmental variables were lower than the optimal value, the distribution probability increased with the increase of the environmental variables. And it worked in the reverse too.

The suitable areas of different types under current climatic scenario
MaxEnt model predicted the potentially suitable areas of D. sinicus included 12 countries or regions under the current environmental conditions. The potential distribution range of the ''straight type'' was between 8 • -30 • N and 73 • -122 • E (Fig. 3A), and the distribution area covered 2234583.37 km 2 (Table 5). Among potential distribution range, the countries with large distribution areas were India (30.86%), China (24.61%), Myanmar (16.72%), and Thailand (15.86%) ( Table 5). In addition, the highly suitable areas, with a cumulative area of 133,072.92 km 2 , were mainly distributed in Myanmar (81.57%), Thailand (12.82%), and Bangladesh (3.94%) ( Table 5). By and large, the suitable areas of the ''straight type'' extended along the longitude and latitude from the existing distribution range, and the extended range along the longitude was much larger than that of the latitude. As for the ''bending type'', its potential distribution range was between 6 • -31 • N and 79 • -109 • E (Fig. 3B), and the highly suitable areas were mainly distributed in southwest China, northern Laos and eastern Myanmar. Different from the ''straight type'', the distribution range of the ''bending type'' mainly displayed a trend of southward expansion along low latitude to the tropics. It was worth noting that the potential distribution areas of the two variants were significantly larger than the existing ranges.

The suitable distribution area of the "straight type" under future climate scenarios
As an excellent timber bamboo species and main type for industrial utilization and cultivation, the ''straight type'' is our concerning focus, and in this study, we predicted the  potentially suitable distribution area of the ''straight type'' under SSP1-2.6 and SSP5-8.5 future climate change scenarios (Fig. 4). Under the SSP1-2.6 scenario in 2050, the highly suitable habitat of the ''straight type'' was mainly in Myanmar and western Thailand. While under SSP5-8.5 in 2050 and SSP1-2.6 in 2090, the highly suitable habitat basically disappeared in western Thailand, but emerged in parts of eastern India. In addition, under SSP5-8.5 in 2090, the suitable habitat, especially for the highly suitable distribution, sharply shrank in eastern India, southern China, Cambodia, Thailand and Myanmar. Compared with the current distribution, three types of the predicted suitable habitats, namely poorly, moderately and total highly suitable regions, exhibited different change trends under the climate change scenarios SSP1-2.6 and SSP5-8.5 in 2050. The total highly suitable regions displayed a decreasing trend, namely 27.42% and 8.19% under SSP1-2.6 and SSP5-8.5 scenarios respectively (Table 6). But the poorly, moderately and total suitable habitat had an increasing trend, and the increase rate of suitable area under low radiation intensity (SSP1-2.6) was higher than that of high radiation intensity (SSP5-8.5).
On the other hand, comparison with the current distribution, almost predicted suitable habitats under two climate scenarios in 2090 exhibited dramatically decreasing trends, except the highly suitable area slightly increased 0.59% under the SSP1-2.6 climate scenario in 2090 (Table 6). Meanwhile, the high radiation intensity scenarios (SSP5-8.5) had a much stronger effect on decreasing suitable habitat than that of low radiation intensity scenarios (SSP1-2.6). The reduction of suitable areas under SSP5-8.5 in 2090, to a large extent, were due to the exclusion of eastern India, southern China (GX, HI, GD, FJ and JX), Cambodia and Thailand (Fig. 4).

Key climate factors determining the survival and distribution of D. sinicus
As the strongest woody bamboo documented in the world, D. sinicus has outstanding performance in producing timber. Meanwhile, because of narrow distribution and provenance scarcity, its germplasm conservation and introduction planning have raised extensive attentions since publication as a new species (Chia & Sun, 1982;Hui, Yang & Du, 2006;Gu et al., 2012). The consequent problem is to clarify the key climate factors determining the survival and distribution of D. sinicus, which is important but yet controversial. Empirically, low temperature in winter (Hui, Yang & Du, 2006) or precipitation (Pu, 2004) were once considered to be the most key climate variable determining the growth and distribution of D. sinicus. In this study, the MaxEnt model detected that the bioclimatic variables of the highest contribution rate to the distribution of two types were different, namely Temperature Annual Range (36.6%) for the ''straight type'' and Isothermality (47.8%) for the ''bending type'' respectively (Table 2), indicating that temperature was probably the most key factor affecting the survival of D. sinicus. This result was similar to opinion of Hui, Yang & Du (2006). It may be due to the biological characteristics of D. sinicus, namely, it is a tropical bamboo species occurred at the edge of the tropics, and the low temperature in winter will seriously threaten its survival (Chia & Sun, 1982;Hui, Yang & Du, 2006).

Suitable distribution range based on the MaxEnt model and introduction practice of D. sinicus
According to the prediction of the MaxEnt model, the current suitable distribution ranges of the two types were different. For the ''straight type'', the suitable distribution range was 8 • -30 • N and 73 • -122 • E with a broad span of east-west distribution, while the potential distribution range of ''bending type'' was 6 • -31 • N and 79 • -109 • E. Compared with the ''straight type'', the highly suitable habitat of ''bending type'' appeared sporadically in south Vietnam, south India and Sri Lanka, This suggested that the distribution range of ''bending type'' exhibited a southward spreading trend to the tropical area, which was similar to the previous research results (Yang et al., 2018). The identical environmental variable affecting the distributions of the two types was Temperature Seasonality, which is the most important variables in the Jackknife test, and the maximum value of Temperature Seasonality in the ''bending type'' was higher than that in the ''straight type'' (Table 4). This might reflect that the ''bending type'' could survive in regions with larger temperature seasonality variation, and had more tropical attribute than ''straight type'' (Chia & Sun, 1982;Hui, Yang & Du, 2006;Yang et al., 2018).
On the aspect of introduction and cultivation planning of D. sinicus, our results also provided some new views. Conventionally, planning for plant introduction and cultivation range was mainly based on similar climate, soil conditions and other factors between the introduction area and the origin habitat (Hui, Yang & Du, 2006). The previous study inferred that the climate conditions of the most suitable area for the ''straight type'' were as follows: (1) the average temperature in the coldest month was ≥13 • C; (2) the days of the daily minimum temperature ≤ 0 • C were not more than 1 day; (3) frost did not happen throughout the year; and (4) the annual precipitation was ≥ 1200 mm (Hui, Yang & Du, 2006). From 2000 to 2003, the scientists introduced and cultivated ''straight type'' at 12 predicted suitable areas in eight counties of Yunnan Province (Hui, Yang & Du, 2006;An & Chen, 2010). After 15 years, the introduced bamboo clumps of D. sinicus survived only at two sites: Ning'er County (23 • 32 N) and Xinping County (24 • 04 N), in which clumps could bear new shoots normally and the culms grew up to 20 cm in diameter. Moreover, approximate 10% introduced clumps bloomed and died in 1-2 years after cultivating. The above results suggested the strict conditions for the introduction and cultivation of D. sinicus. On the other hand, based on the MaxEnt model, we predicted the optimal climate conditions of highly suitable areas for ''straight type'' were: (1) the Temperature Annual Range was <23.9 • C; (2) the Precipitation of Wettest Month was >282 mm and (3) the Annual Precipitation was >1,418 mm, respectively. Compared with the empirical climate variables for highly suitable areas of ''straight type'', our climate conditions based on MaxEnt model were more strict. However, our results were also more consistent with the actual outcomes, indicating our results had a higher reliability.

Change trend on suitable habitat of "straight type" under different climate change scenarios
The adaptation to climate and climate change is vital to plant growth, geographical distribution and biodiversity (Liu et al., 2022;Lv & Wu, 2009;Yi et al., 2017). In the present study, the potential suitable distribution area of ''straight type'', predicted by MaxEnt model, were different under two climate change scenarios (Table 6). Compared with the current distribution, the total suitable habitat of ''straight type'' would increase in 2050 and decrease in 2090 under same radiation intensity. Furthermore, the total suitable habitat would increase slightly (2050) and decrease dramatically (2090) under high radiation intensity (SSP5−8.5) ( Table 6), implying that higher radiation intensity (SSP5−8.5) would limit distribution of the ''straight type''. The previous studies indicated that global warming would increase, fluctuate or decrease the distribution range of species (Thomas et al., 2004;Yuan, Wei & Wang, 2015). Our results also suggested that the impact of climate change on plant distribution might be a long-term process. The slow increase of radiation intensity in a short period would not give rise to significant changes in distribution area of species. This is because plants possess certain capacity of self-regulation and diffusion (Liu et al., 2021b), which result in fluctuated suitable distribution area. In extreme cases, for example, excessive radiation intensity will affect the growth and development of plants, and eventually lead to death. In turn, the distribution area of plants will reduce or disappear (Li, Fan & He, 2020).
In summary, our results confirmed that the ''straight type'' had poorer heat resistance than ''bending type'', which is important to instruct future introduction and cultivation planning of D. sinicus. In addition, we should also realize that the MaxEnt model does not consider the factors such as soil and soil microorganisms at specific sites, so the accuracy of this prediction is limited.

CONCLUSIONS
Based on D. sinicus occurrence records and bioclimatic variables, the current and future suitable habitat of D. sinicus in China and adjacent regions was modeled using MaxEnt model for the first time. The bioclimatic variables of temperature annual range and isothermality were revealed to have crucial effect on D. sinicus distribution. In the next 70 years, the habitat suitability of this woody bamboo may be different with climate change. The prediction of this study is of strategic significance for further germplasm conservation, introduction and cultivation planning of this rare and precious bamboo species.