Richness patterns of endemic and threatened conifers in south-west China: topographic-soil fertility explanation

Understanding the relationships between species richness patterns and environment constitutes a key issue in biogeography and conservation strategies. To our knowledge, this is the first integrative study that incorporates soil and human-influence data into species richness modelling. Our aims were to (a) estimate the richness patterns of four conifers groups (all conifers species, endemics, threatened, and endemic-threatened species) in south-west China, (b) assess the relative importance of environmental predictors (energy, water, climate, topography, and soil) and the human-influence on the conifers richness patterns and (c) identify hotspot ecoregions, nature reserves, or important plant areas as priority conservation areas. Generalized linear models and hierarchical partitioning were used by correlating 8962 distributional records of 97 conifer species with different environmental drivers. Results indicated that central Sichuan, northern Sichuan, northern Yunnan, and the southern areas of the Hengduan mountains were identified as distinct centres of conifers richness in China. Topographic heterogeneity and soil fertility were the strongest drivers of conifer richness patterns, while climate, energy, water, and human drivers were contributed to a lower degree. The identified conifers’ important areas were mostly located outside of the existing nature reserves but inside the ecoregions. Our findings emphasize that incorporating soil data into spatial modelling provides great insights for the conservation of conifers species. We recommend conservationists to use soil variables and other environmental data to generate a comprehensive understanding of the key drivers underlying the patterns of conifer diversity and distribution.


