Elevation‐Dependent Warming Over the Tibetan Plateau From an Ensemble of CORDEX‐EA Regional Climate Simulations

Under the Coordinated Regional Climate Downscaling Experiments‐East Asia (CORDEX‐EA‐II), the outputs from two regional climate models (RCMs) driven by four global climate models (GCMs) are used to investigate the characteristics and possible mechanisms of the projected elevation‐dependent warming (EDW) over the Tibetan Plateau (TP) under the Representative Concentration Pathway emission scenario 8.5 (RCP8.5). Results show that widespread warming over the TP is projected with considerable disagreements in warming intensity and the maximum warming center among RCMs. The largest spread in the surface air temperature (Tas) projections is found above 5,000 m, indicating that a large uncertainty exists over the higher elevations. A marked EDW signal over the TP is simulated under the RCP 8.5 by the multi‐RCM ensemble mean for all seasons, particularly in autumn. Based on the analysis of the surface energy budget, it is found that the surface albedo feedback (SAF) is the primary contributor to EDW and acts as the main source of uncertainty in EDW projections among RCMs. The downward longwave radiation (DLW) is found to be the dominant factor in regulating Tas change over the TP, and its contribution to EDW is model‐dependent. Furthermore, the structure and magnitude of projected EDW are sensitive to the RCM physics and driving GCM, as they can alter the projections of snow cover and albedo, which modulate the simulated SAF and its effect on EDW. Additionally, RegCM4 shows a higher sensitivity to the anthropogenic greenhouse forcing than WRF, evidenced by the larger temperature projections and stronger EDW signal in RegCM4.

tools in reproducing fine-scale processes over complex terrains, with their better resolving heterogeneous topography and solid parameterization physics (Fu et al., 2005;Gutowski et al., 2020). In addition, the RCMs may generate different projections compared to their driving GCMs (Niu et al., 2015;Saini et al., 2015).
According to previous studies, RCMs exhibit an advantage over GCMs in reproducing temperature over the TP with systematic cold biases (Gao, Duan, et al., 2019;Xu et al., 2018a), and they can capture the relationship between warming rate and elevation (Gao et al, , 2018Guo, Yu, & Wang, 2016). Note that the RCM simulations over the TP are extremely sensitive to the model resolution (Lin et al., 2018;Xu et al., 2018a), the driving forces (Gao, Chen, Lettenmaier, et al., 2018), and physical parameterizations applied in the model (Gu et al., 2020;Wang et al., 2016), which lead to different climate projections among RCMs.
For future climate change over the TP, previous RCM studies focus mostly on the seasonal climate projections and find that the magnitude of projected warming is dependent on the RCM and scenario (Sanjay et al., 2017;Wu & Gao, 2020). Nevertheless, only a few use RCMs to investigate future EDW over the TP. Guo, Yu, and Wang (2016), for example, analyzed future temperature projections over the TP using two RCMs at 50 km resolution. They found that enhanced warming is expected to happen in the low-elevation range, and attributed it to the snow changes with elevation. Another study using the Weather Research and Forecasting Model (WRF) at 30 km resolution (Gao, Chen, Lettenmaier, et al., 2018) also suggested a future EDW below 5,000 km over the TP, which can be explained by the surface albedo feedback (SAF), water vapor changes and diabatic heating release. Considering that the character and causes of simulated EDW depend strongly on the adopted RCM configuration (Minder et al., 2018), the future EDW over the TP has not been fully understood and further investigations based on multi-RCM ensemble projections are urgently needed.
To better understand Regional climate system and produce coordinated sets of regional downscaled projections worldwide, the Coordinated Regional Climate Downscaling Experiment (CORDEX) was implemented in 2009 under the auspices of the World Climate Research Program (Giorgi, Jones, & Asrar, 2009). As the East-Asian branch of the CORDEX initiative, the first and second phase of CORDEX-East Asia (CORDEX-EA-I and CORDEX-EA-II) have generated an ensemble of regional climate change projections by downscaling several GCMs using multiple RCMs (Tang et al., 2018;Wang et al., 2019;Zou et al., 2019). Under the framework of the CORDEX-EA-II, we design eight experiments using two RCMs driven by four GCMs to address the following issues: (1) the character of future EDW over the TP by multi-RCM projections; (2) identification of mainly physical processes for future EDW over the TP; and (3) effect of driving GCMs and RCM physics on seasonal temperature projections and EDW.
The study is structured into the following sections. Section 2 introduces the model experiments and observational data as well as the methods used for the analysis. In Section 3, we describe the characteristics of projected temperature changes and EDW over the TP, and discuss the potential mechanism responsible for EDW. The influences of driving GCMs and RCM physics on those projections are also described. Finally, the summary of this study and discussion are provided in Section 4.

