Human influence on the earth's climate has become more and more evident (Lindner et al., 2010). Climate observations have clearly shown that average global temperatures have increased by 0.8 ℃ since 1900 (http://www.sciencedirect.com/science/article/pii/S0378112709006604, Hansen et al., 2005; Hansen et al., 2010). Climate change can shape forest structure and function by altering the frequency, intensity, duration, and timing of fire, drought, introduced species, insect and pathogen outbreaks, hurricanes, windstorms, ice storms, and landslides (Dale et al., 2001). Shifts in vegetation distribution in response to past and current climate changes have previously been described (Peñuelas and Boada, 2003; Vacchiano et al., 2014; Guillera-Arroita et al., 2015). If current climate trends continue or accelerate, major changes to forest management will become necessary (Hamann and Wang, 2006).

The vast subtropical regions of Yunnan extend to mid-elevation areas where the complex topography, which includes hills, basins, river valleys, stone forests and valleys of lime rock areas, is distributed at different elevations. Compared to the subtropical regions in eastern China, subtropical Yunnan has unique features such as lower heat, mild winters and two distinct seasons, dry and wet. Influenced by landform, climate and anthropogenic effects, the present vegetation maintains its own diversity and complexity (KIBCAS, 1994; Zhao et al., 2001). Yunnan's biodiversity is under considerable pressure and subtropical evergreen broad-leaved forests in this area have become increasingly fragmented due to agriculture, logging, planting of economic plants, mining activities and changing environment (Yang et al., 2004; Xu et al., 2005; Li et al., 2007, 2011; Zhou and Grumbine, 2011). With the combined impact of anthropogenic effects and rapid climate change, understanding the compositional patterns of Yunnan's subtropical evergreen broad-leaved forest species and identifying areas of high alpha diversity as well as priority areas for conservation has become more significant for the development of sound conservation policies and their integration into a sustainable land development strategy for Yunnan (Zhang et al., 2012).

Species distribution models (SDMs) have become a fundamental method in biogeography, ecology, biodiversity conservation and natural resources management (Guisan and Thuiller, 2005; Newbold, 2010; Franklin, 2013; Guisan et al., 2013; Guillera-Arroita et al., 2015). SDMs can be combined to model biodiversity at the community level following a ‘predict first, assemble later’ strategy, which has been used to identify threatened plant species hotspots, assess the invisibility of vulnerable native ecosystems and define areas of alpha diversity for conservation planning (Ferrier and Guisan, 2006; Parviainen et al., 2009; De la Estrella et al., 2012 Vorsino et al., 2014). In practice, however, SDMs provide only one point of view to correlate key environmental parameters with species distribution. Even though the models might underrepresent complicated natural ecosystems by neglecting species plasticity, adaptation, time-lag and biological interactions, they can be used as the primary technique for projecting vegetation range shifts, evaluating conservation priorities and assessing reserve designs (Hijmans and Graham, 2006; Gallagher et al., 2013; Duckett et al., 2013).

In this study, we examine how rapid climate change will influence distribution patterns of Yunnan's subtropical evergreen broad-leaved forests in order to develop conservation policies. Our aims are to (1) investigate climate change-induced changes of subtropical evergreen broad-leaved forests in Yunnan; (2) identify areas of current alpha diversity centers for conservation preparation.

1 Materials and methods 1.1 Study area and data

Yunnan province, SW China, is one of the most botanically diverse terrestrial regions on Earth. Located at a transitional zone, Yunnan possesses both tropical Indochina mixed and subtropical East Asian flora, while located between major floristic regions, with the Sino-Japanese region in the east and the Sino-Himalayan in the west (Li and Li, 1997; Myers et al., 2000; Zhu et al., 2006). The region also shows a rich diversity of forest types and has a disproportionate amount of China's overall floristic diversity (51.6%), with over 18, 000 plant species (Wu, 1987; Yang et al., 2004). In Yunnan, Evergreen broad-leaved forests almost extend all over whole tropical and subtropical Yunnan (Wu et al., 1987). According to ecological characteristics of vegetation, Yunnan's evergreen broadleaved forests can be divided into 5 categories: Monsoon evergreen broad-leaved forest (ME), Semi-humid evergreen broad-leaved forest (SH), Mountainous humid evergreen broad-leaved forest (MH), Mountainous mossy evergreen broad-leaved forest (MM) and Summit mossy dwarf forest (SM).