Introduction
Species richness can positively affect community stability and ecosystem functions such as nutrient cycling (Maestre et al 2012), resilience (Oliver et al 2015), productivity, and carbon storage (Liang et al 2016, Chen et al 2018, Liu et al 2018). Also, species richness of the rare, threatened, or endemic plants play a central role in identifying important plant areas (IPAs) that are used to determine the best sites for plant conservation (Anderson 2002).
Understanding the relationship between species richness and the environment is a key issue in ecology and biogeography. Although climate is one of the major factors driving species composition, habitat characteristics are also strongly associated with forest structure and composition (Zellweger et al 2015). In forestry, ensuring species diversity has become one of the important goals which require an understanding of the impacts of disturbance on species diversity and its management concerning other natural drivers (Schmiedinger et al 2012). Multiple drivers like climate, soil conditions, and human influence interact to affect the patterns of species diversity and composition (Ulrich et al 2014). Several hypotheses, such as climatic seasonality, environmental energy, water, and habitat heterogeneity (Currie et al 2004, Shrestha et al 2018, have been suggested to explain the diversity gradient. For example, the water-energy hypothesis suggests that an area with high energy and water availability will promote high species diversity (O'Brien 1998, Francis andCurrie 2003). Human disturbance is known to cause habitat fragmentation. It may have a stronger impact on the loss of species than global warming (Schmiedinger et al 2012, Venter et al 2016. The habitat heterogeneity hypothesis attempts to explain species richness through space niches and diversification (Moeslund et al 2013, Stein et al 2014. In addition, on the basis of geological and topographical characteristics for pertaining biodiversity, geodiversity has an effect on biodiversity at landscape and subnational scales (Hjort et al 2012, Gray 2013, Bailey et al 2017. Gray (2013) defined 'geodiversity' as the variability of geological, soil, and geomorphological characteristics and the physical processes that lead to these characteristics. Furthermore, terminology or definitions for describing the geodiversity has not been standardized (Anderson et al 2015). Geodiversity in this article is defined as the variability of geo-traits of morphological, physical, and chemical components representing three hypothesises of topographic heterogeneity, soil texture type, and soil fertility (Keith 2011, Räsänen et al 2016, Bailey et al 2017, Lausch et al 2019. The approach to compile geodiversity information was highly simplified in the current study, but the previous studies e.g. Hjort and Luoto (2012) and Hjort et al (2012) have described geodiversity as the variability in the Earth's surface materials, geomorphological and hydrological variability. Several studies recommend incorporating geodiversity quantitative data such as topographic heterogeneity or soil nutrients into spatial modelling and conservation studies (Mod et al 2016, Bailey et al 2017, Tukiainen et al 2017.
Several studies have neglected some meaningful variables directly or indirectly affect species distributions and ecosystem functions, such as available soil water, nutrients, carbon, potential evapotranspiration (PET), and solar radiation (SR) Gégout 2008, Mod et al 2016). Studies have either reported negative or positive relationships between soil nutrients and species richness (e.g. Huston 1980, van der Sande et al 2018, Yuan et al 2020), which has resulted in the drivers of plant species richness being controversial. As a result, recent studies recommend the integration of biodiversity, geodiversity, and anthropogenic factors into conservation planning (Tukiainen et al 2017, Xu et al 2017, Tripathi et al 2019. Conifer forests constitute terrestrial biomes present in tropical, subtropical, boreal, and temperate climates worldwide (Olson et al 2001). Globally, more than 50% of the temperate conifer forests are located in Asia and primarily in China (Farjon and Filer 2013, Wan et al 2017, Dakhil et al 2019. Nearly 50% of conifer species are threatened due to anthropogenic activities Zhao 2004, Farjon et al 2019). Woody species have been predicted to decline due to climate change , Dakhil et al 2019. Also, the nature reserve system of China has primarily focused on protecting mammals, while plant diversity and its habitats have not been well captured (Xu et al 2017).
We hypothesized that important conifer areas might be located outside the nature reserve system, which would indicate that they require protection. Therefore, our objectives were to (a) estimate the richness patterns of four conifer groups (all conifer species, endemics, threatened, and endemicthreatened species), (b) determine the relative importance of environmental predictors (environmental energy, climatic seasonality, water availability, and soil factors) and the human-influence on the species richness patterns, (c) determine the set of predictors that best explains the richness patterns, and (d) identify the boundaries of the identified IPAs.

Species distribution data and richness estimation
We downloaded occurrence data from the Global Biodiversity Information Facility (www.gbif.org), Chinese Virtual Herbarium (www.cvh.org.cn), National Specimen Information Infrastructure (www.nsii.org.cn/), and appropriate literature (Ying et al 2004). We cleaned the records by removing synonyms and unresolved names using the conifer database (Farjon et al 2019). Occurrence records from outside the study area and the duplicates were removed using ArcGIS v. 10.3 (ESRI 2014). Finally, we compiled 8962 occurrence records of 97 taxa, including 59 endemics, 63 threatened species, and 41 endemic-threatened species (see appendix A, table S1 (available online at stacks.iop.org/ERL/16/034017/mmedia)). As a conservative approach, we divided the species into four categories (all 97 conifers, endemics, threatened, and endemic-threatened species), and then estimated the richness based on the total number of species per grid cell (50 × 50 km 2 ) using SDM toolbox in ArcGIS v. 10.3 (Brown 2014). Several similar studies have used  (R 2 , %) of the predictors for the species richness patterns of all conifer species, endemic, threatened, and endemic-threatened species in south-west China evaluated by the negative binomial generalized linear model (GLM-NB). Stronger predictors are in bold and non-significant values are marked (ns). All other values are significant at p < 0.05. The positive (+) or negative (−) symbols in brackets indicate the relationship between species richness and the predictor variable.

Environmental and human-influence drivers
We used 23 predictor variables of seven hypotheses (table 1)  We downloaded the Human Influence Index (HII) (https://sedac.ciesin.columbia.edu/data/set/ wildareas-v2-human-influence-index-geographic) at a resolution of 30 arc-seconds, which represents human pressures (Wildlife Conservation Society-WCS 2005). The HII was created from nine layers of global data covering human population density, land use, built-up areas, night-time lights, land use/land cover, roads, coastlines, navigable rivers, and railroads (Xu et al 2014). Eight soil variables (see appendix A, table S4) at depth (0-2 m) with 1 km resolution were obtained from the ISRIC-World database, in which a high density of the sampled points was located in south-west China (ftp://ftp.soilgrids.org/data/aggregated; Hengl et al 2014). Soil nitrogen and phosphorus data were downloaded from http://globalchange.bnu.edu.cn/ research/soil2d.jsp (Shangguan et al 2013). The means of the raster-layer depths (0-2 m) were generated using the raster calculator in ArcGIS v.

Statistical and model analyses
Descriptive statistics (variance, mean, skewness, and kurtosis) and Spearman correlation analyses were applied among the four richness groups using SPSS v. 21.0 (IBM 2012). Given that the species richness values of all groups were abnormally distributed and often over-dispersed with variance-mean ratios greater than 1 (Ver Hoef and Boveng 2007; see appendix A, table S2), we used generalized linear models with negative binomial residuals (GLM-NB) to analyse the relationships between conifer richness and environmental-human variables (figure 1). The GLM-NBs could better compensate for overdispersed data compared to that of the Poissonian residuals, and this method has been widely used for ecological count data (Ver Hoef and Boveng 2007).) For each predictor, we performed a GLM to evaluate the explanatory power (R 2 , %) of each predictor individually, which was calculated according to (Panda et al 2017, Shrestha et al 2018 as follows: The glm.nb() function was used for GLMs in the 'MASS' R package (Venables and Ripley 2013) followed by stepwise regression to explore the combined effect of predictor variables on richness pattern and determine the best group of predictors with the highest explanatory power (Shrestha et al 2018). Due to high correlation, one variable from each predictor category was selected to build a model, and finally, seven predictors were selected (Faraway 2016) (see appendix A, table S4). We used all possible combinations of predictors using four energy, three water, two seasonality, two topographic heterogeneity, one human influence, four soil-texture types, and seven soil fertility variables, which yielded 1344 models (4 × 3 × 2 × 2 × 1 × 4 × 7) for each richness category. We selected the best model based on the lowest Akaike information criterion (AIC) for each richness group. Within each selected model, we calculated the variance inflation factor (VIF) for all predictors using the 'car' package in R v. 3.6.1 to evaluate the significance of multicollinearity (Legendre and Legendre 2012, Fox et al 2016), where collinearity is considered to be significant when VIF > 5. Next, we compared the explanatory power of each predictor individually using GLM, and their combined effect using stepwise regression. Similarly, we used a hierarchical partitioning model to compare the proportion of variance explained by each predictor category (Liu et al 2019). We first conducted a principal component analysis (PCA) with the varimax rotation in SPSS v. 21.0 (IBM 2012) within each of the predictor categories. We then extracted the first axis to represent each category (hypothesis). We standardized the human influence variable because the PCA cannot be applied to a group containing only one variable (Liu et al 2019). Then, we used the PCAextracted first axes to apply hierarchical partitioning analysis using the 'hier.part' R package (Walsh and Mac Nally 2013).

Conifer species richness patterns and hotspot ecoregions
Four distinct centres of conifer richness were identified based on the output maps produced (figure 2(a)). These centres were central Sichuan, northern Sichuan, northern Yunnan, and the southern areas of Hengduan Mountains. The conifer richness ranged from 1 to 30 (figure 2(a)) for all species. Two endemic conifer centres were identified central Sichuan and the northern Qionglai-Minshan ecoregion ( figure 2(b)), in which richness ranged from 1 to 24. Regarding threatened species richness (1-17 species per grid), two centres in central Sichuan and northern Yunnan were identified (figure 2(c)). Furthermore, one richness centre (central Sichuan) for endemic-threatened species was detected with a range from 1 to 14 ( figure 2(d)). The species richness of all categories was highly right-skewed (see appendix A, table S2). Four canters of species richness (proposed nature reserves) were located outside the national nature reserve system (figure 2(e)). Out of 34 non-threatened species, 18 species are endemic. Consequently, the total number of threatened or endemic conifer species is 81 (i.e. 18 (endemic nonthreatened) + 63 (all threatened)). Thus, 81 species endemic or threatened and represent 84% of the total studied conifer species. This percentage (84%) is high enough to represent the map of total species to represent the priority conservation areas for either endemic or threatened conifer species.
The results of the correlation analysis showed moderate to high concordance between conifer richness and the three remaining richness groups (r = 0.68-0.90, appendix A, table S3), implying that the factors driving conifer richness are possibly the same among the different species richness groups.

Relationship between different predictors/hypothesis and richness patterns
The GLMs showed that the explanatory variables of environmental energy, particularly MAT and PET, were stronger predictors than climatic seasonality and water availability, and these showed negative relationships (table 1). The MAT explained 6%-19%, while PET explained 4%-19% of the total variance in species richness. The proportion of variance also indicated that environmental energy determines conifer endemic-threatened richness (figure 3). Environmental energy, climatic seasonality, and water availability showed much lower contributions to the variation in species richness compared to that of topographic and soil predictors ( figure 3).
Of the 23 variables, the explanatory variables of topographic heterogeneity were the strongest predictors of species richness and showed positive relationships with all richness groups (table 1). Most of the soil fertility variables showed positive and stronger relationships than those of soil texture variables (table 1). The contribution of the humaninfluence variable was low for all groups and showed negative relationships among groups (table 1). The proportion of variance using the extracted principal components explained by the individual environmental and human-influence predictors highlighted the significant roles of topographic heterogeneity and soil fertility in determining species richness (figure 3).
Hierarchical partitioning modelling showed similar results to those of the GLM models, suggesting that topographic heterogeneity and soil fertility had the highest independent effects on all coniferrichness groups. In contrast, these variables had the highest joint effect on the richness of endemicthreatened species (see appendix A, figures S1a and S1b). The output of the stepwise GLM combined models showed that three models of the all conifer richness group, were rejected due to multicollinearity (i.e. VIF > 5), while for threatened richness group, only one model with the lowest AIC was rejected due to multicollinearity. For endemic or endemicthreatened richness groups, there was no model rejected because the predictor variables, of the model with the lowest AIC, has no multicollinearity. The combined models produced using stepwise regression (GLM) showed that ELER and soil nitrogen were consistently significant predictors (table 2). The

Conifer richness and conservation priorities
Our results showed that the major canters of conifers richness and endemism were located in Sichuan and Yunnan (Farjon and Filer  Similar studies have used the same spatial resolution as that of our study. They have found that the geographic area of central Sichuan is rich in~1000-1500 woody species per grid cell (50 km 2 ) with a total of 11 405 studied species, while the Hengduan Mountains in northern Yunnan have high woody species diversity (>1500 species per 50 km 2 ; Wang et al 2010). These hotspot areas could be assumed as priority conservation areas for woody species, particularly conifers (Wang et al 2010, Huang et al 2016.

The effects of energy, climatic seasonality, and water availability on conifer richness
The effects of water availability and climatic seasonality were relatively weak compared to those of environmental energy, and this due to the differences among conifer species with regard to the demand for water (Liu et al 2019). Generally, this finding supports the hypothesis that environmental energy limits regional species richness (Hawkins et al 2003, Pandey et al Table 2.
The best combinations of variables for each conifer group in south-west China evaluated using stepwise regression and their coefficients of determination (R 2 ). The best models for each group were selected from 1344 models based on the lowest Akaike information criterion (AIC). The numbers in parentheses are coefficients of respective variables from GLM. The variance inflation factors (VIF) for all predictors were less than 5, indicating insignificant multicollinearity.  The seedlings (early developmental stages of trees) are more susceptible to environmental constraints, with warming or heat stress that decreases the species richness of emerging seedlings (Robinson et al 2018). Also, the adverse effects of energy on conifer richness indicate that endemic-threatened conifers may face a high risk of extinction in the future due to global warming (Ye et al 2015). This indication supports the findings of conifer forest range dynamics under climate change conditions (Dakhil et al 2019). It confirms the suggestion of Rosenblad et al (2019) of species currently being restricted to particular environments, e.g. alpine areas being at risk of extinction due to climate change. It is important to note that the relationships between species richness and climate can be influenced by soil fertility in forest ecosystems . The different combined effects of these drivers support the suggestion that the relative effects of these drivers vary among common and threatened taxa (Shrestha et al 2018, Liu et al 2019.

The effect of topographic heterogeneity and human-influence factors on conifer richness
Our results support strong positive topographic heterogeneity-diversity relationships (see table 1  . Furthermore, historically high rates of plant diversification in south-west China have been due to topographical isolation based on a recent study of molecular phylogeny (Xing and Ree 2017). Hence, heterogeneity is thought to have a significant effect on the structure and dynamics of plant communities (Moeslund et al 2013, Lee andChun 2016). Moreover, topographic heterogeneity can be used as a proxy of habitat heterogeneity through the control of temperature, the distribution of soil water, and nutrient availability, particularly at high elevations, thereby controlling the occurrence and distribution patterns of species (Bisbing et al 2016, Tashi et al 2016). These results support our findings that temperature was negatively correlated with conifer richness. At the same time, topographic heterogeneity and soil fertility were greatly correlated with species richness (see tables 1 and 2). The negative correlation between human influence and conifer richness was comparatively weak and explained by 4%. This finding supports the results of Wang et al (2010), who studied the tree species richness in China. On the other hand, the human influence was comparatively strong for endemicthreatened species (Schmiedinger et al 2012, Venter et al 2016. Anthropogenic activity in south-west China is lower than that of Europe or North America. This may be because most of our study conifers were distributed at higher elevations compared to those of Europe and America (Schmiedinger et al 2012, Threatened Conifers of the World 2019).

The effects of soil texture type and soil fertility on conifer richness
Our findings support the assumption that sitespecific habitat factors (e.g. topography and soil) determine the patterns of species richness in addition to climate at the regional scale (Zellweger et al 2015). Soil texture has played a significant role in shaping plant richness patterns along the elevation gradient in eastern Himalaya (Sharma et al 2019), and this agrees with our findings (see figure 3). Furthermore, Sharma et al (2019) also found that soil type was the second strongest predictor explaining 36.4% of the variance in conifer richness. Soil texture and fertility also explained a notable portion of the variance, and similar findings were found in boreal conifer forests in Europe (Germany) and North America (Canada; Schmiedinger et al 2012).
Our results support the finding that explains the importance of sandy soils in the presence of high soil carbon, which leads to an increase in the waterholding capacity and aggregate stability with a greater number of macropores which are vital to soil functions such as cycling and storing of nutrients, based on their influence on water and air exchange, plant root exploration and habitat for soil biota (Soil Quality Institute 1999). The improvement of macropores and soil stability will increase soil infiltration capacity and enhance root penetration, subsequently boosting germination and growth of the seedlings of various species (Stein et al 2014, Quesada and Lloyd 2016, Yao et al 2019. Moreover, soil water variables contributed to species richness patterns more than climatic water predictors. This result verifies that soil water balance variables performed better than climatic water variables in tree species distribution modelling (Piedallu et al 2013).
Cation exchange capacity is an ecological indicator of soil fertility, explaining the variation in plant species richness in China (Zhang et al 2016), and this confirms our results. Soil fertility has direct and indirect effects on the growth of plant species, potential distributions, and thereby species composition (Mod et al 2016, Quesada and Lloyd 2016, van der Sande et al 2018. In addition, soil physical properties correlated with soil fertility and directly influence forest disturbance levels (Quesada et al 2012). This effect ultimately determines the turnover rates of tree species with differences in stand level floristic composition (Quesada et al 2012, Quesada andLloyd 2016). Also, soil fertility including the nutrient availability, increases the species diversity or the competitive advantage of certain species (Yuan et al 2020). Moreover, the finding that nitrogen and phosphorus positively contributed to the variation in species richness through increasing the soil biota, which increases decomposition. The diverse forms of such nutrients may promote plant diversity through resource partitioning when species differ in their preferences for different forms of these nutrients, leading to rapid forest regeneration (Laliberté et al 2014, van der Putten et al 2016). Soil organic carbon is one of the most important constituents of the soil due to its capacity to affect plant growth as both a source of energy and a trigger for nutrient availability through mineralization. Consequently, this provides aggregate stability and water holding capacity, and consequently, enhance seed germination and ecosystem multifunctionality (Soil Quality Institute 1999, Quesada and Lloyd 2016).

Geodiversity and Conservation implications
Geodiversity components (topography, soil texture and soil fertility) are the main strongest variables shaping the patterns of conifer richness. The soil fertility variables had strong positive effects on the richness patterns, which supports that, rarity may be more related to reduced soil fertility than climatic factors (Ulrich et al 2014, Quesada andLloyd 2016). Thus, areas with high richness and characteristic geodiversity features are key priority areas for conservation and provide greater insights into the distinctive soil conditions that support conifer diversity and ecosystem functioning and should be Plant habitats with high species richness hasn't been well captured by the national nature reserve system in comparison to mammalian habitats (Xu et al 2017). Accordingly, the IPAs with the high richness of endemic or threatened conifers (84% of the total conifer species) were found around giant panda and alpine forest nature reserves (see figure 2(e)) and require the establishment of new nature reserves for their protection, and this agrees with the findings of Xu et al (2017) and Dakhil et al (2019).

Limitations
Broad-scale environment-richness associations are stronger when the grain size increases (Xu et al 2014). On the other hand, some studies have found that landscape attributes may be significant at relatively small scales (e.g. 10-20 km grain sizes (Luoto et al 2007, Reino et al 2013). Moreover, the model contribution of geodiversity was quite similar in the case of the 25 km 2 and 50 km 2 resolutions, providing better model fits than that of the 100-250 km 2 resolution in native species richness (Bailey et al 2017). Hence, future studies should include different resolution scales (i.e. fine 5, or 10 km 2 ; coarse 25, or 50 km 2 ). It is important to note that we observed a suitable number of species in each 50 km 2 grid cell for all groups of conifer richness; hence, if this resolution becomes finer, the number of endemic threatened species in each grid cell will be fewer, and the variation between grid cells will be less so may will not enough to evaluate the relationship with the predictor variables.

Conclusions
Geodiversity, including topographic heterogeneity and soil fertility, was the main driver of conifer richness patterns. The Qionglai-Minshan, Hengduan Mountains, and Nujiang Langcang Gorge ecoregions were the major centres of threatened and endemic conifers in south-west China. Compared to all conifer species. Endemic-threatened conifers were more strongly affected by soil fertility, MAT, and topographic heterogeneity, and this indicates that endemic-threatened conifers may face a high risk of extinction in the future due to climate change and habitat destruction which induced by human activities that may reduce soil fertility. Also, this finding is a promising tool for afforestation planning of the threatened-endemic conifers, as soil fertility increase species diversity or the competitive advantage of certain species. The areas of central Sichuan, with a high richness of endemic-threatened conifers, could be considered to be priority areas for conservation because they were positively correlated with geodiversity, including soil fertility. Carbon, nitrogen, and topographic heterogeneity showed the highest contributions and strongest positive relationships with the richness of all conifer groups. This finding would help in afforestation planning. Climate, energy, and water availability were influenced by topographic-soil fertility. Accordingly, we recommend global researchers to incorporate soil fertility variables into the spatial modelling at community levels, and this will provide greater insights into conservation and management of the ecosystem multifunctionality.
The four proposed nature reserves or areas around the giant panda and alpine forest nature reserves in central Sichuan and northern Yunnan (see figure 2(e)) could be suggested as priority conservation areas. These areas should be protected to rescue all conifers of the endemic or threatened communities and their habitats.
The conifers in the current study may serve as a proxy for other plant or conifer taxa living in similar topo-geographical areas in Europe, America, or Canada. However, human influence should be taken into consideration because this is high in Europe. While the Chinese government has banned logging in south-west China, human disturbances in the region could be the result of urbanization, tourism activities, and agricultural expansions, which promote the growth of invasive species. As such, these human influences should be taken into consideration, along with the recent data of these disturbances.
Finally, logical directions of future research include a spatial conservation assessment of endemicthreatened conifer species under various scenarios of climate change, taking into account the findings of the current study when planning future conservation and ecosystem management strategies.

Data availability statement
The data that support the findings of this study are openly available.
All data that support the findings of this study are included within the article (and any supplementary files).