Evaluation of the surface air temperature over the Tibetan Plateau among different reanalysis datasets

The surface air temperature (SAT) over the Tibetan Plateau (TP) not only affects the physical processes such as local evaporation, snow melting, and glacier ablation, but also has a great impact on the downstream regions and even the global climate change. The development of reanalysis data has gradually overcome the problem of sparse stations over the TP, but there are still some deficiencies. Therefore, the distance between indices of simulation and observation (DISO) method is used to calculate the distance between five reanalysis datasets (ERA5, JRA-55, ERA-Interim, MERRA2, NCEP2) and the CMFD to evaluate the abilities of different reanalysis datasets to capture the SAT over the TP in different seasons. The results indicate that ERA-Interim has a higher ability to reproduce the SAT over the TP in spring and summer, while it is ERA5 in autumn and winter. It should be noted that although the optimal reanalysis has a better performance in capturing the SAT of the TP, there are still a certain degree of deviations in their spatial fields. We further show the spatial deviation fields of SAT over the TP corresponding to the optimal reanalysis data in different seasons, and analyze the possible reasons. The result implies that the SAT deviation field is mainly related to the snow in winter and spring, while the summer SAT deviation field is mainly related to the water vapor, and the autumn is related to both the snow and the water vapor fields. Overall, the quality of reanalysis data needs to be further improved in the future.