Experimental Design and Models
Four GCMs from the CMIP5 ensemble, including Centre National de Recherches Météorologiques Earth system model (CNRM-CM5), the state-of-the-art Earth System Model (EC-EARTH), Geophysical Fluid Dynamics Laboratory (GFDL-ESM2M), and Earth system model (MPI-ESM), are chosen as the driving forces for the historical climate simulations and future climate projections. The detailed information about the GCMs is listed in Table 1.
To generate the TP climate change scenario at high resolution, the RCM version 4 (RegCM4, Giorgi, Coppola, et al., 2012) and WRF model (Skamarock et al., 2008), whose configurations are listed in Table 2, are nested within four GCMs for the historical and future periods under the Representative Concentration Pathway emission scenario 8.5 (RCP8.5). All simulations are integrated over the CORDEX-EA-II domain ( Figure 1) with a uniform horizontal resolution of 25 km. Two periods are chosen to represent the historical  and future (2029-2060) climate, and a two-year spin-up time is applied in the two periods.

Datasets and Methods
The performance of simulated seasonal surface air temperature (T as ) over the TP is evaluated against the gridded temperature dataset with 0.25°×0.25° resolution provided by the National Climate Center of China Meteorological Administration (NCC/CMA). This dataset (hereafter referred to as CN05, Wu and Gao, 2013) is obtained from 2,416 weather monitoring stations using the anomaly interpolation method described in New et al. (2000), and has been widely used to describe the detailed information of regional climate change and validate high-resolution models (Bao et al., 2015;Li, Zhou, et al., 2018a;You et al., 2018). Due to the harsh environment, the in-situ stations over the TP have sparse distribution and are mainly located over the eastern and central TP, limiting its reliability over the TP. Therefore, a satellite-based surface air temperature dataset with a resolution of 1 km (Xu et al., 2018b) is also used in this study. This dataset is produced by machine learning models based on 11 environmental variables derived from Moderate Resolution Imaging Spectroradiometer (MODIS) data, Shuttle Radar Topography Mission (SRTM) digital elevation model (DEM) data, and topographic index data, and is available from the Shallow Water Earth Observation Lab (https://www.shallowwaterlab.com/). However, this dataset (hereafter referred to as MODIS) is available since 2001 and thus only used to assess the climatic distributions.
For the convenience of model evaluation, both simulation and observation data are interpolated into 0.25 ° × 0.25 ° resolutions using bilinear interpolation. Taylor diagrams (Taylor, 2001), which present a concise statistical summary in terms of spatial correlation (CC), centered root mean square (RMSE), and spatial variance, are applied to examine the models' performance. The locally weighted scatterplot smoothing (LOW-NIU ET AL.  ESS) approach is also used in this study to reduce the impact of the uneven distribution of points at different altitudes. Due to distinct seasonal climatic conditions over the TP, the evaluation of models is carried out for the four seasons, namely spring (March-May = MAM), summer (June-August = JJA), autumn (September-November = SON), and winter (December-February = DJF).
A surface energy budget approach proposed by Lu and Cai (2009) is employed to attribute temperature changes to surface energy components. According to the land surface energy balance equation and Stefan-Boltzmann law, the land surface temperature (Ts) change can be written as follows: where σ is the Stefan-Boltzmann constant (5.67 × 10 −8 Wm −2 K −4 );  LW , and  LW represent the surface upward and downward long-wave radiation;  SW and  refer to the downward short-wave radiation and surface albedo; GHF, SH, and LH denote the ground heat flux, surface sensible, and latent heat fluxes, respectively. Additionally, the overbar stands for the historical mean climate and Δ means the climatic differences between the future and historical period. Therefore, the temperature change can be decomposed as: where the six terms from left to right on the right-hand side represent temperature changes induced by the SAF, downward short-wave radiation (DSW), downward long-wave radiation (DLW), ground heat flux (GHF), surface sensible (SH), and latent heat (LH) flux, respectively. Figure 2 presents the distribution of seasonal mean T as from CN05 and the differences between MODIS and CN05 and model biases with respect to CN05. Similar distributions can be found between CN05 and MODIS, with the Qaidam Basin and the southern edge of the TP showing higher temperature and the northwest areas displaying lower temperature (Figures 2a and 2b), which reflects the importance of the topography and natural climatic setting in regulating T as distributions. Compared to CN05, MODIS has higher magnitudes reaching up to 6°C over the northwest part of the TP. All RCMs and their ensemble mean (MME) generally reproduce the observed spatial distribution of T as well, but they always generate cold biases, with the largest biases of −10°C locating over the northwest TP. The underestimation of T as in model simulations over mountainous regions is a well-known deficiency for models (Giorgi, Bi, & Pal, 2004), and maybe related to the weakness in the model physics, such as the stronger surface albedo feedback (Chen et al., 2017), misrepresentations of orographic drag (Lin et al., 2018) and lack of precipitating ice (Lee et al., 2019), etc. Furthermore, the gridded dataset we used in this study is constructed from sparse in-situ stations over the TP, especially for the northwest TP, which may exaggerate model bias (Lucas-Picher et al., 2011).

Evaluation of the Present-Day Simulation
The general performance of RCMs in simulating T as over the TP is summarized in Figure 3. Statistically, all simulations demonstrate a significant correlation with CN05 (larger than 0.8) and have relatively higher spatial variability compared to CN05. Due to different physics in the driving GCMs and RCMs, considerable spread exists among the RCMs, especially during JJA. Note that the simulated temperature distribution is more sensitive to the internal model physics than driving GCMs during the cold season (MAM and DJF), while during the warm season (JJA and SON) it relies on both the driving GCMs and the internal model physics. Additionally, the RegCM4 models simulate better T as distributions than WRF models in MAM and DJF, as higher CC and lower RMSE are found in RegCM4 models. NIU ET AL.

Warming Patterns Under RCP8.5
Widespread warming over the TP is projected by all models in the middle of the 21st century (2031-2060) under RCP8.5. The projected change in T as by MME shows slight seasonal variations, with the largest change of 2.42°C in SON and the lowest of 1.97°C in MAM. However, there is a large disagreement in the seasonal T as projections among the RCMs (Table 3). For example, EC-EARTH_W and MPI_W generate the greatest warming in JJA, while other RCMs in SON. Even driven by the same GCM, MPI_W shows the lowest change in SON when the largest change occurring in MPI_R. Furthermore, the distribution of T as change follows the topography over the TP in most models' projections, with large inter-model variability for the magnitude and location of the warming center ( Figure 4). For instance, the maximum warming is located over the central of TP in MPI_W's projections while over the northwest TP in other models during JJA. Those disparities in the T as projections among the RCMs underline that large uncertainties exist in T as projections over the TP. Such disagreements arise from the differences in large-scale circulation produced by GCM and the discrepancy in representations of mesoscale processes in RCMs.
To quantify the uncertainty in the T as projections over the TP, the spread of multi-RCM, represented by the standard deviation of the eight RCMs, is shown in Figure 5. The spread in T as projections is mainly between 0.2°C and 0.5°C, with the largest above 1.5°C lying over the southeast TP during MAM and DJF and over the southwest TP during JJA and SON. Comparatively, the T as projections among RCMs are more consistent in cold seasons when large scale circulation dominates than in warm seasons when the meso-microscale convective systems always occur, which is consistent with Wu and Gao (2020). Interestingly, the LOWESS fit suggests a close relationship between the spread and elevation for all seasons. The spread in JJA and DJF remains constant up to about 4,500-5,000 m and then increases rapidly. While for other seasons, it increases slowly until 4,000 m, then decreases from 4,000 to 5,000 m, and finally increases sharply above 5,000 m. Such disagreement in temperature projections over the high altitudes mainly derives from the differences in treatment of snow cover and accompanied SAF among models (Qu & Hall, 2014;Winter et al., 2017). Therefore, reliable future projections over mountainous regions require further improvement of the components affecting the SAF in the model, such as surface albedo, snow cover and vegetation characteristics, etc. Additionally, the internal variability of each experiment, affected by the RCM physics and driving GCM, increases with elevation under increasing greenhouse-gas concentrations, which in turn leads to increased uncertainty at high-elevation regions (Dimri et al., 2018).

EDW in the Future
A clear EDW signal over the TP during recent decades has been identified by previous studies based on the observations (Guo, Sun, et al., 2019;Yao et al., 2019) and simulations (Palazzi, Filippi, & von Hardenberg, 2017;Rangwala et al., 2009), which accelerates the rate of change in cryosphere system, mountain ecosystems and biodiversity (Pepin et al., 2015). Therefore, whether this signal will continue in the future urgently requires further investigation.
As is shown in Figure 6, a marked EDW signal under RCP8.5 is projected in the MME for all seasons with the strongest occurring in SON, corroborating previous modeling study (Palazzi, Mortarini, et al., 2019). However, the relationship between the temperature change and elevation de-NIU ET AL.   pends on the season. The warming rate in MAM tends to increase with height up to 4,500 m followed by a slight decrease thereafter. Similar features in warming rate can be found for JJA and SON, but with higher magnitudes and the maximum around 5,000-5,500 m. In DJF, the EDW signal is expected to appear below 4,500 m. Like the glacier, snow and permafrost are mainly located above 4,500 m over the TP (Guo, Sun, et al., 2019), a larger projected warming rate over the higher elevation regions may result in the glacier melting and permafrost thawing as well as the changes in the fraction of solid and liquid precipitation, thereby altering river flow regimes the streamflow patterns and associated water resources.
The spatial variability of the T as projections by MME is almost less in the cold season for the entire altitudinal range compared to that in the warm season ( Figure 6), again demonstrating the importance of the meso-microscale convective systems in temperature projections. Interestingly, the largest spatial variability about 0.4°C lies around 3,500-4,000 m and the lowest occurs around 5,500 m in DJF. While in JJA and SON, the spatial variability tends to increase as the elevation increases, with the maximum reaching up to 0.5°C above 5,500 m. Those diverse responses in the grid points within the same elevation bins to global NIU ET AL.
10.1029/2020JD033997 7 of 17 warming may be related to the local factors, such as land-use type, topography, slope, etc (Pepin et al., 2015). Moreover, larger spatial variability over higher elevations may be attributed to a wide range of altitudinal variation within short distances over there, which results in a sharp gradient of vertical changes in the lapse rate (Barry, 2008). Besides, the cloud condensation can increase the diabatic processes in the mid and high troposphere and hence the temperature variation at high altitudes (Ohmura, 2012).
Looking at the individual RCM projections shown in Figure 6, it is noted that most of those RCMs project a clear EDW signal over the TP for all seasons under the RCP8.5. Despite that, there exists disagreement in both structure and magnitude of EDW among the RCMs due to the RCM physics. For example, though downscaled from the same GCM, MPI_R exhibits a consistent enhanced warming rate up to about 5,000-5,500 m in SON, while the EDW signal in MPI_W can only be found below 4,000 m. Although similar structures could be found between GFDL_R and GFDL_W, different magnitudes are seen between them, NIU ET AL.
10.1029/2020JD033997 8 of 17 particularly in the high elevations during the warm season. Such disagreements appear to be explained by different RCM configurations, such as the land surface model (Minder et al., 2018), microphysics scheme , etc., which affect the simulation of snow cover and specific humidity and ultimately EDW. Also, the structure and magnitude of simulated EDW are sensitive to the driving GCMs, as the spread in EDW among WRF (RegCM4) models remains very high. Climate warming tends to vary with elevation in the troposphere, and such variations influenced by the driving GCMs through the lateral boundary may contribute somewhat to EDW (Minder et al., 2018;Rupp et al., 2017;Walton et al., 2016). Comparatively, the inter-model discrepancies in EDW among RegCM4 models are larger than that in WRF models, indicating greater sensitivity of RegCM4 to the anthropogenic greenhouse forcing.

Mechanisms Influencing EDW
The EDW has been proved over mountain regions of the world (Dimri et al., 2018;Pepin et al., 2015), and is attributed to several factors, such as SAF (Guo et al., 2016;Walton et al., 2016), lapse rates in the troposphere (Minder et al., 2018;Rupp et al., 2017), cloud feedback (Duan & Wu, 2006), longwave water vapor feedback (Rangwala et al., 2013), etc. As a detailed analysis based on regional energy budgets could potentially identify the main processes responsible for EDW (Pepin et al., 2015), Equation 2 is used to quantify the contributions of each energy budget component to EDW over the TP. Considering the pronounced EDW and associated large spread among multi-RCM in SON, we will hereafter discuss in detail for SON only.
As is shown in Figure 7, the projected EDW can be largely explained by the SAF during SON, because of the strikingly similar profiles between the SAF-induced temperature changes and EDW, supporting the dominant role of SAF in shaping EDW over the TP (Rangwala et al., 2013;Palazzi, Mortarini, et al., 2019). Under the background of global warming, there is an altitude dependence of losing snow cover over the TP (Figure 8), which leads to elevation dependent reduction in surface albedo and increase in absorbed solar radiation and ultimately EDW via SAF (Minder et al., 2018). The warming rate due to SAF almost increases with elevation during SON, with the maximum warming about 1-6.4°C around 5,000-5,500 m ( Figure 7b). Besides that, the change in albedo can alter the partitioning of latent heat flux (LH) and sensible heat flux (SH) at the surface. Therefore, it is interesting to find that the temperature change due to SH follows that caused by SAF, but with opposite sign and weaker magnitude (Figures 7b and 7f). The canceling out of those two terms exerts a positive contribution to the EDW over the TP under RCP8.5.
NIU ET AL.
10.1029/2020JD033997 9 of 17 Figure 6. Relationship between the projected change in T as (unit:°C) and elevations over the TP in different regional climate models and their MME. Vertical whiskers show spatial variability of the T as changes within each bin.
The importance of DLW in shaping EDW has been documented by Rangwala et al. (2013) and Palazzi, Filippi, and von Hardenberg (2017). Instead, the contribution of DLW to EDW is under debate based on our simulations, although the correlation between the DLW-induced temperature changes and T as projections is significant (r = 0.954) in MME. For example, there is no clear relationship between the DLW-induced temperature change and elevation in some of the models (such as CNRM_W, EC-EARTH_W, and GFDL_R), while MPI_W and GFDL_W project an EDW due to DLW at certain altitudes ( Figure 7d). The effect of downward shortwave radiation (DSW) influenced by cloud cover is also model dependent since its contribution to EDW is negligible in WRF models and is positive in most RegCM4 models (Figure 7c). For the ground heat flux (GHF) and LH, their changes under RCP8.5 lead to small temperature variations along elevations, with temperature change differing by less than 0.3°C between elevation bins.
NIU ET AL.

10.1029/2020JD033997
10 of 17 Overall, the projected EDW over the TP during SON is primarily caused by the SAF, and its contribution is somewhat canceled out by SH. The effect of DLW is model dependent, as DLW appears to be the second role in shaping EDW in some of RCMs while contributes little in other RCMs. Similar conclusions can be drawn for the other seasons, although the vertical structure and magnitudes of each energy budget are different ( Figures S1-S3). Note that the contribution of DLW to EDW in EC-EARTH_R is larger than that of SAF during DJF, implying the relative role of different driving mechanisms may be dependent on model configurations.

Influencing Factors of Temperature Projection and EDW
To identify the factors that dominate the temperature projections and inter-model differences, the contribution of each energy component to temperature projections and their inter-model standard deviation are presented in Figure 9 and S4. Here, the contribution is calculated as the ratio of temperature change induced by the energy component to the T as change. The contribution of each component varies within the elevation bins and models. The warming pattern over the TP is driven primarily by DLW, with its contribution decreasing with elevation in most of RCMs. The essential role of DLW in the surface warming amplification over the TP during SON has been verified by Gao, Duan, et al. (2019) for the present climate using two reanalysis datasets. The SAF also makes a large contribution to the TP warming under the RCP8.5, and its contribution increases with elevations for most of the models. Contrary to SAF and DLW, changes in LH have a cooling effect on the TP, with smaller contributions over the higher altitudes than the lower altitudes. While for other components, their contributions depend on model and elevation. For example, positive contributions of DSW to T as change are simulated by CNRM_R below 4,500 m, while negative contributions are found in CNRM_W.
Interestingly, the inter-model standard deviation of DSW shows a similar magnitude to that of DLW below 4,500 m ( Figure S4), although the contributions of DSW are far less than that of DLW. DLW together with DSW acts as the main contributors to the inter-model differences below 3,500 m, which indicates that the spread in the cloud and water vapor among RCMs seems to explain the uncertainty in temperature projections over the lower altitudes. Besides, the SAF strongly enhances the inter-model spread of temperature projections under RCP8.5 over the regions above 3,500 m, which is partly counteracted by the SH (Figure S4).
NIU ET AL.

10.1029/2020JD033997
11 of 17 When comparing the EDW projections among the RCM, it is found that the disagreements in EDW during SON are primarily caused by the SAF (Figure 7). The spread in SAF rises with elevation and reaches its maximum up to 8°C around 5,000-5,500 m. This may derive from the uncertainty in projections of snow cover (Figure 8), which is found to be sensitive to the choice of RCM and the driving GCM (McCrary & Mearns, 2019). Moreover, considerable spread in the projected temperature due to DSW and DLW is also found below 4,500 m during SON, and their ranges reach about 1.5°C with small variation along elevation. These imply that the spread in DSW and DLW among the RCMs contributes little to the difference in EDW. Interestingly, the spread in individual partial temperature changes among RegCM4 models is almost larger than that among WRF models.

Effect of Driving GCM and RCM Physics on Temperature Projections
As discussed above, the driving GCM and RCM physics exert a profound influence on the temperature projections over the TP, in terms of magnitude, spatial distribution, and EDW. To understand to what extent the two factors affect the projections, we compare the temperature change and the elevational gradients of T as change between the RCMs and their driving GCMs ( Figure 10) and between the two models driving by the same GCM (Figure 11), respectively. Here, the elevational gradients of T as change are expressed as the slope of the linear regression between T as changes and elevation.
Compared to the driving GCMs, RCMs substantially modulate the projected temperature changes over the TP, with the relationship exhibiting seasonal dependence. The warming rate is always enhanced in JJA and SON but decreased in DJF and MAM after downscaling. The strength of temperature response in RCM does not always follow that in GCM. For example, EC-EARTH project a larger warming rate than GFDL in MAM while the RCMs driven by EC-EARTH show smaller temperature change than those driven by GFDL. Note that the spread in the seasonal temperature projections is much larger in the RCMs than the GCMs except MAM. This implies that the temperature projections over the TP may be more sensitive to the RCM physics and less influenced by the driving GCM. McCrary and Mearns (2019)   RCP8.5 (Figure 10b), indicating the essential role of the driving GCM in projected EDW. However, the nested RCMs always amplify the EDW signals in GCMs or even change their signs. Additionally, the spread in slopes is larger among RCMs than GCMs. Those results suggest the importance of RCM physics in the EDW.
The RegCM4 ensemble average (MMA_R) generates a greater temperature response than the WRF ensemble average (MMA_W) for all seasons, particularly in SON (Figure 11a), which indicates stronger sensitivity of RegCM4 to the anthropogenic greenhouse forcing than that of WRF. As described above and shown in Figure S5, there is little agreement on the vertical profile of T as changes between RegCM4 and WRF. For example, the projected EDW in DJF occurs below 4,000 m in MMA_W while it can reach 5,000 m in MMA_R. The spatial variability of T as change is always larger in MMA_R compared to that in MMA_W, especially during JJA and SON. Similar to the T as projection, MMA_R also presents stronger EDW signals than MMA_W (Figure 11b). Driven by MPI, the elevational gradients of warming rates in SON could reach 0.7°C/ km in RegCM4 while WRF generates negative elevational gradients of warming. Such notable differences in EDW between RegCM4 and WRF may be related to their different representation of land surface processes, which affect the simulated snow cover and then modulate the simulated SAF and its effect on EDW (Minder et al., 2018). All those results demonstrate that the projected EDW is highly sensitive to the choice of RCM.
NIU ET AL.

10.1029/2020JD033997
13 of 17 Figure 10. Comparison of projected temperature changes (a, unit:°C) and elevational gradients of projected temperature change (b, unit:°C/km) in the driving global climate models and nested regional climate models. Figure 11. Comparison of projected temperature changes (a, unit:°C) and elevational gradients of projected temperature change (b, unit:°C/km) in two regional climate models and their ensemble average (MMA). The number stands for the driving global climate model.

Conclusions and Discussion
Under global warming, warming rates in mountain regions often depend on elevation (Hock et al., 2020).
Observations show that the dependence can be either positive or negative and the relationship between the warming rate and elevations may not be linear (Palazzi, Mortarini, et al., 2019;You et al., 2020). There is growing evidence of EDW (positive elevational gradients of warming) over the TP in recent decades (Guo, Sun, et al., 2019;Liu et al., 2009), and whether this feature will continue to exist in the future gains wide attention. This study focuses on the character and causes of future EDW over the TP under the RCP8.5 scenario using two RCMs driven by four GCMs, and highlights the effect of driving GCM and RCM physics on those projections.
For the historical climate, large cold biases over the TP exist in all RCMs and their ensemble mean (MME), although they could capture the spatial distribution of observed surface air temperature (T as ). Such cold biases may be attributed to the driving GCMs (Gao et al., 2016), deficiencies in the model (Meng et al., 2018), and observational uncertainties . Additionally, the location of stations lying in valleys and lack of high-elevation stations may also contribute to the model cold bias . Comparatively, the internal model physics has more influence on the T as distribution than the driving GCMs during MAM and DJF, while both of them play an important role during other seasons in our experiments.
Widespread warming over the TP is projected by MME under RCP 8.5 for 2031-2060, with the largest change of 2.4°C in SON and the lowest of 2.0°C in MAM. Large disagreements in the warming intensity and maximum warming center exist among the RCMs' projections, and those projections are more consistent with each other in cold seasons (MAM and DJF) than in warm seasons (JJA and SON), which is consistent with findings of Wu and Gao (2020). The largest spread in temperature projections is found above 5,000 m, which can to a large extent be explained by the difference in surface albedo feedback (SAF) among the RCMs. Further investigations employing a larger number of multimodel simulations should be carried out to obtain reliable climate projections.
Consistent with previous studies using the CMIP5 model ensemble (Palazzi et al., 2017(Palazzi et al., , 2019Rangwala et al., 2013), a clear EDW signal over the TP is projected under the RCP8.5 by MME, with its magnitude and structure varying with seasons. Larger spatial variability of T as projections is found at higher elevations in MME during most seasons, suggesting different temperature responses to anthropogenic greenhouse forcing at higher elevations. Furthermore, the structure and magnitude of projected EDW are sensitive to the RCM physics and driving GCM, implying that large uncertainty exists in the projected EDW over the TP.
The projected EDW in SON seems to be primarily caused by the changes in SAF, because of the strong correspondence profiles between the SAF-induced temperature changes and EDW. The dominant role of SAF has been confirmed by sensitivity experiments in which EDW is eliminated without an active SAF (Minder et al., 2018). The SAF also acts as the main source of uncertainty in EDW projections among RCMs. The different choice of RCM or driving GCM results in diverse projections of snow cover and albedo, which can modulate the simulated SAF and its effect on EDW. The downward longwave radiation (DLW) is found to be the dominant factor in regulating T as change over the TP, and its contribution decreases with elevation. However, the relationship between the DLW and EDW is dependent on the model used, which is different from previous studies that emphasize the significant correlation between DLW and EDW (Palazzi, Filippi, & von Hardenberg, 2017;Rangwala et al., 2016). Additionally, the spread in DLW is large with small vertical variations, implying that DLW contributes little to the different EDW among RCMs.
When assessing the influence of driving GCMs and RCM physics on the temperature projections over the TP, it is found that the spread in the seasonal temperature projections is always larger in the RCMs than the GCMs, implying the temperature projections over the TP may be more sensitive to the RCM physics and less influenced by the driving GCM. A similar conclusion for the snow projections is also drawn by McCrary and Mearns (2019) based on the RCM-GCM ensemble. Our study shows that the EDW signal in the nested RCMs mostly follows that in the driving GCM, indicating the essential role of the driving GCM in projected EDW. However, the amplified EDW signal in RCMs and the disagreement in magnitude and sign of EDW among RCMs indicate the projected EDW is sensitive to RCM physics. Additionally, compared to WRF, RegCM4 shows higher sensitivity to the anthropogenic greenhouse forcing, evidenced by the larger temperature projections and stronger EDW signal in RegCM4.
Our study confirms the EDW signal for T as is expected to continue in the future over the TP, which will exert a profound influence on the high-altitude cryosphere system and mountain ecosystems. Another related question is whether the EDW will occur for the daily maximum (T max ) or minimum (T min ) temperature. As suggested in Figures S6 and S7, MME projects EDW in both T max and T min for all seasons, particularly strong in JJA for T max and in SON for T min . Different patterns between the T max and T min are found in DJF when the warming rate increases with elevation for T max while reaches its maximum at middle elevations for T min . Such robust EDW over the TP in the future is consistent with the previous projections from GCMs (Palazzi, Mortarini, et al., 2019;Rangwala et al., 2016;You et al., 2019) and RCMs (Gao, Chen, Lettenmaier, et al., 2018;Guo, Sun, et al., 2019). Compared to T max , T min exhibits a stronger tendency toward EDW, and shows a larger spread in simulated EDW among the RCMs ( Figure S8), indicating that separate mechanisms may be at work during the day and night. However, the temporal resolution (six-hour) in our study is not high enough to identify the main processes leading to EDW for T max and T min . Pepin et al. (2015) find that the specific temperature response depends on soil moisture, as the response will be more intense in T max (rather than T min ) when the increased surface shortwave absorption is balanced by increases in SH (rather than LH). Also, the SAF is expected to have a more dominating influence on T max than T min (Liu et al., 2009;Rangwala et al., 2010), and the near surface specific humidity plays the primary role in shaping EDW for T min (Rangwala et al., 2016). Note that the signs of the slope for T min in RCMs are more consistent with those in driving GCMs compared to T max ( Figure S8).
Considering the T as projection and EDW are substantially moderated by the RCM physics and driving GCMs, our conclusions apply specifically to the model settings considered in this study. Note that the surface energy budget approach used to diagnose the energetic contributors to temperature change is based on the prerequisite that these contributors are independent of each other. However, this may not true. For example, previous studies often consider DLW as an independent forcing on the surface energy budget that drives surface warming (Burt et al., 2016;Gao, Duan, et al., 2019). Zeppetello et al. (2019) stress that DLW is tightly coupled to surface temperature, and it cannot be regarded as an independent component. Therefore, new approaches are still required to improve such analyses. Additionally, Pepin et al. (2015) emphasize that models can also be used to perform sensitivity experiments to quantify the role of specific variables to EDW, and the climate response to greenhouse-gas forcing in complex terrain can only be adequately resolved by models with a resolution of 5 km or less. Now, a variety of experiments using the convection-permitting models have been (Gao, Chen, & Jiang, 2020;Li, Furtado, et al., 2021;Zhou et al., 2021) and will be performed under the new framework of CORDEX-FPS-CPTP (Convection-Permitting Third Pole). Those results could strengthen recent findings from the CORDEX-EA project and provide more robust climate projections over the TP as well as the EDW.

Data Availability Statement
All data used in this study are freely accessible. The satellite-based surface air temperature datasets are available from https://zenodo.org/record/1250228#.