《1. Introduction》

1. Introduction

As the third pole of our planet, the Tibetan Plateau in China is the highest and most extensive alpine region in the world [1]. It is also a hotspot of climate change, with an observed warming rate (Fig. 1(c)) twice that of the global average (0.25–0.27 °C per ten years) [2]. Its exposure to high levels of solar radiation distinguishes the Tibetan Plateau from other high-latitude regions with a similarly cold climate. The plateau is home to a population of three million people, half of which live on agriculture [3]. The cropping system across Tibet is also quite unique. The traditional staple food of Tibetan people is highland barley (Hordeum vulgare L., or ‘‘Qingke” in Chinese), which is the only crop that can grow over 4000 m above sea level; it accounts for about 60% of the cropgrowing area and about 70% of the cereal production within the Tibet Autonomous Region (TAR) in China. The cultivation and production of highland barley also nurture the unique Tibetan culture. A simultaneous increase in yield and population in the TAR during the past three decades (Fig. 1) has significantly intensified both food supply and demand. The tightening demand–supply balance, which is due to the leveling of yield increment in the recent decade (Fig. 1), along with growing evidence of globally widespread negative warming impacts on crop yield [4–6], should have raised concerns regarding the impacts of climate change on Tibetan crop yield. However, a lack of evidence on the impacts of rapid warming on Tibet limits the current understanding of Tibetan crop yield response to warmer climate. In this research, using crop yield statistics (1985–2015) and historical climate datasets, we investigate the yield–climate relationships in the TAR and their evolution during the past three decades.

《Fig. 1》

Fig. 1. (a) Highland barley in the Tibetan Autonomous Region (TAR) in China; (b) a time series of human population and crop yield in the TAR; and (c) a time series of anomalies in growing-season temperature and precipitation for barley-growing areas of the TAR. GST: growing-season temperature; GSP: growing-season precipitation.

《2. Methods》

2. Methods

《2.1. Datasets》

2.1. Datasets

Time series of growing area and production for cereals in the TAR between 1985 and 2015 were obtained from the Tibet Statistical Yearbook. The dominant crop grow over the TAR is highland barley, which accounts for ~60% of crop growing area and ~70% of cereal production. Although the data for highland barley were only available for 1985–1994, the variability of cereal yield across the TAR generally reflects the variations in highland barley (R2 = 0.91, P < 0.01; Fig. S1 in Appendix A). Therefore, our analysis on cereal yield variability mostly reflects that of barley yield. The barley-growing area was obtained from the MIRCA2000 dataset [7] (Appendix A Fig. S2). Monthly climate data (including mean/ maximum/minimum air temperature, precipitation, and solar radiation) were obtained from 0.1° gridded climate data from the Chinese Academy of Sciences of China Meteorological Forcing Dataset (CMFD) [8]. The growing season was defined as May to August, and the growing-season temperature (GST), growing-season precipitation (GSP), and growing-season radiation (GSR) for each grid of barley-growing area were calculated by means of area-weighted averaging of the corresponding monthly data. The annual mean GST, GSP, and GSR across the entire TAR were obtained by weighting each grid cell (0.5° × 0.5° grids) according to its cereal growing area [9]. Temperature data for the 1.5 and 2.0 °C above-preindustrial-level warming scenarios were derived from the biascorrected climate change projection by Institut Pierre–Simon Laplace Coupled Model, version 5, coupled with NEMO, low resolution (IPSL-CM5A-LR) [10], in which simulations extend from the pre-industrial period to the end of the 21st century. Temperature changes across the TAR under different warming scenarios were also obtained from IPSL-CM5A-LR.

《2.2. Data analyses》

2.2. Data analyses

To analyze the relationship between yield and climate, we first detrended all the time series based on a common approach of first difference (that is, year-to-year changes) [11,12]. The use of first differences minimizes the influence of slowly changing factors such as crop management, technology advances, and rising atmospheric carbon dioxide (CO2). Simple correlations between the detrended time series were used to analyze the relationships between yield and climate during 1985–2015. To minimize the confounding impacts of co-varying variables, we also performed partial correlations between yield and one climatic variable while statistically controlling for the other two variables. Temporary changes in temperature–yield relationships were examined by repeating both the ordinary and partial correlations using a moving time window of 15 years during 1985–2015. The robustness of the time window length choice was tested by repeating the same analyses using different time window lengths ranging from 10 to 17 years (Appendix A Fig. S3). Correlation analyses are generally more reliable with longer time series. We present the results from the 15-year window in the main text because it is the longest time window that ensures that the data points in the first and last time windows are fully independent. The sensitivity of yield to temperature (ST) was obtained from a multiple regression between the yield and the GST, GSP, and GSR. Regressions between ST and GST were performed in order to explore possible reasons for the temporal changes in ST. Using the regression model, we also extrapolated ST into different climate change scenarios (1.5 and 2.0 °C) with 95% confidence intervals.