Introduction
The Tibetan Plateau (TP) with an average elevation of more than 4,000 m is close to the middle of the troposphere (Qiu, 2008;Yao et al., 2012;Kuang and Jiao, 2016), which locates the first level of terrain in China and also called as "Asian Water Tower" (Xu et al., 2008;Xu et al., 2019;Xue et al., 2021). Covering an area of about 2.5 million square kilometers, the TP is crisscrossed by mountains and dotted with lakes and known as a typical region of the multi-sphere interactions in the Earth's climate system (Lei et al., 2013;Lei et al., 2017;Yao et al., 2015;. Due to its special geographical location and topographic distribution, small environmental disturbances and changes can result in the "butterfly effect" of the climate system and ecological environment in this region, thus the TP is also the sensitive and the amplifying area of the global climate change You et al., 2020;Tang et al., 2022;Wang et al., 2022).
In the context of global warming, the surface air temperature (SAT) over the TP has increased substantially (Duan and Xiao, 2015;Duan et al., 2020;You et al., 2021;Zhou and Zhang, 2021), and the resultant thermal forcing changes have a significant impact on regional and global climate. For example, the increase of the spring SAT over the TP will influence the surface sensible heat flux (SSHX), and the change of the SSHX will further affect the downstream East Asian summer monsoon precipitation and North Pacific sea surface temperature (Duan et al., 2011;Wang et al., 2014;Sun et al., 2019;Liu et al., 2020); meanwhile, the warming of the TP can enhance the intensity and frequency of the warmest day and night over it (You et al., 2008;; on the other hand, drastic changes have also occurred in its cryosphere, such as the accelerated melting of snow and glaciers, which reduces the surface albedo and forms a positive feedback to promote the TP warming Li et al., 2017;Dai et al., 2017;Xu et al., 2017;Zhang et al., 2021). All these changes may accelerate the pace of the TP to become a global tipping point (Liu et al., 2023) However, due to the complex terrain characteristic of the TP, most conventional meteorological observation stations are located in the valleys and relatively low-altitude areas in the central and eastern plateau, and there are few stations in its western region (Duan et al., 2018;Duan et al., 2022a), where there exist a large number of glaciers and permafrost (Yang et al., 2019;. Therefore, the sparse and uneven distribution of observational stations is the primary difficulty in accurately understanding the climate change on the TP.
The gridded reanalysis datasets are generated by assimilating the observational data obtained by ground stations, satellite remote sensing, aircraft and ships into the numerical model, which has the advantages of long time and wide coverage area (Dee et al., 2011;Rienecker et al., 2011;Kobayashi et al., 2015;Hersbach et al., 2020), providing us important database for understanding atmospheric circulation and global and regional climate change (Simmons et al., 2004;Yi et al., 2011;Lei et al., 2017;Tang et al., 2022). The development of reanalysis data overcomes the problem of sparse stations and has been widely used in the study of climate change on the TP Sheng et al., 2021;. Nonetheless, due to differences of the observational data, numericalmodels and assimilation methods used in the different reanalysis datasets, the results obtained from different reanalysis datasets may be different particularly in mountain areas Peng et al., 2021;Yang et al., 2021). Therefore, the purpose of this paper is to use the DISO (distance between indices of simulation and observation) method (Hu et al., 2019; to evaluate the ability of capturing the SAT over the TP among different reanalysis datasets in different seasons, so as to lay a foundation for better research on the temperature change of the TP.

Research data
The observational data used in this paper are monthly average surface air temperature named as China Meteorological Forcing Dataset (CMFD) and daily average snow depth in China provided by A Big Earth Data Platform for Three Poles. The spatial resolution of CMFD data is 0.1°× 0.1°and the period covers from 1979 to 2018. IGURE 1 (A) Trends (°C/decade), (B) climatology (°C), and (C) interannual variability (°C) of surface air temperature (SAT) over the Tibetan Plateau (TP) among the different datasets for four seasons. The color of red, light blue, blue, green, orange, and yellow represents the observational data CMFD, ERA5, JRA-55, ERA-Interim, MERRA2, and NCEP2 data, respectively.

Frontiers in Environmental Science frontiersin.org
This dataset is merged into remote sensing products, reanalysis datasets, and meteorological station data , and has a good performance over the TP (Lun et al., 2021;Peng et al., 2022), which can be regarded as observational data. In addition, the daily average snow depth in China, with the spatial resolution of 0.25°× 0.25°and a time range of 1979-2019 (Che et al., 2008;Dai et al., 2015), is widely used over the TP (Duan and Xiao, 2015;Bao and You, 2019). The domain of the TP is defined as the region above 2500 m within the range of 25°-40°N and 70°-105°E. The atmospheric reanalysis datasets include ERA-Interim (Dee et al., 2011) andERA5 (Hersbach et al., 2020) provided by the European Center for Medium-Range Weather Forecast (ECWMF). The ERA5 dataset is the fifth-generation reanalysis product developed by ECMWF. Compared with ERA-Interim data, the ERA5 data has the advantages of higher spatial-temporal resolution and longer covered time. The variables including the monthly average surface air temperature, snow depth and threedimensional specific humidity from 1979-2019 provided by ERA-Interim and ERA5, with the spatial resolution is 1°× 1°and 0.25°× 0.25°, respectively. Furthermore, we also use the monthly mean surface air temperature datasets from JRA-55 (Kobayashi et al., 2015), MERRA2 (Rienecker et al., 2011;Gelaro et al., 2017) and NCEP2 (Kalnay et al., 1996), which are provided by the Japan Meteorological Agency (JMA), National Aeronautics and Space Administration (NASA) and National Centers for Environmental Prediction (NCEP), respectively. The coverage time of the JRA-55 and NCEP2 datasets are 1979-2019, while MERRA2 is 1980-2019 and the spatial resolutions are 1.25°× 1.25°, 2.5°× 2.5°, and 0.5°× 0.625°, respectively. The detailed information about the reanalysis datasets is listed in Table 1.
Considering the consistency of data coverage time, the comparison period between the SAT of the reanalysis datasets and the CMFD is 1980-2018. All the observational and reanalysis data are uniformly interpolated to the spatial resolution of 0.5°× 0.5°.

Methods
The DISO (distance between indices of simulation and observation) method is a means of evaluating data quality. Its essence is to calculate the distance between the simulated and observed data. The smaller the distance, the higher the quality of simulated data (Hu et al., 2019;Hu et al., 2022;. The advantage of this method is that it is simple, convenient and flexible, reflecting the idea of Da Dao Zhi Jian . Researchers can choose the corresponding evaluation index according to their own needs, and the calculation of DISO value can be extended to multi-dimensional space. Therefore, the method is able to evaluate the data comprehensively and can be applied to many disciplines, such as meteorological and disease studies Hu et al., 2020;Deng et al., 2021;Wang et al., 2021;Xu et al., 2022;Yin et al., 2022).
In this paper, three metrics (trend, climatology, and interannual variability) are selected to calculate the deviations between the reanalysis datasets and the CMFD by the DISO method. Given

FIGURE 2
The ability of different reanalysis datasets to capture the SAT over the TP in winter (A), spring (B), summer (C) and autumn (D). The smaller the DISO value, the better quality of the data.

Frontiers in Environmental Science
frontiersin.org

FIGURE 3
The spring SAT spatial distribution in CMFD and ERA-Interim. (A-C) represent the CMFD temperature trend field, the ERA-Interim temperature trend field, and the difference field (°C/decade) between the ERA-Interim and CMFD, respectively; (D-F) are similar to (A-C) but represent the climatology fields (°C); (G-I) are similar to (A-C) but represent interannual variability fields (°C). The TP domain is outlined by the black curve with height above 2,500 m.

FIGURE 4
Spatial distribution of the spring snow depth over the TP in observational and ERA-Interim data. (A) the observational trend (cm/decade) of spring snow depth; (B) is similar to (A) but for the ERA-Interim data (cm/decade); (C) the difference of snow depth trend field (cm/decade) between ERA-Interim and observation data; (D-F) are similar to (A-C) but represent the climatology fields (cm). The TP domain is outlined by the black curve with height above 2,500 m.
Frontiers in Environmental Science frontiersin.org 04

FIGURE 5
The summer SAT spatial distribution in CMFD and ERA-Interim. (A-C) represent the CMFD temperature trend field, the ERA-Interim temperature trend field, and the difference field (°C/decade) between the ERA-Interim and CMFD, respectively; (D-F) are similar to (A-C) but represent the climatology fields (°C); (G-I) are similar to (A-C) but represent interannual variability fields (°C). The TP domain is outlined by the black curve with height above 2,500 m.

FIGURE 6
The summer water vapor fields (g/kg/decade) integrated from 600 hPa to 200 hPa over the TP in ERA-Interim (A) and ERA5 (B), respectively; (C,D) are similar to (A,B) but represent the autumn water vapor fields. The TP domain is outlined by the black curve with height above 2,500 m. Stippled areas indicate regions exceeding the 90% statistical confidence level.

Frontiers in Environmental Science
frontiersin.org 05 that the different magnitudes of the deviations, the deviation of each metric should be normalized before calculating the DISO value. The equation for calculating the DISO value is as follows: where |Re − OBS| represents the absolute value of the deviation between the reanalysis data and observation. The smaller the DISO value, the smaller the difference between the reanalysis data and the observation, and the higher the data quality.
Taking the year 1980 as an example, the spring mean value is defined as the average of March, April, and May, the summer mean value is defined as the average of June, July, and August, the autumn mean value is defined as the average of September, October, and November and the winter mean value is defined as the average of December 1980, January and February 1981. Figure 1A implies that the TP has the most intense warming of 0.7°C/decade in winter, followed by 0.55°C/decade in autumn, and about 0.3°C/decade in spring and summer. Although the reanalysis datasets capture the characteristics of severe warming over the TP in winter, the amplitude of the warming trend is underestimated relative to the observation. The winter warming trend of 0.5°C per decade in ERA5 is the closest to the observation, and the deviation is the largest for the MERRA2 data, which is only 0.27°C/decade. In autumn, the warming trends displayed by the five sets of reanalysis datasets are still weak. The result of ERA5 is the closest to the CMFD, while the NCEP2 has the strongest bias. In spring and summer, except for NCEP2, the deviation of the warming trends presented by other reanalysis data have weakened compared with the counterpart in autumn and winter, and are closer to the observation ( Figure 1A). In terms of climatology ( Figure 1B), the CMFD shows that the winter temperature over the TP is about −13°C, and the ERA5 and NCEP2 have cold biases, while other reanalysis data have warm biases. Compared with the results of the trend, the reanalysis datasets present a smaller deviation of the winter temperature climatology. In summer, except for JRA-55, the climatology presented by the other reanalysis data are almost consistent with the observed values. In spring and autumn, the ERA5 and NCEP2 have weaker cold biases, and the ERA-Interim and MERRA2 have weaker warm biases, while the JRA-55 has stronger warm biases. For the interannual variability, the CMFD data displays the strongest interannual variability in winter, and all the reanalysis datasets reproduce this feature, among which the

FIGURE 7
The autumn SAT spatial distribution in CMFD and ERA5. (A-C) represent the CMFD temperature trend field, the ERA5 temperature trend field, and the difference field (°C/decade) between the ERA5 and CMFD, respectively; (D-F) are similar to (A-C) but represent the climatology fields (°C); (G-I) are similar to (A-C) but represent interannual variability fields (°C). The TP domain is outlined by the black curve with height above 2,500 m.

Frontiers in Environmental Science
frontiersin.org ERA5 and NCEP2 are the closest to observation, while the MERRA2 is the worst. In other seasons, the results of the ERA5 exhibit stronger interannual variability relative to the observations, and the JRA-55, MERRA2, and ERA-Interim exhibit weaker interannual variability, while the interannual variability given by the NCEP data is stronger in spring and weaker in summer and autumn ( Figure 1C). Clearly, the five sets of reanalysis data have different ability in reproducing the trend, climatology and interannual variability of SAT over the TP in different seasons. Therefore, we select three factors to comprehensively evaluate the ability of the reanalysis datasets to capture the SAT over the TP in different seasons by the DISO method. According to the DISO values, the ability of the reanalysis data to reproduce the winter temperature over the TP from strong to weak is: ERA5, NCEP2, JRA-55, ERA-Interim, MERRA2 (Figure 2A). The ranking is ERA-Interim, MERRA2, ERA5, JRA-55, NCEP2 in spring ( Figure 2B), while it is ERA-Interim, MERRA2, ERA5, NCEP2, JRA-55 in summer ( Figure 2C) and ERA5, ERA-Interim, NCEP2, MERRA2, JRA-55 in autumn ( Figure 2D).
Generally, the ERA-Interim has the highest quality among these reanalysis data in spring and summer, while it is ERA5 in autumn and winter. It is worth noting that although the capturing ability of the ERA-Interim in spring and summer and the ERA5 in autumn and winter is better than other reanalysis data, the DISO values are about 0.3 and 0.8 respectively, indicating that there still exist some deviations with the observation even for the optimal reanalysis data. Figure 3 shows the spatial distribution of spring surface air temperature over the TP from the CMFD and ERA-Interim. The ERA-Interim shows that the spring warming is weaker over the most parts of the TP ( Figure 3B) relative to the observation ( Figure 3A), and the bias field also reflects this feature ( Figure 3C). For the climatology, the SAT fields of the ERA-Interim and the observation are generally consistent ( Figure 3D; Figure 3E), while the deviation field implies that the ERA-Interim has a warm bias in most areas of the TP ( Figure 3F). Figure 3I demonstrates that the interannual variability given by ERA-Interim on the TP is weaker compared with the observation ( Figure 3G). Considering the ERA5 is the improved version of the ERA-Interim, we also further analyzed the SAT bias field of the ERA5 dataset. Supplementary Figure S1C reveals that the spring SAT trend deviation field in ERA5 is similar to that in Figure 3C, and there is a cold bias over the most parts of the TP. The difference is that the ERA5 data has a cold bias for climatology (Supplementary Figure S1F) and the strong interannual variability (Supplementary Figure S1I) over the most areas of the TP. The previous study pointed out that the temperature change over the TP in spring was mainly related to the surface albedo (Gao et al., 2019). Figure 4 and Supplementary Figure S2 present the variation of snow depth for the observations, ERA-Interim and ERA5. It can be found that compared with observations, the snow depth decline trends captured by the ERAinterim and ERA5 are weaker ( Figure 4C and Supplementary Figure  S2C). However, the ERA-Interim has less snow ( Figure 4F) and the ERA5 has more snow (Supplementary Figure S2F) than the observation ( Figure 4D) in terms of climatology. Therefore, the deviation of the spring snow depth could explain the deviation field of the spring SAT to a certain extent.

FIGURE 8
Spatial trends and time series of autumn snow depth over the TP from the observational and ERA5 data. (A) The observational trend (cm/decade) of autumn snow depth over the TP; (B) is similar to (A) but for the ERA5 data (cm/decade); (C) difference in snow depth trend field (cm/decade) between ERA5 and observation data; (D) time series of autumn snow depth changes (cm) in the observation and ERA5 data. The TP domain is outlined by the black curve with height above 2,500 m.

Frontiers in Environmental Science
frontiersin.org

FIGURE 9
The winter SAT spatial distribution in CMFD and ERA5. (A-C) represent the CMFD temperature trend field, the ERA5 temperature trend field, and the difference field (°C/decade) between the ERA5 and CMFD, respectively; (D-F) are similar to (A-C) but represent the climatology fields (°C); (G-I) are similar to (A-C) but represent interannual variability fields (°C). The TP domain is outlined by the black curve with height above 2,500 m.

FIGURE 10
Spatial trends and time series of winter snow depth over the TP in observational and ERA5 data. (A) The observational trend (cm/decade) of winter snow depth over the TP; (B) is similar to (A) but for the ERA5 data (cm/decade); (C) difference in snow depth trend field (cm/decade) between ERA5 and observation data; (D) time series of winter snow depth changes (cm) in the observation and ERA5 data. The TP domain is outlined by the black curve with height above 2,500 m.
Frontiers in Environmental Science frontiersin.org 08 Figure 5 is the spatial distribution of summer SAT over the TP from the observations and ERA-Interim. The ERA-Interim shows a weaker warming over the TP relative to the CMFD in summer, and the deviation field also reflects this feature (Figures 5A-C). For the climatology of summer SAT, the ERA-Interim and observation are relatively consistent, while the bias field indicates that ERA-Interim has a warm bias over the eastern region of the TP and a cold bias over the western region of the TP (Figures 5D-F). Moreover, Figure 5I demonstrates that the summer SAT presented by the ERA-Interim displays weaker interannual variability over the TP. Supplementary Figure S3 implies that the results of ERA5 are basically consistent with the ERA-Interim, but the summer DISO value of ERA5 is about 0.7. This is because the bias of the ERA5 summer SAT is stronger than that of the ERA-interim in terms of the trend and interannual variability. The previous studies concluded that the SAT changes are mainly related to water vapor over the TP in summer (Gao et al., 2019;. Therefore, we calculate the ERA-Interim and ERA5 summer water vapor trend fields. Both the ERA-Interim and ERA5 show a significant increase of water vapor in summer, which corresponds to the warming of the TP (Figures 6A,B). This result also reflects that the temperature deviation over the TP in summer may be mainly related to water vapor that induces atmospheric downward long-wave radiation.
In autumn, we present the spatial distribution of SAT over the TP from the observations and ERA5. Figure 7B shows that the plateau warming in autumn given by ERA5 is weaker than the observation ( Figure 7A), and the deviation field also displays this feature ( Figure 7C). As for the climatology, the spatial distribution of ERA5 and observation is relatively consistent (Figures 7D,E), while the bias field reveals that the ERA5 data overall presents a cold bias ( Figure 7F). Moreover, the interannual variability of ERA5 over the TP is stronger than the observation for the autumn SAT ( Figures  7G,I). Previous study indicated that the temperature changes over the TP are mainly related to water vapor and snow in autumn (Gao et al., 2019). Figure 6D reveals that the water vapor of ERA5 shows a significant increase in autumn, which indicates that it has a certain contribution to the autumn warming of the TP. However, Figures 8A-D imply that the ERA5 reproduces a weaker snow decline trend and has more snow than the observation. The above analysis demonstrates that the autumn SAT deviation field in ERA5 may be mainly related to its snow deviation and water vapor field.
Although the ERA5 presents the characteristics of the winter warming over the TP, the overall trend is weak (Figures 9A-C). As for the climatology of winter SAT, the spatial distribution of the ERA5 and observation are relatively consistent ( Figures 9D,E), while the deviation field shows that ERA5 has a cold bias in the southeast of the plateau and a warm bias in the northwest region ( Figure 9F), and the overall is the cold bias ( Figure 2B). The interannual variability captured by the ERA5 is stronger in the southeastern plateau and weaker in the northwest region relative to the observation ( Figures 9G-I). The previous study demonstrated that the winter temperature change over the TP is mainly related to the surface albedo (Gao et al., 2019). Figures 10A-D reveal that the ERA5 reproduces a weaker snow trend and exhibits more snow than the observation. Hence the bias of winter snow in ERA5 data contributes to its temperature bias field to a certain extent.

Conclusion and discussion
This paper selects three indicators including the trend, climatology, and interannual variability to comprehensively assess the ability of the five sets of reanalysis data to capture the SAT over the TP in different seasons by the DISO method. Results indicate that in spring and summer, the DISO values of ERA-interim are the smallest, which are 0.3 and 0.2, respectively. However, the DISO values of ERA5 are the smallest in autumn and winter, which are 0.8 and 0.7, respectively. Although the ERA-Interim and ERA5 are better than other reanalysis data in the corresponding seasons, there are still some deviations from the observed data. Therefore, we further analyze the spatial differences in the SAT field between the Frontiers in Environmental Science frontiersin.org optimal and observed data in different seasons. The results show that the spring warming trend of the TP captured by the ERA-Interim is weaker than the CMFD, and there is a warm deviation in most regions in terms of climatology, which may be mainly related to the snow deviation existed in the ERA-Interim. In summer, the warming trend presented by the ERA-Interim over the TP is weak, and there is a warm deviation in the eastern part of the plateau and a cold deviation in the western part of the plateau for the climatology, which may be mainly related to the water vapor field in this data. Compared with the CMFD, both the warming trend and climatology captured by the ERA5 in autumn are cold deviations, which may be mainly related to the water vapor field and the larger snow bias existed in the data. The winter warming trend displayed by the ERA5 is weaker than the CMFD, while the climatological bias field indicates that the ERA5 has a cold bias in the southeast of the plateau and a warm bias in the northwest region. The results show that the deviation of winter snow depth has a certain contribution to the deviation of the surface air temperature in ERA5.
In addition, the deviations of the SAT exist in the different reanalysis datasets over the TP are essentially due to the differences in the numerical models used by the development institution, assimilated observational data, assimilation methods, and parameterization schemes. Therefore, it is necessary to conduct a more detailed analysis for the deviation reasons existing in these reanalysis datasets to promote the improvement of the data in the future.
Besides the drastic warming and the strongest interannual variability over the TP in winter as shown in Figure 1, we also find that the winter SAT of the TP has an interdecadal transition from a negative phase to a positive phase in 1999 ( Figures 11A-C). The further analysis demonstrates that the winter sea ice concentration (SIC) over the Barents-Kara Sea also has an interdecadal transition from a positive phase to a negative phase in 1999 ( Figures 11B,C). The recent studies have concluded that the sea ice loss over the Barents-Kara Sea can aggravate the warming of the TP in winter (Duan et al., 2022b) and explaining 61.4% of the winter cooling over China in the past 2 decades . Furthermore, it shows that the interdecadal variability of the winter SIC over the Barents-Kara Sea is an important factor affecting the interdecadal transition of the cold surge path over East Asia . The previous studies also reveal that the interannual variability of the winter SIC over the Barents-Kara Sea has a great impact on the Asian climate system (Wu et al., 2016;Mori et al., 2019;Li et al., 2020;Sun et al., 2022). Therefore, will the winter Barents-Kara SIC loss also have an impact on the SAT of the TP on the interdecadal and interannual time scale? Further exploration of this interesting scientific question is needed in the future.

Data availability statement
The datasets presented in this study can be found in online repositories. The surface air temperature observational data of the China Meteorological Forcing Dataset (CMFD) on 0.1°× 0.1°spatial resolution used in this study are openly available at http://poles.tpdc. ac.cn/en/data/8028b944-daaa-4511-8769-965612652c49/.

Author contributions
YP wrote the initial paper and carried out most of the data analysis. AD revised the paper. CZ helped in the analysis of the data. BT and XZ checked the paper and proposed amendments. All authors contributed to the paper and approved the submitted version.

Conflict of interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher's note
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.