1. Introduction

Cushion-forming plants are considered ecosystem engineers in alpine and arctic ecosystems that create suitable micro-habitats for less stress-tolerant plants to inhabit (Arroyo et al., 2003; Körner, 2003; Badano et al., 2006; Cavieres et al., 2005, 2007; Chen et al., 2015b, 2019). As a result, cushion plants facilitate higher plant diversity and above-ground productivity at both local and regional scale (Arroyo et al., 2003; Butterfield et al., 2013; Cavieres et al. 2014, 2016; Chen et al., 2015a; Kikvidze et al., 2015; Gavini et al., 2020). Cushion plants also contribute to alpine arthropod diversity and community dynamics (Molina-Montenegro et al., 2006; Molenda et al., 2012; Reid and Lortie, 2012; Chen et al., 2021) by constructing and sustaining the structures of plant-pollinator networks in alpine ecosystems (e.g. Losapio et al., 2019), which is critical for alpine ecosystem function and sustainability (Badano et al., 2006). Because cushion plants are keystone species in alpine ecosystems, understanding their current patterns of diversity and geographical distribution allows us to predict future community dynamics of the entire alpine ecosystem and develop effective conservation strategies (Anthelme et al., 2014; Dolezal et al., 2019).

Cushion plants are distributed on mountains across the Earth. However, few studies have aimed to assess their patterns of diversity and geographical distribution (but see Aubert et al., 2014; Boucher et al., 2016). Aubert et al. (2014) revised the worldwide catalogue of cushion plants and confirmed that temperate Asia has 37.1% (487 species) of all species, the highest number of cushion plants in the world. Boucher et al. (2016) used this revised catalogue of cushion plant to confirm that the Himalayas are a center for cushion plant diversity and found that coldness was more important than dryness in explaining the presence of cushion plants. Furthermore, adaptive capacities of cushion plant have been shown to be species-specific (Badano and Cavieres, 2006; Chen et al., 2015a, b). Cushion plant typology differs by compactness, from loose dwarf/mat to densely compact cushions (Aubert et al., 2014). These variations suggest that cushion plant morphology may facilitate the ability of plants to adapt to different environmental stressors. However, previous studies have yet to identify the factors that drive such morphological and adaptive radiations.

The Qinghai-Tibet Plateau (QTP) is a cushion plant diversity hotspot that hosts over 110 cushion plant species (Aubert et al., 2014; Boucher et al., 2016; Li et al., 1987; Huang and Wang, 1991). The unique topology of the QTP, with closely clustered, huge mountains separated by deep valleys, has generated diverse climatic conditions across the entire QTP region, with relatively warm and wet conditions in the southern and southeastern parts and cold and dry conditions in the northern and northwestern parts (The comprehensive scientific expedition to the Qinghai-Xizang plateau, 1983; Chen et al., 2017; also see results in this study). These diverse habitats have in turn produced high plant diversity (more than 12, 000 seed plants, Zhang et al., 2016). In addition, the QTP has undergone massive climatic changes due to several uplift events, which have likely affected the distribution of plant species and the evolution of plant community structure (Wen et al., 2014 and references therein; Spicer et al., 2020). Thus, the patterns of cushion plant diversity and their distribution might differ across the entire QTP, especially because cushion plants are strongly associated with cold and dry conditions (Körner, 2003; Aubert et al., 2014; Boucher et al., 2016). Although previous studies have attempted to characterize the flora and distribution of cushion plants on the QTP, most of these studies were conducted at local or regional scales (Li et al., 1985, 1987; Huang and Wang, 1991; Huang, 1994). In addition, the QTP is very sensitive to climate change (Liu and Chen, 2000) and the continually increasing temperature in alpine regions due to global warming has already significantly affected the species composition, distribution and community dynamics of this alpine ecosystem (Lenoir et al., 2008; Qin et al., 2009; Elsen and Tingley, 2015; Liang et al., 2018). Understanding which specific climatic factors affect cushion plant communities is essential for planning suitable conservation strategies, as cushion plants are especially sensitive to climate change (Neuner et al., 2000; Cranston et al., 2015; Chen et al., 2020).