Each category is composed of a characteristic set of species chosen from the Vegetation of Yunnan (Wu et al., 1987). Using the table of classification for Yunnan's vegetation, we identified four important families in evergreen broad-leaved forests: Fagaceae, Magnoliaceae, Theaceae and Lauraceae. Because SDMs with too few occurrences are considered less accurate, only valid species with at least 15 unique presences were adopted to prevent the generation of low-performance models (Stockwell and Peterson, 2002; Hernandez et al., 2006; Wisz et al., 2008). In this study, more than fifty species were dominant in Yunnan's evergreen broadleaved forest. However, only fifty-five forest woody species, including eighteen constructive species and thirty-seven companion species, were qualified because of their reliable status in vegetation classification and their dominant position in comparison with other non-preferred tree species which only had extremely small records ( < 15). Every chosen species was assigned to one of five categories above in order to analyze environment conditions more individually and model species distributions more precisely as these five categories are greatly distinct in community constitution and structure as well as having very dissimilar terrains and climate specificity. The presence data of 55 forest woody species (Table 1) were provided by some herbariums in China, including 11 institutions in Chinese Academy of Sciences system and 21 national universities, totally 13, 038 records (Fig. 1). All these records were also stored in Chinese Virtual Herbarium (http://www.cvh.org.cn/) or National Specimen Information Infrastructure (http://www.nsii.org.cn/) and are generally accurate, although a few identification errors occurred. We were finally able to produce 2405 geo-referenced collections which were effective and workable to continue operating. Presence data were scored at 5 arc min grid cells (ca. 10 × 10 km), because this spatial resolution is able to match the resolution for environmental data well since spatial errors in the geo-referenced records cannot be ignored with an overly high resolution.

To successfully model the distributions of 55 woody species for current environmental conditions, we used 19 bioclimatic variables (average for 1950-2000), as well as mean elevation, population density, 8 land cover or land use variables and 7 soil quality variables which are all at 5 arc min grid cells level. Current bioclimatic variables were used from WorldClim v 1.4 dataset (http://www.worldclim.org/) and land cover or land use variables, and soil quality variables were extracted from the Harmonized World Soil Database v 1.2 (http://www.fao.org/soils-portal/en/) (Hijmans et al., 2005; Fischer et al., 2008). To model distributions for future conditions only bioclimatic variables were used, as we considered soil conditions and land use beyond prediction. Therefore, 19 bioclimatic variables at 5 arc min resolution were collected from the WorldClim v 1.4 dataset for future climate conditions in 2070 (average for 2061-2080). We adopted ACCESS1-0 general circulation model under IPCC-CMPI5 RCP4.6. The ACCESS1-0 model provides the best performance simulating the climatology of atmospheric general circulations in East Asia and reproducing the historical inter-annual variability and the consistency during twenty-first century projections (Tian, 2013; Perez et al., 2014).

To avoid multi-collinearity of variables which can result in model over-fitting, highly correlated environmental predictors were removed by Pearson's pairwise correlation analyses in R v 3.2.2 (https://cran.r-project.org/) when the correlation coefficient > |0.70| (Graham, 2003; Pearson et al., 2007). Because 55 species could not share totally equal environmental conditions, we generally did the Pearson's pairwise correlation analyses based on the distributions of five classified vegetation categories. Five sets of variables were produced for both current climate conditions and future climate conditions (Supplement1).

1.2 Species distribution modeling and testing

BIOMOD2 v 3.1-64 was used to create SDMs (Thuiller et al., 2009). BIOMOD2 is a freeware, open source, package, which can efficiently generate ensemble forecasting of species distributions, and was implemented in R (Thuiller et al., 2013). This technique has been shown to greatly improve the accuracy of predictions over single-algorithm approaches and has also been widely applied in biogeography, invasion biology and conservation biology (Marmion et al., 2009). Six modeling techniques implemented in BIOMOD2 were employed in fitting and averaging the predictions: generalized linear model (GLM), generalized additive model (GBM), artificial neural networks (ANN), multivariate adaptive regression splines (MARS), Random Forest (RF) and Maximum entropy model (MaxEnt). These six models were chosen based on computation requirements and ability to evaluate response curves. The same set of responses, predictors and scenarios were used within modeling strategies. Modeling specifications were as follows: (1) for GLM, Bayesian Information Criteria (BIC) was applied for the stepwise selection procedure; (2) for GBM, distribution parameter was Bernoulli. For other modeling techniques, we used the default settings, as these were optimized for SDMs (Thuiller et al., 2009). All models were fit using a maximum of 100 iterations. Models were built using 80% of the occurrences and pseudo-absence data, and the remaining 20% were used for evaluating predictions. Three partitions were produced by randomly splitting datasets for every species. True skill statistics (TSS) and the area under the receiver operating characteristic curves (AUC) were used as measures of model performance. TSS takes both omission and commission errors into consideration and ranges from -1 to 1 (Allouche et al., 2006). Values ranging from 0.2 to 0.5 were considered poor, from 0.6 to 0.8 useful, and larger than 0.8 were good to excellent (Coetzee et al., 2009). As for AUC, values are poor when in the range 0.5-0.7, fair when in the range 0.7-0.9 and excellent when greater than 0.9 (Swets, 1988). An ensemble forecasting framework was applied by averaging all the projections from every single model using TSS and AUC (Araújo and New, 2007). To build the prediction, only the models with AUC values above 0.75 were kept and averaged (Thuiller et al., 2013).

1.3 Spatial pattern of alpha and beta diversity

Stacked species distribution models (S-SDMs) allow us to analyze the potential influence of rapid climate change on biological communities (Guisan and Thuiller, 2005; Guisan and Rahbek, 2011). The process by which SDMs are combined into community-level models is often referred to as ‘stacking’ (Ferrier and Guisan, 2006). We projected S-SDMs over the study area in the current period as well as one future climate scenario. Because one main aim of our research is to estimate alpha diversity and beta diversity, it should be noted whether binary conversion of S-SDMs predictions need be applied (Guillera-Arroita et al., 2015). Recent studies that argue against discretization of S-SDM outputs for specific applications have shown that continuous outputs of SSDMs provide richer information than discrete outputs and therefore discretization of S-SDM results is harmful in the majority of applications (Calabrese et al., 2014; Lawson et al., 2014). Since thresholds will lead, quite generally, to biased results, in our study thresholds were not used to generate species richness, hereafter referred to as alpha diversity. The process of stacking SDMs represents a summation of the occurrence probabilities of every grid cell (Calabrese et al., 2014). The graphs of alpha diversity for five vegetation categories as well as total vegetation in both current and future climate scenarios were mapped. As for beta diversity, a 5 arc min spatial resolution's moving window performed using a ninepixel window size was operated to calculate beta diversity (Ochoa-Ochoa et al., 2012). Beta diversity was formed based on Whittaker's classic formula, β=γ/α, where γ (gamma diversity) means the total number of species within one moving window andaimplies the arithmetic mean of alpha diversity of all the pixels within the window. Whittaker's formula is not only straightforward to interpret, but also not sensitive to contrasting values of gamma diversity (Srivastava, 1999). We also obtained potential change of woody species diversity in study areas by overlapping their current and future distribution patterns.

2 Results 2.1 Potential distribution of 55 woody species

We present both the logistic probability of presence and the binary (e.g. suitable/unsuitable) output of the ensemble models for the 55 species over Yunnan based on ten sets of predictors: five sets for current predictions and an additional five for future predictions (Supplement 2). All SDMs were performed well, with individual models producing AUC mean values between 0.74 and 0.89 and TSS mean values between 0.47 and 0.74 within six models (Supplement 3). Only individual models having an AUC score higher than 0.75 were chosen to build ensemble models.

2.2 Potential distribution of alpha diversity, beta diversity and change rate

For this study, we obtained the geographic distribution of alpha diversity for Yunnan's evergreen broad-leaved forests by using the stacked models of all species. In order to see the predicted results directly, the maps were transformed to binary outputs, which were set with a series of thresholds (Fig. 2). We identified regions where alpha diversity remained at high levels across the maps based on a series of thresholds. Also, we generated projected spatial patterns of Yunnan's evergreen broad-leaved forests between the present and 2070 based on the stacked models, including current alpha diversity (Fig. 3a), alpha diversity in 2070 (Fig. 3b), change rate of alpha diversity (Fig. 3c), beta diversity (Fig. 3d), beta diversity in 2070 (Fig. 3e) and change rate of beta diversity (Fig. 3f).

After showing potential geographic transformation, we found great variation on alpha diversity and beta diversity of the evergreen broad-leaved forests, which were projected to experience decrease, increase or standstill by 2070 relative to the current forests. For alpha diversity, 35% of the current distribution of evergreen broad-leaved forests may increase in 2070, 51% may stay stable and 14% may decline. However, beta diversity may not enlarge notably under future climate scenarios. We found that 36% of areas would decrease, 58% remain stable and only 6% will increase (Fig. 4).

3 Discussion 3.1 The centers of alpha diversity

Alpha diversity is a basic measurement of regional and community diversity and forms the foundation of some ecological models of community structure (MacArthur and Wilson, 1967; Connell, 1978; Magurran, 2013; Stevens, 1989). It is practical to use alpha diversity for basic comparisons among sites as well as analysis on the saturation of local communities colonized from regional source pools (Cornell, 1999). In this study, we found that the major areas of alpha diversity were evenly distributed between western Yunnan, southeastern Yunnan and central Yunnan (Figs. 2 and 3). These alpha diversity centers showed a marked mountainous character, falling almost completely within main mountain ranges and highlands (López-Pujol et al., 2011). In western Yunnan, the major centers were situated in great valleys, hills and basins. These centers seemed comparatively fragmented, which corresponded roughly to crisscrossing mountains and river valleys. Situated from west to east, Gaoligong Mountain, Nu Mountain, and Yun mountains are accompanied by the Dulong, Nujiang, Lanang, and Jinsha rivers. In central Yunnan, alpha diversity centers are colonized with composite units of karst landforms, low mountains, hills, basins and lakes in the broad razing level. Richness centers of southwest Yunnan were situated in karst plateau and karst valleys where mountains, hills and basins intricately overlap. Climate dynamics in these centers are strongly and easily influenced by complex topography. For example, elevation range contributes greatly to micro-climate diversity in these richness centers, where massifs with more than four climate zones at similar latitudes are common. (Yang et al., 2004). Due to unique geographic development, climate variation, and relative isolation, each richness center has been able to maintain moderately high alpha diversity.

3.2 Change of alpha diversity and beta diversity

Future changes of forest distribution are closely related to the rapid changes of both temperature and precipitation (Li et al., 2015). Based on our research, results from different diversity parameters for estimating the influence of rapid climate change can be greatly dissimilar. Alpha diversity, which could possibly increase under future climate change for RCP 4.5 in this region (Figs. 3 and 4), means that the future climate could give species that belong to evergreen broad-leaved forests more suitable spaces. Our models predicted that the regions roughly corresponding to Yulong Mountain, Yaoshan Mountain, and Wumengshan Mountain on Yunnan's northern frontier will be areas where alpha diversity increases greatly in the future. The warmer climate may make forest species invade northern communities where deciduous forests are dominant. However, few studies have identified the northeast part of Yunnan as a conservation priority (Zhang et al., 2012) and more attention aimed at increasing alpha diversity should be focused on this area. Meanwhile, the change of alpha diversity in other domains should also be noted. For example, alpha diversity may slightly increase in western Yunnan, including the longitudinal valleys of Hengduan Mountains and the transitional zone between the Tibetan Plateau and Yunnan-Guizhou Plateau. As the Hengduan Mountain range stretches from north to south, cold-tolerant plant species from the central Tibetan plateau may be able to expand their distribution south, forming a cold-tolerant plateau biogeocoenosis. This is the lowest latitude where representatives of the Palaearctic and Holarctic realms may be found in the Northern Hemisphere (Chen, 1998). The increase in alpha diversity in western Yunnan could be a complex issue, due to special geographical features and effects on susceptible micro-climate conditions. Even though alpha diversity might increase, the major zones where evergreen broad-leaved forests exist are stable, such as Honghe, Wenshan and Kunming. In general, the distribution of alpha diversity would not move drastically and only the transformation in the marginal areas of the evergreen broad-leaved forest would seem obvious. On the contrary, the reduction in beta diversity in most of the current forests represents a general pattern of homogenization in vegetation communities. Two opposite trends for alpha diversity and beta diversity raise a question about how rapid climate change plays a role in forest structure. Global warming will elevate average temperatures and this effect may lead to transformations of ecosystems on earth. Particularly, the change could be drastic when it occurs at biodiversity hotspots. Yunnan is a wellknown biodiversity center, not only because of the extremely rich abundance in species diversity, but also because of a high ratio of endemic, rare and endangered species.

3.3 Five vegetation categories in the evergreen broad-leaved forest

Evergreen broad-leaved forests of Yunnan are one kind of complex community which contains five subtypes. In general, all subtypes could contribute to the overall alpha and beta diversity of the forest ecosystem and are obviously not equally weighted. For example, monsoon evergreen broad-leaved forests are mainly located in southeastern Yunnan where a large number of palaeoendemic species exist (López-Pujol et al., 2011). Many palaeoendemic species which show disjunctive distributions are the key species to the evergreen broad-leaved forest (Kruckeberg, 2004). These species were thought to be unable to recover their original range as a result of the constriction of their specialized habitats and the loss of genetic variability (Stebbins and Major, 1965; Kruckeberg and Rabinowitz, 1985). Generally, it is true that alpha diversity in southeastern Yunnan shows decline in the future. Meanwhile, beta diversity also shows falling in the coming period. As for semihumid evergreen broad-leaved forest, it occupies central Yunnan mostly and may extend its range under future climate conditions. Climate warming could provide semi-humid evergreen broadleaved forest with more suitable conditions for maintenance and expansion due to its preference for the drier climate in central Yunnan, compared with monsoon evergreen broad-leaved forest in the south of Yunnan.

Consequently, alpha diversity may increase in central Yunnan, as a result of new habitats produced by warming. However, beta diversity in central Yunnan will decline significantly, representing a general pattern of homogenization in semi-humid evergreen broad-leaved forest. Regarding the remaining three mountainous vegetation categories, mountainous humid evergreen broad-leaved forest, mountainous mossy evergreen broad-leaved forest and summit mossy dwarf forest, they are distributed largely on the top of the subtropical mid-mountains and high cloudy mountains greatly correlated with mountain distribution patterns in Yunnan. For instance, these mountainous forests chiefly in west Yunnan could widen their ranges. Increasing precipitation would represent a more and less beneficial alteration for mountainous forests because of their preference for humid conditions. To sum up, rapid climate warming may gradually cause changes in the water-energy dynamics that have an impact on the species diversity and habitat diversity of the mountainous forests.

When we focused on the diversity of forest communities, surprisingly, alpha diversity and beta diversity represented distinct distribution patterns respectively. This finding supports future conservation of evergreen broad-leaved forest in Yunnan, highlighting the need to consider fully the problem of vegetation homogenization caused by alterations in water-energy dynamics. Rapid climate change may gradually have severe effects on evergreen broad-leaved forest. However, we should note that rapid change may only partially contribute to the true vegetation distribution. In practice, forest communities largely depend on biological interactions and processes, such as inter-specific competition, predation on different trophic levels and migration activity. These changes may also have a profound effect on the structure and function of forest communities. For instance, the homogenization of forest communities may be generated from the interplay between environmental and biological interaction. If this occurs, it would substantially alter not only the environment but human society. If we carefully take habitat and species diversity into consideration, the work on forest protection, artificial afforestation and related programs will have a beneficial and sustainable development in future.

Funding

Biodiversity Conservation Research Project of Yunnan Environmental Protection Department (Y430112261).

Acknowledgments

The authors thank Ming-Gang Zhang and Katharina Filz for suggestions about problem of multicollinearity and thank Damien Georges for suggestions about modeling.

Appendix A. Supplementary data

Supplementary data related to this article can be found at http://dx.doi.org/10.1016/j.pld.2016.04.003.