Warming-induced shrubline advance stalled by moisture limitation on the Tibetan Plateau

1 –––––––––––––––––––––––––––––––––––––––– © 2021 The Authors. Ecography published by John Wiley & Sons Ltd on behalf of Nordic Society Oikos This is an open access article under the terms of the Creative Commons Attribution License, which permits use, distribution and reproduction in any medium, provided the original work is properly cited. Subject Editor: Jian Zhang Editor-in-Chief: Miguel Araújo Accepted 30 July 2021 44: 1–11, 2021 doi: 10.1111/ecog.05845 44 1–11


Introduction
More than 50% of the Earth's terrestrial surface is covered by treeless vegetation communities (Friedl et al. 2010). Among these biomes, alpine regions provide habitats for endangered and endemic plant and animal species, support the livelihood of many people and are highly relevant for surface energy budgets, hydrological and biogeochemical cycles (Metcalfe et al. 2018). A comprehensive understanding of the spatiotemporal dynamics of alpine woody plants is therefore essential for understanding and sustainably managing wildlife biodiversity, trophic cascades, ecosystem services 2 and assessing vegetation/climate feedbacks (Pearson et al. 2013). Forecasting the range dynamics of alpine plants under a warming climate, however, is challenging, and current projections are highly uncertain due to the complex interplay of various determinants of range (IPCC 2014, Ford andHilleRisLambers 2020).
Rates of warming are more than twice as high at high altitudes than for the rest of the world (Pepin et al. 2015). Many montane and alpine plant species have accordingly already started to shift their seasonal activities, abundance and geographic ranges in response to climatic warming (Case and Duncan 2014, IPCC 2014, Vandvik et al. 2020, Bader et al. 2021. The expansion of shrublines to higher altitudes or more northern latitudes is one of the most striking examples of the impacts of global warming (Dial et al. 2016, Myers-Smith andHik 2018). Several studies of tundra have reported the northward expansion of shrublines into the Arctic (Myers-Smith and Hik 2018). The thermal limitations of woody plants at the limits of their latitudinal ranges are comparable to those at high altitudes (Randin et al. 2013), so latitudinal and alpine shrublines should share similarities in their demographic and range dynamics. Evidence for the advance of alpine shrublines in high-altitude regions, however, is limited (Lu et al. 2021).
Shrub encroachment under global change is occurring across large parts of alpine and arctic biomes (Hallinger et al. 2010, Criado et al. 2020. Such land-cover transitions, however, are spatially heterogeneous and temporally unstable (Brodie et al. 2019). For example, the rate of shrubline advance in southern Alaska closely matches warming trends (Dial et al. 2016), whereas shrublines in central Alaska are advancing more slowly (Brodie et al. 2019). The biological responses of shrublines to climatic warming vary with time, in part due to adaptation (Myers-Smith et al. 2020). These responses demonstrating disequilibrium between shrubline shifts and climatic warming at a broad spatial scale or over different time periods suggest that factors other than temperature (e.g., seasonal availability of moisture) could slow or even inhibit shrubline dynamics driven by warming (Olofsson et al. 2009, Frost and Epstein 2014, Grytnes et al. 2014. Thus, undisturbed shrubline dynamics would be increasingly codetermined by moisture at cold and dry sites, because the climate is expected to be dry under continued warming if precipitation does not increase proportionally . Alpine shrublines subjected to cold and dry climatic conditions are common on the Tibetan Plateau (TP), where both air temperature and atmospheric aridity have increased substantially in recent decades (Yao 2019, Yuan et al. 2019. The upper boundaries of alpine vegetation across the TP represent a natural laboratory for assessing the effects of global change on alpine plants (Lu et al. 2019a). Subnival ecosystems below permanent snowlines and above treelines play an increasingly important role in regulating hydrology, nutrient turnover, carbon cycling and biodiversity as snowlines and glaciers retreat (Anderson et al. 2020). Salix oritrepha Schneid. (hereafter Salix) is a widespread alpine deciduous shrub constituting the upper shrublines across the TP (Fang et al. 2011). Importantly, seed dispersal does not limit the upward expansion of Salix shrubs (personal field observation), so this species represents a model plant for monitoring range shifts in response to climate change. The structure and spatiotemporal patterns of TP Salix shrublines and their long-term responses to changes in temperature and precipitation, however, remain unknown.
We addressed this knowledge gap using an in-depth survey across 24 plots encompassing the Salix shrubline ecotone along a latitudinal gradient of 30-38°N on the eastern TP (96-103°E) (Fig. 1). The 24 Salix shrubline plots are little affected by local human disturbances (Supporting information) due to poor road access, long walking distances, the presence of topographic barriers (e.g., narrow valleys or river crossings) and a broad altitudinal range of shrublands or forests with dense plant communities. The dynamics of these undisturbed shrubline sites on the TP can thus be directly linked to climate variability. We reconstructed spatiotemporal shrubline dynamics using dendroecological techniques. This study 1) reconstructed the variations of Salix shrub recruitment and shrubline position in recent 80 years on the TP and 2) evaluated the relative influences of temperature, precipitation and moisture variables on Salix shrubline dynamics. The growth of Salix shrubs above the treeline is assumed to be mainly limited by low temperature (Lu et al. 2019b), so our hypothesis was that infilling and shrubline advances of Salix shrubs would be driven by climatic warming but that the rate of these processes would be modulated by warming-induced moisture stress.

Study area
A network of 24 shrubline plots was established on the eastern TP to reveal the large-scale changes of alpine willow (Salix) shrublines in western China (Fig. 1). The southern and northern parts of this region is influenced by the Indian and East Asian monsoons, respectively (Yao et al. 2012). The study region is characterized by cold and continental climate. The warmest and coldest months are July (mean temperature of ca 10°C) and January (<−10°C), respectively. The significant warming trend in summer and enhanced annual and summer evaporative demand since 1950 were detected in the study region ( Fig. 1 and Supporting information). Increased summer drought events were also found during the past decades in this region (Supporting information). In particular, a negative trend in summer Palmer drought severity index (PDSI) series over last decade was detected in a substantial number of sites (Supporting information).

Field sampling
Salix is an erect deciduous shrub belonging to the Salicaceae family, ranging from 3000 to 4900 m on the north-facing 4 slopes of the eastern TP (Fang et al. 2011). Field observations, however, indicated that isolated Salix shrub could reach an elevation of 4960 m in Basu County on the eastern TP, implying that seed dispersal did not limit the advance of Salix shrublines. Salix shrubs mainly grow on shady and semi-shady slopes in alpine regions, are the dominant woody species above the treeline and form the alpine shrubline in the high mountains on the TP. Crushed stones are commonly seen in the substrate with some alpine plants dominating the understory (Androsace tapete Maxim., Caragana jubata (Pall.) Poir., Rhodiola fastigiata (Hook. f. et Thoms.) S. H. Fu, Artemisia dalai-lamae Krasch.).
The Salix shrubline ecotone forms a narrow zone of transition varying from 40 to 135 m. A total of 24 rectangular plots 30 m wide and lengths varying between 50 and 135 m were sampled between July 2018 and August 2020 at a total of 10 subregions in the main alpine regions across the eastern TP (Fig. 1, Supporting information). Each plot was established in a topographically uniform area of the ecotone. The longer sides of the plots paralleled the steepest slopes and encompassed the current shrublines (Wang et al. 2015). A GPS device was used to measure the upper and lower elevations of each plot. The lower-left corner of each plot facing up the slope was set as the origin (0, 0) of a local coordinate system. Several variables were recorded for each Salix individual in the field: location (x and y coordinates), maximum height and basal diameter. The diameter of the thickest basal stem of each individual was measured using an electronic caliper (±0.01 mm). Shrub height (±10 cm) was measured using a measuring tape. None of the plots contained dead shrubs.
We collected cross sections about 3 cm thick at the base of the stems of all shrubs in each plot. Samples of wood were processed following standard dendrochronological procedures, including air-drying, sanding and the visual cross-dating of cores under a stereomicroscope (Stokes and Smiley 1996). Salix stems usually have circular cross sections (Lu et al. 2019b). Series of annual ring widths from mature individuals can be cross-dated well if mean series correlations exceed 0.32. The average series correlations of 0.34-0.60 for the cross sections of mature Salix shrub evidenced the reliable cross-dating or synchronization of the ring-width series.

Data analysis
We identified the age of each shrub individual from the annual ring counts on the 3 cm thick discs at the stem base (Lu et al. 2019a). The static age structure in each plot, which represented the tradeoff between shrub survival and mortality, was used to explore the dynamics of recruitment in shrubline plots (Myers-Smith and Hik 2018, Lu et al. 2019a). Recruitment will increase when shrub establishment increases more than mortality under climate change. Earlier analyses of the dynamics of treelines and shrublines have supported the notion that tree or shrub mortality would be low under favorable climate change (Camarero and Gutierrez 2004, Elliott 2012, Wang et al. 2015, Liang et al. 2016). Recent studies have used annual-resolution recruitment series from tree-ring data to analyze tree population dynamics in subalpine forests and treeline (Andrus et al. 2018, Bailey et al. 2021. To date, collecting annual-ring data is the only choice for the determination of the germination age of trees and shrubs in the studied remote regions. In the eastern TP, Salix shrub forms well-defined annual rings and can be used in dendroecological analyse (Lu et al. 2019b). We sampled all Salix shrubs located in each plot (Supporting information), so we could reliably reconstruct the yearly variations of shrub recruitment. Pearson correlation coefficients between shrub recruitment and mean minimum summer temperature in the 24 plots were used to evaluate the sensitivity of recruitment to temperature along the moisture gradient across the TP. We also used Pearson correlation coefficients between shrub recruitment and mean summer PDSI to evaluate the sensitivity of recruitment to drought along the gradient of soilmoisture content (depth of 10 cm).
Climatic records since the 1950s were available for most of the study sites on the TP. Data for monthly temperature and precipitation retrieved from gridded (0.5 × 0.5°) data of the Climate Research Unit (CRU) (Harris et al. 2014), however, were used for the sites where meteorological records were not available. A previous study evidenced that CRU data was highly correlated with the local data during the 1951-2000 period in China (Wen et al. 2006). In particular, the CRU data well captured the main features in the annual changes of temperature and precipitation (Wen et al. 2006). Thus, CRU data is considered as reliable climatic data source in the studies of shrubline-climate relationships on the TP. Estimated data for soil-moisture content  at each site (0.25 × 0.25°) was obtained from the NASA GLDAS-2 (<http:// disc.gsfc.nasa.gov>). The monthly PDSI, a dryness index calculated from temperature and precipitation, at each site was obtained from CRU gridded data (0.5 × 0.5° grids) (Harris et al. 2014). Atmospheric vapor pressure deficit (VPD, an indicator of evaporative demand)  at each site was calculated from CRU TS4.03 gridded data (0.5 × 0.5° grids). These climatic data sets were used to assess the relationships between shrub recruitment and climate change.
The series of shrub recruitments for each plot were converted to yearly shrub recruitment per hectare. Mean recruitment across the eastern TP was calculated by averaging recruitment for all plots. Monthly mean NDVI (Normalized Difference Vegetation Index) between 1982 and 2015 at each site retrieved from the GIMMS NDVI data set (0.083° resolution) (Pinzon and Tucker 2014) was averaged to generate a mean summer (June-August) NDVI across the eastern TP. Correlations between mean shrub recruitment and mean summer NDVI were calculated at the plateau scale (n = 24). If mean shrub recruitment and mean summer NDVI were significantly positively correlated, variations of shrub recruitment contributed to the greening of the eastern TP.
The sequential algorithm developed by Rodionov was used to detect abrupt changes in the spatiotemporal patterns of shrub recruitment (Rodionov 2004). This method is driven by data, in which an a priori hypothesis for the timing of shifts is not needed, and sequential t tests are performed on time-series data to identify changes in shifts (Rodionov 2004). An abrupt change was identified if the cumulative sum of normalized deviations from the average of a potentially new pattern of recruitment differed significantly from the average of the previous pattern of recruitment. Shifts in the means of the recruitment series were identified using a four-year cut length at p < 0.05. This method has also been successfully used to elucidate abrupt changes in recruitment data at treelines in the Rocky Mountains in the USA (Elliott 2012).
A dendroecological analysis indicated that the age of Salix shrubs 1 m tall varied between 10 and 30 years. We thus defined shrubline position as the uppermost elevation of shrubs, regardless of height. This definition is similar to that of previous studies of shrublines on the TP and in arctic tundra (Wang et al. 2015, Myers-Smith andHik 2018). Altitudinal variations of shrubline position in each plot based on this definition were reconstructed for the last 30-81 years over a 10-year interval using the locations and age structures of all shrubs.
Data for monthly temperature and precipitation, soilmoisture variables and spring and summer PDSI, which has negative values for dry conditions and positive values for wet conditions, near the 24 plots (Supporting information) were obtained from meteorological records, the CRU gridded database and the NASA GLDAS-2. Changes of mean minimum spring and summer temperatures, spring and summer precipitation and their changes, spring and summer PDSI, and mean seasonal VPD near each shrubline site were calculated using these three climatic data sets (Fig. 1, Supporting information).
Pearson correlation analysis was used to identify the relationships between the magnitude of shrubline shifts in the 24 plots during the past decades and various predictor variables across the TP. The potential predictor variables of shrubline migration in recent decades were: elevation of the current shrubline, mean slope and aspect of the shrubline ecotone, changes of mean minimum spring and summer temperatures, spring and summer precipitation and their changes, mean soil-moisture content in summer and the non-growing season (spring, autumn and winter), mean spring and summer PDSI, VPD in spring, summer, autumn and winter and the linear distances between the current shrubline positions and nearby mountain peaks. The decadal advance rates of shrublines were calculated and compared to warming velocity in recent decades. We tested if advance rate of shrublines before 2010 was higher than after 2010 using one-way ANOVA.
Structural equation models (SEMs) are frequently used to resolve complex multivariate relationships among interrelated variables using a form of path analysis (Lefcheck 2016). We used SEMs to identify the impacts of summer warming, summer PDSI, seasonal VPD and seasonal soilmoisture content at a depth of 10 cm on the advance of the shrublines before 2010 . Similarly, SEMs were used to link these factors and shrubline stability in the last decade (2011)(2012)(2013)(2014)(2015)(2016)(2017)(2018)(2019)(2020). The goodness of fit of the SEMs was evaluated using a chi-square test (p) and the root mean square error of approximation (RMSEA) in the selected SEMs. If p > 0.05 and RMSEA < 0.05 in the selected SEM, then the model fits the complex relationships among the observational data well (Lefcheck 2016). In such cases, the full model incorporated all the variables mentioned above may not perform better than a simpler model in which nonsignificant parameters were dropped. Effective recruitment (binomial variable: 1, 0) in each plot was defined as 1 if the shrubline in a plot shifted upward by at least half of the mean magnitude of shrubline advance in all plots. If not, the effective recruitment was set to 0. Effective recruitment was used to explore the links between shrub recruitment and shrubline shifts in the SEMs for the following reasons. First, the data for effective recruitment were significantly and positively correlated with the recruitment data during the advance of the shrublines in the past decades (r = 0.46 and p < 0.05 before 2010; r = 0.70 and p < 0.001 after 2010; n = 24 in both cases). Second, the data for effective recruitment were significantly correlated with shrubline dynamics of the past decades and several climatic variables (p < 0.05). Shrub recruitment, without considering the magnitude of the shifts, was not a significant predictor of shrubline dynamics. Similarly, tree recruitment, disregarding the magnitude of treeline shifts, did not predict the advance of treelines on the TP during the last century (Liang et al. 2016).

Magnitude of shrubline advances
We reconstructed the spatiotemporal variations in shrubline position in the past decades based on the ages of individual shrubs sampled in the 24 plots across the TP (Supporting information). The magnitude of Salix shrubline advances under substantial summer warming, infrequent seasonal drought and enhanced seasonal evaporative demand in recent decades (Fig. 1, Supporting information) varied between 0.10 and 59.3 m in elevation (mean = 18.24 m, n = 24) (Supporting information). The rates of advance of the Salix shrublines throughout the study period ranged from 0.0 to 12.2 m decade −1 (mean = 4.4 m decade −1 ), and the mean rate of summer warming was 0.27°C decade −1 . Based on an altitudinal lapse rate of 0.6°C per 100 m, the mean and maximum rates of upward shifts were equivalent to 0.026 and 0.073°C decade −1 , respectively. Based on the observed warming and approximate lapse rates, climatic envelope would theoretically permit an upward shift of the shrublines of 45 m decade −1 , which is much faster than the reconstructed shrubline advances. The mean upward shift of the Salix shrublines before and after 2010 was 5.16 and 2.34 m decade −1 , respectively. If we excluded two abnormal plots (HSX1 and HSX3), the mean advance rate of shrublines before 2010 was significant higher than after 2010 (One-way ANOVA: F = 7.24, p = 0.014, n = 22).

Shrub recruitment and climate
The maximum age of Salix shrubs along the TP latitudinal gradient varies from 30 to 81 years. Shrub recruitment tended to increase over time in each plot ( Fig. 2a-j). On average, recruitment in each plot peaked in 2010, but since then recruitment had been decreasing or was absent. Recruitment prior to 2010 was positively correlated with mean minimum summer temperature in 91.6% of the plots (Fig. 1,  Supporting information). In contrast, recruitment after 2010 has mainly been associated with moisture variables in all cases (spring and summer precipitation, the spring and summer PDSI, and seasonal mean VPD) (Supporting information). Mean minimum summer temperature and seasonal moisture variables thus had the strongest relationship with shrub recruitment before and after 2010, respectively. Spatial patterns in the sensitivity of shrub recruitment to climate along the moisture gradient differed between the two time periods. Before 2010, the sensitivity of shrub recruitment to mean minimum summer temperature tended to be positively correlated with the gradient in summer soilmoisture content along the TP (r = 0.24, p = 0.25; Fig. 3a), whereas the sensitivity of shrub recruitment to summer PDSI had a weaker tendency to be correlated with the summer soil-moisture content along the TP (r = 0.09, p = 0.65; Fig. 3b). In contrast, the sensitivity of shrub recruitment to mean minimum summer temperature after 2010 tended to be negatively correlated with the gradient of summer soilmoisture variables on the TP (r = −0.25, p = 0.24; Fig. 3a), whereas the sensitivity of shrub recruitment to summer PDSI after 2010 was significantly and positively correlated with the gradient of summer soil-moisture variables on the TP, with drought sensitivity higher at the wetter than the drier sites (r = 0.54, p = 0.005; Fig. 3b).

Drivers of shrubline shifts
Changes of mean minimum summer temperature (∆SUT) were the most important driver of the shrubline shifts in recent decades, whereas the moisture variables were also important. Among different climatic and non-climatic factors, ∆SUT and mean summer PDSI (SUPDSI) were positively and significantly correlated with the magnitudes of the shrubline advances in recent decades (∆SUT: r = 0.50; SUPDSI: r = 0.44; p < 0.05 in the two cases) (Supporting information).
Both mean shrub recruitment and mean summer NDVI at the TP scale increased and decreased significantly before and after 2010, respectively (Fig. 2k). Mean shrub density (individuals ha −1 ) was significantly and positively correlated with mean summer NDVI for 1982-2015 (r = 0.35, p = 0.043).

Forecasting shrubline shifts
We used SEMs to identify the impacts of potential key drivers (summer warming, mean summer PDSI, mean seasonal VPD and seasonal soil-moisture content) on Salix shrubline shifts before 2010 (1939-2010) and after 2010 (2011-2020). The full models incorporating all the variables showed poorer performances (before 2010: p = 0.0001 and RMSEA = 0.411; after 2010: p = 0.009 and RMSEA = 0.22) (Supporting information). Thus, we used the simpler models in which non-significant parameters were dropped (Fig. 4). The first simpler SEM (p = 0.618 and RMSEA = 0.001) indicated that summer warming positively affected shrub recruitment before 2010 and the advance of the Salix shrublines, whereas mean VPD during the non-growing season (spring, autumn and winter) negatively affected shrub recruitment before 2010 (Fig. 4a). The second simpler SEM (p = 0.768, RMSEA = 0.001) indicated that summer warming negatively affected shrub recruitment after 2010, while summer soil-moisture content positively affected shrub recruitment (Fig. 4b). Both of the two SEMs showed that shrub recruitment positively affected shrubline shifts before and after 2010 (Fig. 4). The third (p = 0.359, RMSEA = 0.032) and fourth simpler SEM (p = 0.872, RMSEA = 0.001) indicated that summer PDSI positively affected shrub recruitment and shrubline stationarity after 2010, whereas early growing season VPD (particularly June) negatively affected shrub recruitment (Supporting information). Mean summer VPD was negatively and significantly associated with mean summer soilmoisture content since 1954 near each study plot (p < 0.05) (Supporting information).

Discussion
Processes of recruitment (e.g., germination and seedling survival and establishment) of woody plants are strongly limited at their altitudinal or latitudinal limits of distribution by temperature during the growing season (Körner 2012). These limitations have been demonstrated in both observational studies and manipulative experiments. For example, sexual reproduction in shrubs near the shrubline increased substantially after 12 years of experimental warming (Klady et al. 2011). The alleviation of cold stress under a warmer climate will thus likely increase shrub recruitment by enhancing reproduction at the shrubline if soil moisture is high enough. Indeed, our study found increased recruitment before 2010 in all plots across the TP, regardless of topographical or latitudinal differences. This increase in recruitment indicated the strong impact of large-scale climatic forcing. Changes in minimum summer temperature among the climatic variables were correlated most strongly with the recruitment data, suggesting that summer warming was the main driver of shrub densification at these alpine shrublines. Our results are consistent with studies of the Pan-Arctic region (Tape et al. 2006, Rundqvist et al. 2011, Formica et al. 2014, Büntgen et al. 2015 and suggest that increased shrub recruitment could lead to extensive shrub expansion across the TP. Shrub recruitment tended to decrease after 2010, in strong contrast to this upward shift of the shrubline during the 20th century, likely due to the increasing drought stress during the warm season, which is generally creating less favorable moisture conditions. A sharp increase in land surface VPD and occurrence of drought events with rising temperature in the last decade on the TP suggested that precipitation conditions may not compensate for increasing evaporative demand induced by warming, thereby resulting in moisture stress (Ding et al. 2018, Wu et al. 2019, Yuan et al. 2019. These results could also be partly due to episodic masting, as observed in alpine treelines and arctic tundra (Elliott 2012, Körner 2012, Büntgen et al. 2015, Elliott et al. 2020. Extreme climatic events such as frosts could also inhibit the recruitment of these woody plants (Büntgen et al. 2015), but a lack of data from in situ ecological monitoring has prevented the quantitative study of such relationships at TP shrublines.
The sensitivity of shrub recruitment to temperature before 2010 increased with summer soil-moisture content, indicating that warming has increased shrub recruitment when moisture variable during the growing season was not limiting. The sensitivity of shrub recruitment to drought before 2010 varied little as a function of the summer soil-moisture gradient, suggesting that drought played a minor role in pre-2010 recruitment dynamics. The sensitivity of shrub recruitment to temperature after 2010 decreased with summer soil-moisture content, indicating that warming has decreased shrub recruitment when seasonal moisture variable might be limiting. The sensitivity of shrub recruitment to drought after 2010 increased with summer moisture gradients, suggesting that the responses of recruitment to seasonal drought stress tended to be higher at the wetter than the drier sites. For instance, relatively high drought sensitivity of recruitment in wetter sites such as HSX1-4 were observed in comparison with low drought sensitivity of recruitment in drier sites such as BD1-2. Indeed, previous studies have reported a higher sensitivity of woody plant growth to drought conditions at wet than dry sites (Cavin and Jump 2016). The plastic responses to drought of alpine plants could be higher in dry sites than in wet sites where a soil-moisture reduction may lead to strong impacts on ecological processes. The sensitivity of shrub recruitment to the climate in our study varied along the moisture gradient on the TP, indicating that spatial variations of seasonal moisture conditions modulated the population dynamics of these alpine shrublines under climate change.
Salix shrub densification and greening trends at the TP scale during the 20th century stopped or reversed after 2010. A greening trend in the 1980s and 1990s has been previously reported (Shen et al. 2015, Zhu et al. 2016 and suggests that infilling with alpine shrubs may play a role in the widespread greening on the TP (Shen et al. 2016). In parallel, evidence is accumulating that shrub densification is the main driver of vegetation greening across circumpolar regions (Myers-Smith et al. 2020). Our study identified abrupt changes in Salix shrub recruitment around 2010, which illustrated the importance of climatic and ecological thresholds in demographic processes across montane vegetation communities (Elliott 2012). Ample evidence suggests that sudden and potentially severe ecological effects of climate change can occur when temperature means and variability exceed a physiological inflection point (Trisos et al. 2020). The trend of decreasing Salix shrub densification and vegetation greenness since 2010 on the TP further indicated that alpine ecosystems may have already passed their climatic optima, that is further climatic warming would not promote recruitment without a simultaneous increase of the moisture availability (Liu et al. 2013, Lu et al. 2019a, Yuan et al. 2019. Accelerated warming would release plants from cold limitations but induce water shortage (Dolezal et al. 2020). The warmer and drier climate has already resulted in the recruitment decline of alpine shrubs since the 1990s in the central Himalayan Mountains (Sigdel et al. 2021). In parallel, owing to the hotter drought conditions, recent studies found no seedling establishment at treelines in the last decade along a 400 km latitudinal gradient in the Southern Rockies of USA, signifying that conditions are now beyond the climatic optimum for successful seedling establishment above treeline moving forward (Bailey et al. 2021). Remotely sensed data have provided evidence of this increasing drought stress in the last decade (Shen et al. 2015), and previous findings have also indicated that vegetation activity above 30°N is decoupling from climatic warming, with stronger associations with moisture conditions (Piao et al. 2014).
Our results indicated that the rate of advance of Salix shrublines did not match summer warming. The mean rate of advance of spruce and fir treelines on the TP is 2.9 m decade −1 (Liang et al. 2016), which is only 62% of the rate of advance of the Salix shrublines. The mean rate of advance of birch and fir treelines in the central Himalayas is 1.1 m decade −1 (Sigdel et al. 2018), which is also slower than that of the Salix shrublines. Our results are also comparable to those for different alpine plants in the European Alps, with rates of elevational shifts of 1-4 m decade −1 (Walther et al. 2002). The Salix shrublines are therefore more sensitive to climatic warming than are alpine treelines on the TP. We have demonstrated that Salix shrublines on the TP can migrate uphill to track climatic warming, supporting our initial hypothesis. Our results are consistent with previous studies that summer warming has already latitudinally expanded Salix shrublines in the circumpolar tundra and alpine areas of Eurasia and North America (Rundqvist et al. 2011, Formica et al. 2014, Myers-Smith and Hik 2018. The shrubline positions of evergreen shrub species (e.g., Juniperus pingii var. wilsonii (Rehder) Silba) on the central TP have remained relatively stable over the last two centuries (Wang et al. 2015), in contrast to the upward advance of deciduous Salix shrublines on the eastern TP. Such a sharp contrast could be partially attributed to their different functional traits associated with competition for light, resource acquisition and conservation (Mekonnen et al. 2018). A smaller foliar mass per unit area and less dense canopies in deciduous than evergreen shrubs lead to a greater interception of light for a given investment in foliar carbon (Mekonnen et al. 2018). The superior competitiveness for light, higher rates of growth and rapid decomposition and nutrient cycling of deciduous shrubs may allow Salix to rapidly colonize newly opening habitats, unless hampered by natural disturbances or moisture stress (Mekonnen et al. 2018, Myers-Smith andHik 2018).
Our fitted models provided insight into the direct and indirect relationships between hydrothermal factors and shifts in Salix shrublines during two periods (1939-2010 and 2011-2020). Summer warming during 1939-2010 likely stimulated shrub recruitment by alleviating cold stress. Increased shrub recruitment due to the combined effects of summer warming and higher availability of soil moisture on the TP eventually led to upward advances of the Salix shrublines between 1939 and 2010. However, reduced shrub recruitment due to enhanced evaporative demand during the non-growing season may slow the rates of shrubline advance, thus largely explaining the mismatch between the rates of Salix shrubline advance and warming velocity in the past decades. This finding was partly consistent with result of a meta-analysis of treeline shifts in the Northern Hemisphere where autumn moisture conditions regulated the rates of alpine treeline advances under climatic warming (Lu et al. 2021). The negative effects of summer warming and evaporative demand during the early growing season and the positive effect of summer PDSI and summer soil-moisture content on shrub recruitment during the 2010s suggests that a critical thermal threshold was reached and passed, introducing a sudden shift in conditions (Elliott 2012). Due to warming-induced moisture limitation in the growing season, Salix shrublines on the TP remained relatively stable in the recent decade. Warming-recruitment links varied from positive before 2010 to negative coupling in the last decade, while the negative impact of evaporative demand on shrub recruitment shifted from the non-growing season before 2010 to the early growing season in 2010s, indicating that a hotter-drier climate could have strengthened the coupling of land-atmosphere interactions and ecological processes of alpine plants.
In summary, this field study is the first to quantify the dynamics of Salix shrublines along a wide latitudinal gradient across the TP spanning ca 900 km. The altitudinal expansion of Salix shrublines into previously shrubless areas in recent decades is among the most striking ecological changes on the TP. Shrub infilling and shrubline advance induced by summer warming was evident before 2010, but moisture limitation, which acted as an important environmental modulator of changes to alpine vegetation, has recently slowed the rate of Salix shrub recruitment and shrubline advance under climatic warming. The dependence of the dynamics of alpine shrublines on moisture highlights the importance of in situ monitoring of soil-moisture variables and of controlled experiments and plant-water relations in studies of alpine ecosystems (Buchwal et al. 2020). Salix shrublines under scenarios of future increased warming and aridity are at risk of downhill shifts due to the negative impacts of drought stress on key physiological, demographic and community-scale processes (Wernberg et al. 2016). Such changes will have important implications for alpine ecosystem services, vegetation/climate feedbacks and bush management. Climate/shrubline relationships should thus be better represented in the next generation of Earth system models because previous models have inadequately represented the effects of temperature and moisture limitations on vegetation activity, particularly in mountainous regions. Additionally, tree-ring-based reconstructions of heatwaves and soil moisture for the past 260 years reveal an abrupt shift to hotter and drier climate in inner East Asia . Buchwal et al. (2020) likewise reported the declines in the greenness and production of tundra vegetation over a considerable part of the terrestrial Arctic in response to warming and drying associated moisture limitation. These findings combined with our results indicate that warming and drying as a large-scale phenomenon have led to an abrupt shift in climate (Zhang and Fang 2020) and would have enhanced a land-atmosphere-ecosystems coupling, which warrants further studies.