《3. Results》

3. Results

《3.1. Crop yield variations are primarily driven by temperature anomalies》

3.1. Crop yield variations are primarily driven by temperature anomalies

First, we analyzed the interannual variations of the yield (Y), growing-season temperature (GST), precipitation (GSP), and incoming solar radiation (GSR) over the past three decades (Fig. 2), obtained through detrending with the first difference (see Section 2) [11,12]. Across the entire period of 1985–2015, large anomalies in yield exhibited a strong anti-phase with those of temperature (Fig. 2(a)). Correlation analyses indicate that Y has a significant negative correlation with GST when controlling for GSP and GSR (RY–GST = –0.37, P = 0.05; Fig. 2(b)). On average, a one-degree warmer temperature anomaly results in a yield loss of 0.11 t·ha-1 . On the other hand, Y also appears to correlate with GSP and GSR (Figs. 2(e) and (f)); however, these correlations disappear after controlling for other climatic factors in partial correlation analyses (P > 0.05; Figs. 2(e) and (f)). Together, the results suggest that GST is more important than GSP and GSR in driving variations of yield in the TAR, which is reasonable, as Tibetan croplands are exposed to abundant solar radiation (mean annual solar radiation ~7.2 × 109 J·m–2 ·a–1) and are widely managed with irrigation [13].

《3.2. Crop yield appears to be more sensitive to warmer temperature》

3.2. Crop yield appears to be more sensitive to warmer temperature

As yield variations are primarily driven by temperature anomalies (Fig. 2), which rose rapidly over the past three decades (Fig. 1), the next question is whether the response of crop yield to warmer temperature has changed over the same period. To answer this question, we analyzed the relationship between Y and GST (Fig. 3(a)) over 15-year moving time windows. As Fig. 3 shows, the partial correlation between Y and GST (RY–GST) was not statistically significant during the late 1980s and 1990s, but became significantly negative during 2000s (RY–GST < –0.60, P < 0.05). This change in the correlation coefficient between Y and GST is also robust when GSP and GSR are controlled for. Since standard deviations of Y and GST do not change significantly across time windows (Appendix A Fig. S4), this increasing negative correlation between Y and GST should come from enhanced negative responses of crop yield to warmer temperature. Indeed, multiple regression of Y against GST, GSP, and GSR indicates that a 1 °C increase in GST would reduce yield by (0.13 ± 0.20) t·ha-1 in the first 15-year time window; however, the yield loss responding to the same amount of GST increase would almost double ((–0.22 ± 0.14) t·ha-1 ) during 2001–2015.

《Fig. 2》

Fig. 2. (a, c, e) Time series and (b, d, f) relationships between first differences of yield and (a, b) GST, (c, d) GSP, and (e, f) GSR in the TAR. Rsimple and Rpartial represent ordinary and partial correlation coefficients, respectively.

We further performed two additional analyses to assess the robustness of the change in the partial correlation coefficient between Y and GST (RY–GST) and the doubling of the apparent ST variations over the past three decades. First, we examined whether data from a few extreme years might have led to the observed changes in RY–GST and ST. By performing 500-time boot-strapping analyses, we found that RY–GST robustly changed from –0.42 ± 0.21 in the first 15-year period to –0.70 ± 0.17 in the last 15-year period. The difference between the RY–GST of the two periods is significant (P < 0.01). Similarly, ST significantly (P < 0.01) increased from (–0.13 ± 0.07) t·ha-1 ·°C–1 during the first 15-year period to (–0.22 ± 0.06) t·ha-1 ·°C–1 in the last 15-year period (Appendix A Fig. S5). Thus, the increasing negative RY–GST and ST are not caused by the data from a few extreme years. Second, we tested whether the observed changes in RY–GST and ST were artefacts of selected time window lengths. We performed the same moving window analyses with a time window length varying from 10 to 17 years. The results show that, regardless of the time window length, RY–GST changes from statistically nonsignificant (P > 0.10) in the first time window to significantly negative (RY–GST < –0.68, P < 0.05) in the last time window, except for 10-year time windows, in which RY–GST is marginally significant in the last time window (RY–GST = –0.62, P < 0.10). The average increment in ST between the first and the last time windows across different time window lengths is 100% ± 49%, ranging from 55% (17-year time windows) to 177% (12-year time windows) (Appendix A Fig. S3). Thus, both the increasing negative RY–GST and the almost doubling ofST over the past three decades are robust to the choice of time window length.