Phylogenetic approaches are frequently used to measure biodiversity and characterize the assembly of evolutionary history in geographical space (Forest et al., 2007; Thornhill et al., 2016, 2017; Zhang et al., 2021b). Additionally, phylogenetic approaches can reveal how climatic features filter or accumulate (phylogenetic clustering or overdispersion) plant species under specific climatic conditions (Webb, 2000; Swenson et al., 2012; Aldana et al., 2017; Liu et al., 2019). Generally, strong climatic oscillations induce phylogenetic clustering, whereas milder climatic changes induce higher phylogenetic diversity (PD) and phylogenetic overdispersion (Feng et al., 2014). Although cushion plants are among the best adapted life-forms, especially in the alpine ecosystem of the QTP where they play vital roles in the local ecosystem functions (Li et al., 1985; Yang et al., 2010; Butterfield et al., 2013; Liu, 2014; Chen et al., 2015a, b, 2019; Cavieres et al., 2016), phylogenetic community structures of cushion species under varying environments on the entire QTP have yet to be characterized. Consequently, the evolutionary history of cushion plant assemblages and the factors that have driven this evolution across geographic space remain unclear.

In this study, we assess patterns of total cushion plant diversity and phylogenetic community structure to answer the following questions: (1) What are the current patterns of cushion plant diversity and their distribution on the QTP? (2) What potential climatic variables drove the current patterns of cushion plant diversity and distribution? (3) What factors drive differences in distribution of cushion groups with distinct typology?

2. Materials and methods 2.1. Data sources for cushion plant distribution