《3.3. An alarm from the emerging negative impacts of warming on the TAR》

3.3. An alarm from the emerging negative impacts of warming on the TAR

The magnitude of ST (–2%·°C–1 to –4%·°C–1 ) in the TAR during the past three decades is less than the yield loss of global barley to temperature increase (–9%·°C–1) [12]. However, the emerging negative impacts of warmer temperatures on Tibetan crop yield are particularly alarming, given previous assumptions that crop yield at high latitudes and altitudes may more or less benefit from warmer temperatures [14]. The findings may also imply that the contemporary benefits of warming on crop yield in some wet and cold regions [15] may soon disappear in the near future. Furthermore, the shift of ST from nonsignificant in the 1980s and 1990s to significantly negative in the 2000s can be associated with the rapid rise in GST (Fig. 1(c)). Figs. 3(c) and (d) show that RY–GST is strongly correlated with GST (R = –0.78, P < 0.001). This relationship does not vary much (R = –0.58, P = 0.03) when GSP and GSR are statistically controlled for. Although GSP also shows a large increasing trend (3 mm·a–1, P < 0.05; Fig. 1(c); insignificant trend of GSR in Appendix A Fig. S6), the changes in GSP and GSR are not significantly correlated with those in RY–GST in partial correlation analyses (Appendix A Fig. S7).

《Fig. 3》

Fig. 3. Temporary changes in correlation coefficients between first differences of yield and GST, and yield sensitivity to GST in the TAR. (a) Temporary changes in correlation coefficients (Rsimple and Rpartial represent ordinary and partial correlation coefficients, respectively; dashed dots represent 10% significance levels); (b) temporary changes in the GST sensitivity of crop yield (ST represents the sensitivity derived from multiple linear regressions); (c) the relationship between Rpartial and GST; (d) the relationship between ST and GST.

《4. Discussion》

4. Discussion

We propose two possible mechanisms that may explain our findings. First, given the exponential relationship between saturated water pressure and temperature [16], the same amount of temperature increase could induce a greater atmospheric water deficit in a warmer climate. Such a greater water vapor pressure deficit (VPD) can reduce stomatal conductance and thus photosynthesis [17], leading to declining crop productivity. Distinguishing direct warming impacts from indirect impacts by modifying the atmospheric water demand is not easy. However, partial correlation analyses between Y and variations in the growing-season average maximum daily temperature (Tmax) and minimum daily temperature (Tmin) showed that the partial correlation between Y and Tmin was not significant in most recent decades, and the partial correlation coefficients remained relatively stable across the study period (Fig. S8). On the contrary, the partial correlation between Y and Tmax was significant (P < 0.05) in the recent decades, and the negative partial correlation between Y and Tmax was found to be strengthening over time (Fig. S8), which explained the observed temporal changes in RY–GST. The more dominant role of daytime temperature (Tmax) in comparison with nighttime temperature (Tmin) in the yield–temperature relationship indicates that temperature impacts on photosynthetic processes, rather than on respiratory processes, are what drive the change in RY–GST, which is consistent with our first hypothesis that warming-induced higher VPD stresses crop productivity in the TAR.

Second, although the mean GST is rather low in the TAR, summer (July and August) daytime temperatures can still be quite high (> 25 °C) (Appendix A Fig. S9). A recent study of the cardinal temperature thresholds indicated that, even for the Tibetan croplands with their relatively lower altitude, the optimum temperature for photosynthesis is still less than 25 °C [18]. This finding indicates that there is significant temperature stress beyond the optimum temperature for barley development [19], which could negatively affect the yield [20]. In fact, summer in Tibet coincides with highland barley’s reproductive growth period, which is known to be the period that is particularly sensitive to heat stress [21–23]. Rising temperatures result in an increase in both the intensity and frequency of hot days, which may lead to stronger yield decline [24–26]. Nonetheless, information on the reproductive growth of Tibetan highland barley is still very limited, which hinders us from narrowing down the exact point stress that is predominantly responsible for the negative warming impacts. Future studies should enhance the monitoring of phenological and growth indicators of Tibetan highland barley, and manipulative warming experiments on different growing periods are encouraged in order to further understand the mechanisms driving the increasing negative yield response to warmer temperatures.

Deducing from the increasing yield sensitivity to GST under a warmer climate over the past three decades, we expect a stronger negative response of crop yield variations to temperature variations in the even warmer future. By extrapolating the historical relationship between GST and ST to projected climate change, we found that ST may change from (–0.19 ± 0.04) t·ha-1 ·°C–1 , that is, (4.3%·°C–1 ± 0.9%·°C–1 ) under contemporary climate, to (–0.33 ± 0.10) t·ha-1 ·°C–1 (–7.5%·°C–1 ± 2.3%·°C–1 ) under the 1.5 °C above pre-industrial levels scenario, and to (–0.51 ± 0.18) t·ha-1 ·°C–1 (–11.6%·°C–1 ± 4.1%·°C–1 ) under the 2.0 °C warming scenario(Fig. 4). This means that, even if the climate target set by the Paris Agreement (the 1.5 °C warming scenario) [27] can be achieved, the sensitivity of Tibetan crop yield to temperature will still increase to almost twice that of the past three decades. Previous studies assuming an unchanged ST over time [4,28] may have drastically underestimate future warming impacts on crop yield. With continuous warming, the vulnerability of Tibetan crop yield to the projected more frequent extreme heat events will also increase, putting food security and the unique culture of the Tibetan people in danger.

《Fig. 4》

Fig. 4. Predicted ST in the TAR under different warming scenarios. The left and right parts represent the historical and predicted future relationships between ST and GST. The grey and red areas represent the 95% confidence intervals for historical and future estimates of ST, respectively. The two black vertical lines represent the scenarios of 1.5 °C and 2.0 °C warming above pre-industrial levels, respectively.

《5. Conclusions and future perspectives》

5. Conclusions and future perspectives

To summarize, we found emerging negative responses of crop yield to temperature change in the TAR during the past three decades. The apparent sensitivity of yield to GST approximately doubled. This finding is in contrast to previous studies suggesting weakening impacts of temperature variations on maize yield across the United States [29], and implies that increasing precipitation or atmospheric CO2 may not mitigate negative warming impacts, at least in the TAR. Our analyses also call into question the often-used assumption of constant crop yield sensitivity to temperature in predicting future crop yield in response to climate change. Contemporary literature has largely focused on the ‘‘big four” crops; however, our results provide new insights into a crop system that is both biologically and culturally unique. While our understanding of how the unique crop system in the TAR may respond to climate change still contains large uncertainties—particularly when compared with process modeling and other approaches—it provides an additional line of evidence for constraining crop model projections [30]. Our finding highlights the urgency of further experiment and modeling efforts, in order to ensure regional/global food security and the lifestyle of the Tibetan people across the world’s third pole.

《Acknowledgments》

Acknowledgments

This work was supported by the Second Tibetan Plateau Scientific Expedition and Research Program (2019QZKK0405) and the National Natural Science Foundation of China project Basic Science Center for Tibetan Plateau Earth System (41988101).

《Authors’ contribution》

Authors’ contribution

Shilong Piao and Xuhui Wang designed the research. Chuang Zhao performed the statistical analyses. Tsechoe Dorji, Chuang Zhao, Shilong Piao, Xuhui Wang, and Baohua Liu drafted the manuscript. All authors contributed to interpretations of the results and discussions of the contents.

《Compliance with ethics guidelines》

Compliance with ethics guidelines

Tsechoe Dorji, Shilong Piao, Xuhui Wang, Chuang Zhao, Baohua Liu, Anping Chen, Shiping Wang, and Tao Wang declare that they have no conflict of interest or financial conflicts to disclose.

《Appendix A. Supplementary data》

Appendix A. Supplementary data

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