We used an updated worldwide catalogue of cushion plant species (Aubert et al., 2014; http://www.cushionplants.eu/) to identify and gather data on 123 cushion-forming plant species in the Qinghai-Tibet Plateau (QTP). The database was checked for recent additions, which were added to our species list. We also used the database to categorize plants according to their level of compactness. Of the total cushion plant species present in the QTP, 92 species have densely compact morphology (either dome-shaped or flat) and 31 species have "loose morphology" (either dome-shaped or mat) (Table S1). The densely compact cushion species used in this study are distributed in 22 genera and 12 families, most belonging to Saxifraga (Saxifragaceae), Arenaria (Caryophyllaceae) and Androsace (Primulaceae). The loose cushion species belong to 19 genera and15 families, most belonging to Eritrichium (Boraginaceae), Oxytropis (Fabaceae) and Acantholimon (Plumbaginaceae).

County-level distribution records of cushion plant were transferred into a 0.5-degree × 0.5-degree resolution following the protocol adopted by Zhang et al. (2021a). Distribution records were collected from national and regional flora records, Chinese Virtual Herbarium (http://www.cvh.ac.cn), Zhang et al. (2016), and the National Specimen Information Infrastructure (NSII: http://www.nsii.org.cn).

2.2. Bioclimatic variables

Cushion plants are strongly associated with cold temperatures and dry conditions. Thus, in this study, we examined whether, temperature and precipitation play a role in driving cushion plant distribution. For each parameter, we selected four sub-parameters that directly reflect the relevant conditions. Mean annual temperature and annual precipitation are climatic variables frequently used to represent energy and water conditions (e.g., Xu et al., 2013). For temperature, we selected mean annual temperature (Bio1), temperature seasonality (Bio4), mean temperature of warmest quarter (Bio10) and mean temperature of coldest quarter (Bio11). Temperature seasonality is defined as the standard deviation of monthly mean temperature and precipitation seasonality is defined as the coefficient of variation of monthly precipitation both can reflect the climatic instability (Shrestha et al., 2017). For precipitation, we selected annual precipitation (Bio12), precipitation seasonality (Bio15), precipitation of wettest quarter (Bio16) and precipitation of driest quarter (Bio17). Importantly, mean temperatures of warmest or coldest quarters, mean precipitation of wettest or driest quarters can reflect extreme or limiting environmental conditions, which are assumed to affect plant distribution and survival (Liu et al., 2020). For each bioclimatic variable of a species, we calculated the mean climatic values of all the distribution grid cells for this species. Bioclimatic variables were obtained from the WorldClim database v.1.4 (Hijmans et al., 2005).

2.3. Statistical analysis

For phylogenetic analyses, we downloaded cushion plant ITS sequences (ITS1-5.8S-ITS2) from the NCBI database (https://www.ncbi.nlm.nih.gov/). For species without sequences, we used ITS sequences of closely related species (GenBank accessions, Table S2). Sequence alignments and model selection were performed according to Zhang et al. (2019). A forced bifid branch phylogenetic tree without time calibration was constructed using BEAST v.1.8.4 (Drummond and Rambaut, 2007). We conducted an uncorrelated relaxed clock and birth-death process and an independent run of Markov chain Monte Carlo (MCMC), totaling 100 million generations, sampling every 5000 generations. TRACER v.1.6 (Drummond and Rambaut, 2007) was used to test the effective sample size (ESS) of the log file and ESS > 200 was regarded as having less influence by autocorrelation. The maximum clade consensus tree was generated in TREE ANNOTATOR v.1.8.4 (Drummond and Rambaut, 2007) and the first 25% of trees were discarded as burn-in. We used BIODIVERSE v.2.0 (Laffan et al., 2010) to calculate species richness (SR), phylogenetic diversity (PD), net related index (NRI) and the categorical analysis of neo- and paleo-endemism (CANAPE). NRI values reflect the phylogenetic structure in a grid cell (e.g., phylogenetic clustering or overdispersion) (Webb et al., 2002). CANAPE can identify these grid cells as evolutionary "museums" or "cradles", which shows important conservation values (e.g., Zhang et al., 2021b). Blomberg's K was calculated to estimate the phylogenetic signal using the package "phylosignal" (Keck et al., 2016) in R 3.5.1 (R Core Team, 2018). A K value close to 1 indicates that there is a certain degree of phylogenetic signal or of conservatism during the evolutionary process, a K > 1 indicates that traits are conserved, and K close to 0 indicates that evolution tends to be random.

To account for spatial autocorrelation, spatial simultaneous autoregressive error (SAR) models were run using MuMIn (Bartoń, 2019) and spdep (Bivand et al., 2018) packages. Moreover, we conducted ordinary least squares (OLS) linear regressions. OLS and SAR were only used to explore bivariate relationships between spatial indices and each variable due to the collinearity of the variables. In order to contain collinear variables, we conducted a partial regression analysis using the "vegan" package (Oksanen et al., 2015). All the bioclimatic variables were divided into four groups: energy, water, seasonality and extreme environmental factors. We conducted the same procedures on compact and loose cushion plants, respectively, to investigate the effects of climatic variables on different cushion structures (compact and loose).

3. Results 3.1. Cushion plant richness and distribution patterns

Cushion plant richness was negatively correlated with latitude (Fig. A1), indicating that geographical attributes affect the distribution of cushion plants. The highest total cushion species richness was found in the southern and southeastern QTP area, especially in the south-central Hengduan Mountains; in contrast, from south to north and southeast to northwest, total cushion plant richness declined (Fig. 1a). The pattern of cushion plant phylogenetic diversity (PD) was similar to that of total cushion plant richness; however, PD showed a different pattern than that of total cushion plant richness (Fig. 1b). Specifically, cushion plants showed higher PD in the southern part of the QTP, with the south-central Hengduan Mountains harboring the highest PD, which progressively decreased to the northern and northwestern parts (Fig. 1b).

Consistent with total cushion plant richness, cushion plants with compact morphology showed higher diversity in the southern part of the QTP than in other areas (Fig. 2a). In contrast, the loose cushion plants showed a relatively random diversity pattern, with the south-central Hengduan Mountains and the central and northern QTP harboring the highest diversity (Fig. 2b).

3.2. Cushion plant phylogenetic community structure

Cushion plants from the eastern and southwestern parts of the QTP showed higher NRI (most > 0), indicating that the phylogenetic community structure is more clustered in these areas (Fig. 1c). In contrast, in other areas, especially in the central and northwestern parts of the QTP, the NRI values were relative lower (most < 0), indicating that the phylogenetic community structure in these regions is frequently random or overdispersed (Fig. 1c). Generally, higher degrees of phylogenetic clustering were observed in the eastern and southwestern cushion plant communities than communities in other areas of the whole QTP.

Additionally, CANAPE analysis identified only one area of paleo-endemism for cushion species: the northeastern QTP. This finding indicates this area might act as an evolutionary museum for cushion plants (Table S3). Neo-endemism centers were mainly found in the southern QTP, whereas very few were found in the western, northern and northwestern QTP, indicating the southern QTP might act as an evolutionary cradle for cushion plants. Mixed- and super-endemism centers were identified in the southeastern and northwestern QTP and scattered centers were distributed in the northeastern and southwestern QTP, indicating these areas might act as both the above roles (Fig. 1d; Table S3).

3.3. Drivers of cushion plant diversity and distribution

Overall, mean annual temperature, annual precipitation and precipitation of the wettest quarter progressively decrease from the southeastern to the northwestern parts of the QTP (Fig. 3a, e, g); whereas the temperature seasonality shows the opposite trend, with the northeastern QTP having the lowest seasonality and northwestern QTP having the highest seasonality (Fig. 3b). The highest mean temperature of the warmest quarter was found around the QTP, especially in the southeastern part, i.e., the Hengduan Mountains subregion (Fig. 3c). The highest precipitation seasonality was found in the south-central QTP (Fig. 3f), while the highest precipitation of driest quarter was found in the southeastern and southwestern parts of the QTP (Fig. 3h).

All the climatic variables represented by annual energy and water trends, seasonality and extreme environmental factors had significant effects on patterns of cushion plant diversity (Table 1). According to the SAR model, all climatic variables negatively affected cushion plant richness and phylogenetic diversity (PD), except for precipitation seasonality (Bio15), which positively affected diversity (Table 1). These results indicate that poor environmental conditions (e.g., dry and cold) promoted cushion plant diversity. However, most climatic variables had no effect on cushion plant phylogenetic community structure, with the exception of precipitation seasonality, which showed significantly negative effects, and temperature seasonality and mean temperature of warmest quarter, which showed marginally significant effects (Table 2). These results indicate that the intensity of climate fluctuations did not strongly filter the cushion plant community assemblages. Partial regressions indicate that seasonality and extreme environmental factors accounted for more variance than energy and water in explaining cushion plant richness and PD patterns (Fig. 4a and b), but these factors showed low explanatory power in NRI (Fig. 4c). Climatic variables produced no phylogenetic signal (all K ~ 0), indicating that the effect of climatic variables was randomly distributed across cushion plant phylogeny (Fig. 5; Fig. A2). Specifically, the distribution of all cushion plants can be predicted by coldness (Bio1), water condition (Bio12), seasonality (Bio15) and extreme climate (Bio10, Bio11 and Bio16), indicating that cushion life-form is effectively adapted to these climatic factors. However, the ability to adapt to most of these climatic factors is stronger in compact cushion plants than in loose cushion plants (Table S4).

4. Discussion

This study demonstrates that the patterns of cushion plant diversity and their distribution on the Qinghai-Tibet Plateau (QTP) are driven by environmental stressors, including temperature and water conditions. These findings increase our understanding of the evolution and formation of alpine cushion plant diversity in this region and establish an important foundation for future research and conservation efforts.

4.1. Cushion plant diversity and distribution

Total cushion plant richness is highest in the southern and southeastern parts of the QTP, especially in the south-central Hengduan Mountains subregion, and decreases from south to north and from southeast to northwest (Fig. 1a). This pattern indicates that the Hengduan Mountains subregion and its adjacent areas are important diversity centers of cushion plants. Such results are consistent with those from other plant taxa and functional groups (Wu, 1988; Wu et al., 1995; Yu et al., 2017; Zhang et al., 2017; Shrestha et al., 2017), including Rhododendron (Yu et al., 2017; Shrestha et al., 2017), Saussurea (Xu et al., 2019; Zhang et al., 2021a) and many other groups (Zhang et al., 2017, 2021b). Multiple factors have been proposed to explain the high plant diversity in these mountainous areas (e.g., Sun et al., 2002a, b; Wen et al., 2014; Zhang et al., 2016; Xing and Ree, 2017; Shrestha et al., 2017; Yu et al., 2019). One explanation may be that the rapid uplift of the QTP has led to increased topographical complexity and drastic disparities in climatic conditions on the plateau. For instance, the eastern part of the QTP, where the Hengduan Mountains are located, receives heavy rainfall in the summer, whereas precipitation decreases along the remainder of the QTP (An et al., 2001; Shrestha et al., 2017). These regional differences in precipitation are likely to have driven plant speciation and diversification (Sun et al., 2002a, b). Such mechanisms could, solely or jointly, explain the establishment of cushion plants in this region.

Our results also showed that the general trends of most climatic features are in accordance with the patterns of cushion plant diversity (Figs. 1a, b and 3). For instance, climatic seasonality (i.e., temperature and precipitation seasonality) and extreme climatic variables (i.e., temperature of warmest/coldest and precipitation of wettest/driest seasons) were highest in areas with the greatest cushion plant richness and diversity. The SAR results showed that all climatic variable examined, with the exception of temperature seasonality on total richness and annual precipitation and precipitation of wettest quarter on phylogenetic diversity (PD), had significant impacts on cushion plant total richness and PD (Table 1). These findings indicate that climatic conditions significantly affected the current patterns of cushion plant diversity, and that cushion plants are more tolerant to large climatic fluctuations and extreme climatic events (e.g., Anthelme et al., 2014) similar to those induced by the geological events that characterize the history of QTP.

Cushion plants with different typology (i.e., compact vs. loose) have different patterns of diversity on the QTP (Fig. 2). Specifically, compact cushion plants are concentrated in the southern part of the QTP, whereas loose cushion plants showed a relatively random distribution (Fig. 2). Although previous studies have indicated how cushion plants adapt to cold and drought conditions on the QTP (Yang et al., 2010, 2017; Liu, 2014; Chen et al., 2015b, 2019), these studies focused on compact cushion plants, failing to address which factors drive the distribution of different cushion plant typologies. Our results show that both compact and loose cushion plants are well adapted to most climatic variables examined here (Table S4). However, the adaptive capacity of compact cushion plants is stronger than that of loose cushion plants. Specifically, compact cushion plant diversity is more in accordance with the patterns of temperature seasonality (Bio4), annual precipitation (Bio12) and precipitation in wettest quarter (Bio16) (Fig. 2, Fig. 3); whereas, loose cushion plant diversity was not correlated with climatic factors (Fig. 2, Fig. 3). These findings may indicate that capacity to collect and store water from rainfall (and thus regulate their temperature) is greater in compact cushion plants than in loose cushion plants, which allows compact cushion plants to resist larger temperature fluctuations (i.e., seasonality). Thus, these climatic factors appear to jointly promote the diversification of compact cushion plants in the southern and southeastern QTP, whereas the microclimates induced by topographical features of the QTP may play a role in driving the current, random patterns of loose cushion plant diversity and their distributions. However, further evidence is needed to confirm these conclusions.

Previous studies have suggested that harsh environmental conditions lead to phylogenetic clustering in plants, whereas mild conditions lead to phylogenetic overdispersion (e.g., Kluge and Kessler, 2011; Li et al., 2014; Zu et al., 2019). Phylogenetic clustering also indicates that plants in harsh environmental conditions show niche conservatism (Webb et al., 2002), in which habitat filtering only allows some lineages to survive harsh environments (Li et al., 2014). In contrast, phylogenetic overdispersion indicates that negative interactions (e.g., competition) are important in plant assembly in such communities (Graves and Gotelli, 1993; Webb et al., 2002). Climatic variables showed no evidence of phylogenetic signal, indicating that the adaptation of cushion plants of different taxa to climate is relatively consistent across phylogeny. However, NRI showed relatively higher values (> 0) in the southern and eastern parts than in the remaining areas of the QTP (Fig. 1c), indicating that cushion plants still show spatially phylogenetic clustering in these areas. Spatial phylogenetic clustering indicates that environmental stressors may filter certain cushion plants. As a result, closely related species with the capacity to adapt to the same environmental stresses persist in these environments. This finding is consistent with previous studies that showed communities of alpine plants are phylogenetically clustered likely due to environmental filtering at high elevations (e.g., low temperature and precipitation) (Li et al., 2014; Zu et al., 2019). Phylogenetic clustering of cushion plant communities in the southern and eastern parts of the QTP may be explained by the extremely high topographical complexity of the region, which is characterized by huge mountains and deep valleys (Muellner-Riehl, 2019). Such topographical features can, in combination with changes in climatic features between various micro-topographies, promote species accumulation by offering more niche space, preventing extinction and providing increased opportunities for allopatric speciation (e.g., Kadereit et al., 2004; Linder et al., 2014; Wallis et al., 2016; Shrestha et al., 2017). In the central-northern and northwestern QTP, which have a less complex topography (i.e., plateau terrace, basins or semiarid deserts) and relatively uniform climate features (i.e., arid and/or semiarid) (The comprehensive scientific expedition to the Qinghai-Xizang plateau, 1983; Chen et al., 2017; Muellner-Riehl, 2019), cushion plants are distributed randomly or compete with each other, as indicated by random or overdispersed community phylogenetic structures (Fig. 1c). Relatively uniform climates might lead to high homogeneity in niche requirements of cushion plants, which in turn may lead to neutral or even highly competitive interspecific interactions.

4.2. Future studies on cushion plant functions and conservation on the QTP

Climate change has decreased biodiversity globally (Raven and Wackernagel, 2020 and references therein). Because mountain ecosystems, particularly those in the QTP, are located on the top of the earth's surface, they are extremely sensitive to climate changes (Körner, 2003; Liu and Chen, 2000; Lenoir et al., 2008; Dullinger et al., 2012). In fact, some cushion plant species are already threatened at high altitudes by surrounding vegetation (Chen et al., 2020). Cushion plants play a key role as foundation species in mountain ecosystems (e.g., Molenda et al., 2012; Cavieres et al., 2014, 2016; Kikvidze et al., 2015; Yang et al., 2010; Chen et al., 2015a, b, 2019; Gavini et al., 2020). The extinction of cushion plants (as key species) in the community could lead to local extinction or even cause, through cascading effects, the collapse of the whole species networks (Memmott et al., 2004; Fortuna and Bascompte, 2006; Losapio and Schöb, 2017). Therefore, preserving cushion plant communities, especially in centers of diversity and endemism (Fig. 1), is very important for the conservation of alpine diversity and the long-term sustainability of mountain ecosystems.

This study has primarily revealed the current cushion plant diversity and their distribution on the QTP. However, a number of unknows remain. Importantly, the evolutionary history of cushion plants on the QTP remains unclear. Understanding where and how cushion plants evolved may be a requisite for predicting where and how they will fare in the future. An additional major concern for conservationists is that most cushion plant centers of diversity are outside of established nature reserves (Fig. 1d; also see Zhang et al., 2021b). In addition, although cushion plants are well known to positively influence plant and arthropod diversity and community dynamics (Molenda et al., 2012; Liczner and Lortie, 2014; Kikvidze et al., 2015; Cavieres et al., 2016; Gavini et al., 2020; Chen et al., 2021), little is understood about how they contribute to species networks. Finally, understanding how specific cushion plant communities recruited, degraded and sustained under specific harsh environments and particular areas (like the paleo-endemism and neo-endemism centers) could be very important not only for predicting community dynamics but also for establishing suitable conservation programs.

Author contributions

JG Chen and H Sun planned the study. YZ Zhang, LS Qian and XF Chen collected and organized the data. YZ Zhang and LS Qian analysed the data with input from JG Chen and L Sun. JG Chen wrote the first version of the manuscript with input from YZ Zhang and LS Qian into the method section. All authors made efforts to improve the manuscript.

Declaration of competing interest

No conflict of interest exits with the submission of this manuscript, and it is approved by all authors for publication.

Acknowledgement

This study was supported by grants from the Second Tibetan Plateau Scientific Expedition and Research (STEP) program (2019QZKK0502), the Strategic Priority Research Program of Chinese Academy of Sciences (XDA20050203), the Key Projects of the Joint Fund of the National Natural Science Foundation of China (U1802232) and the Yunnan Applied Basic Research Project (202001AT070060). None of the contents of this manuscript have previously appeared online.

Appendix A. Supplementary data

Supplementary data to this article can be found online at https://doi.org/10.1016/j.pld.2021.09.001.