Evaluation of CMIP6 Models and Multi-Model Ensemble for Extreme Precipitation over Arid Central Asia

<jats:p>Simulated historical extreme precipitation is evaluated for Coupled Model Intercomparison Project Phase 6 (CMIP6) models using precipitation indices defined by the Expert Team on Climate Change Detection and Indices (ETCCDI). The indices of 33 Global Circulation Models (GCMs) are evaluated against corresponding indices with observations from the Global Climate Center Precipitation Dataset (GPCC V2020) over five sub-regions across Arid Central Asia (ACA), using the Taylor diagram, interannual variability skill score (IVS) and comprehensive rating index (MR). Moreover, we compare four multi-model ensemble approaches: arithmetic average multi-model ensemble (AMME), median multi-model ensemble (MME), pattern performance-based multi-model ensemble (MM-PERF) and independence weighted mean (IWM). The results show that CMIP6 models have a certain ability to simulate the spatial distribution of extreme precipitation in ACA and the best ability to simulate simple daily intensity (SDII), but it is difficult to capture the spatial bias of consecutive wet days (CWD). Almost all models represent different degrees of wet bias in the southern Xinjiang (SX). Most GCMs are generally able to capture extreme precipitation trends, but to reproduce the performance of interannual variability for heavy precipitation days (R10mm), SDII and CWD need to be improved. The four multi-model ensemble methods can reduce the internal system bias and variability within individual models and outperform individual models in capturing the spatial and temporal variability of extreme precipitation. However, significant uncertainties remain in the simulation of extreme precipitation indices in SX and Tianshan Mountain (TM). Comparatively, IWM simulations of extreme precipitation in the ACA and its sub-regions are more reliable. The results of this study can provide a reference for the application of GCMs in ACA and sub-regions and can also reduce the uncertainty and increase the reliability of future climate change projections through the optimal multi-model ensemble method.</jats:p>

section_content: According to the World Meteorological Organization report, the global average temperature is 1. 11 ± 0. 13 • C warmer than the pre-industrial level [1]. The most recent seven years, 2015-2021, were the seven warmest years on record. With global warming, the water cycle accelerates, causing an uneven spatial and temporal distribution of global water resources [2]. The frequency, intensity and duration of extreme climate events (e. g., extreme precipitation, drought, hurricanes, etc. ) triggered by global warming have significantly increased, with significant impacts on society and ecosystems [3] [4] [5] [6]. Therefore, it is an important reference value and of practical significance to evaluate the simulation ability of extreme precipitation for the scientific projection of future climate change and the formulation of disaster prevention and mitigation policies. 

Recently, the International Coupling Model Comparison Project has been developed into the sixth phase (CMIP6, see Table 1 for nomenclature, similarly hereinafter). Compared with CMIP5, CMIP6 has noticeably improved in terms of spatial resolution, physical parameterization and composition [7] [8] [9]. However, the different sensitivities and dynamical frameworks of different climate models respond to climate change, resulting in widely varying simulation results of seasonal variability, interannual variability and trend changes of regional precipitation [10] [11] [12] [13] [14] To effectively improve the simulation capability of GCMs, a lot of studies have been conducted by domestic and foreign scholars, such as deviation correction [15] and downscaling [16, 17]. In addition, studies show that the multi-model ensemble can reduce the uncertainty in simulation, and the results are better than those of individual models [18] [19] [20]. The most widely used method is the unweighted scheme, which is simple to calculate, but it ignores the differences between models, and the simulation effect is not as good as that of the weighted scheme [10, 21]. Bishop and Abramowitz [22] developed a new independence weighted mean method (IWM), which has achieved good simulation results in several studies. For example, Bai et al. [23] reported that IWM monitors historical trends in extreme climate indices better than individual models and MME, with lower root-mean-square errors. Zarrin et al. [24] also noted that IWM can significantly reduce uncertainty and outperformed other models in terms of seasonal variations. 

Arid Central Asia is extremely sensitive to global climate change due to its special geographical location as well as climatic conditions [25] [26] [27]. Over the past 60 years, the average annual temperature and extreme temperature in Central Asia (CA) have shown a significant warming trend, with a much higher rate of increase than that of the global and northern hemisphere [28]. Extreme precipitation has also shown an increasing trend with significant spatial heterogeneity [15, [29] [30] [31]. Many scholars attribute this climate change to the internal variability of the climate system, such as the Southern Oscillation, the North Atlantic Oscillation and the Asian subtropical jet [32] [33] [34]. 

In recent years, extensive studies have been conducted on evaluating and projecting extreme climate change over ACA. The prediction studies of temperature and temperature converge, with a significant warming trend expected for the 21st century [20, [35] [36] [37], while the simulation and prediction of extreme precipitation are more uncertain [9, 38]. Li et al. [27] showed that almost all models (HighResMIP) overestimated the precipitation frequency in CA, with a greater overestimation of the precipitation frequency and amount in mountainous areas. Zhang and Wang [39] reported that extreme precipitation in CA is expected to increase approximately linearly with increasing global warming. It has also been projected that seasonal extreme precipitation in southern CA follows an opposite trend than that in northern CA [40]. 

Nevertheless, it can be found that most of the existing studies analyze the simulation ability of the study region as a whole and fail to clarify the influence of topography on the model simulation ability. Since the topography of Central Asia is a complex and significant spatial heterogeneity of climate change, it is actually unscientific to ignore the regional differences and discuss the applicability of the models. In addition, the historical simulation and projection of climate change in ACA are still mostly based on the unweighted scheme, and the research on the weighted scheme still needs to be improved. In view of the above, this paper divided the study area according to topography and precipitation. In addition, both equal-weight and non-equal-weight methods were used to average the multi-model ensembles so as to analyze the performance of individual models (multi-model ensemble) in terms of extreme precipitation in each region. 

Due to the sparse site and poor temporal persistence of observation data in ACA, the application of grid data makes up for these shortcomings. Previous studies have shown that GPCC is the most accurate in the global gridded precipitation dataset for average and extreme precipitation in Central Asia, outperforming other precipitation products on both the daily and annual scales [41] [42] [43] [44]. As a result, it is widely used for precipitation and extreme precipitation studies in Central Asia [39, 40, 45]. The structure of this paper is as follows: The second section introduces the reference dataset and the CMIP6 models dataset used in this paper as well as the research methods. Section 3 displays the performance evaluation analysis of GCMs and the multi-mode ensemble average simulation. The discussion and the conclusions are provided in Sections 4 and 5, respectively.
section_name: Introduction

section_name: Materials and Methods

section_content: The study region of this paper is Arid Central Asia, which comprises five Central Asian countries (Kazakhstan, Kyrgyzstan, Tajikistan, Turkmenistan and Uzbekistan) and Xinjiang, China. The topography of this region is complex, with the terrain being high in the east and low in the west, and the elevations range from -226 m to 8491 m (Figure 1 ). Located in the hinterland of Eurasia, it is one of the largest arid and most semi-arid regions on Earth. Most of the region has a typical continental arid and semi-arid climate, which is very sensitive to the regional climate. Precipitation is scarce and is greater in mountainous areas than in basins. The daily variation (range) in temperature is relatively large (around 20. 0-30. 0 • C), while the annual change (fluctuation) is relatively small. From the 1960s to 2000s, the average annual temperature increased significantly by 0. 30 • C/10 a. In recent years, the rate of increase has been increasing, especially in the spring [46]. Over the past hundred years, the average near-surface wind speed has been high in the arid northwest plain of Central Asia (>4. 2 m/s), and the mean near-surface wind speed in the southeastern mountainous area has been low (<3. 6 m/s) [47]. The decline rate of wind speed from 1978 to 2014 was close to the global decline rate (<-0. 03 m/s/10a) [48]. The maximum value of monthly evapotranspiration occurs in May, and the minimum value occurs in October, with waters and bare soil taking the lead as main ET contributors [49]. 

Located in the hinterland of Eurasia, it is one of the largest arid and most semi-arid regions on Earth. Most of the region has a typical continental arid and semi-arid climate, which is very sensitive to the regional climate. Precipitation is scarce and is greater in mountainous areas than in basins. The daily variation (range) in temperature is relatively large (around 20. 0-30. 0 °C), while the annual change (fluctuation) is relatively small. From the 1960s to 2000s, the average annual temperature increased significantly by 0. 30 °C/10 a. In recent years, the rate of increase has been increasing, especially in the spring [46]. Over the past hundred years, the average near-surface wind speed has been high in the arid northwest plain of Central Asia (>4. 2 m/s), and the mean near-surface wind speed in the southeastern mountainous area has been low (<3. 6 m/s) [47]. The decline rate of wind speed from 1978 to 2014 was close to the global decline rate (<-0. 03 m/s/10a) [48]. The maximum value of monthly evapotranspiration occurs in May, and the minimum value occurs in October, with waters and bare soil taking the lead as main ET contributors [49]. Arid Central Asia is far away from the ocean. Precipitation is influenced by largescale circulation in addition to complex topography, which is an important factor contributing to the spatial heterogeneity of precipitation in this region [31]. To facilitate spatial analysis, this paper refers to the study of Guo et al. [50], which divided Arid Central Asia into five sub-regions based on the topographic condition and climate zones: Northern Kazakhstan (NK), the Central Desert (CD), the Tianshan Mountain (TM), Northern Xinjiang (NX) and Southern Xinjiang (SX). Arid Central Asia is far away from the ocean. Precipitation is influenced by large-scale circulation in addition to complex topography, which is an important factor contributing to the spatial heterogeneity of precipitation in this region [31]. To facilitate spatial analysis, this paper refers to the study of Guo et al. [50], which divided Arid Central Asia into five sub-regions based on the topographic condition and climate zones: Northern Kazakhstan (NK), the Central Desert (CD), the Tianshan Mountain (TM), Northern Xinjiang (NX) and Southern Xinjiang (SX).
section_name: Study Region

section_content: The reference dataset used in this study is from GPCC (Full Data Daily Version 2020, GPCC V2020) a daily precipitation grid data product. It is available from https://opendata. dwd. de/climate_environment/GPCC/html/fulldata-daily_v2020_doi_download. html (accessed on 5 November 2021). The dataset is a global land surface precipitation grid dataset obtained by the World Climate Research Program (WCRP) and based on observed precipitation data from nearly 85,000 stations worldwide, with a spatial resolution of 1 • × 1 •, covering the time period 1982-2019 [51]. 

The model data used in this paper are the daily precipitation dataset output from 33 GCMs in CMIP6, considering only the first run of each model (r1i1p1f1) and the availability of data. The data can be downloaded from the open website (https://esgf-node. llnl. gov/search/cmip6/) (accessed on 15 October 2021); the basic information of the selected CMIP6 is shown in Table 1. In order to facilitate subsequent comparison, the bilinear interpolation method was used to uniformly interpolate all GCMs into the same standard grid as GPCC, i. e., 1 • × 1 •. Considering the historical period covered by GPCC (1982 GPCC ( -2019) ) and CMIP6 (1850-2014), the research period of this paper is 1982-2014.
section_name: Data

section_content: Based on the extreme indices defined by the Expert Team on Climate Change Detection and Indices (ETCCDI, http://etccdi. pacificclimate. org/indices. shtml) (accessed on 20 September 2021), combined with the characteristics of the study region, eight extreme precipitation indices were selected in terms of precipitation intensity, precipitation frequency, precipitation duration and percentile-based thresholds to investigate extreme precipitation over ACA (Table 2 ). All indices are calculated on an annual time scale. For model comparison, the spatial percentage bias (PBIAS) between the eight precipitation indices simulated by 33 GCMs and the precipitation indices observed by GPCC was calculated (Equation (1)), and the ability of individual models to reproduce the ACA extreme precipitation mean fields was evaluated by using the Taylor diagram. The Taylor diagram also takes into account the spatial correlation coefficient (CC), centralized root-mean-square error (RMSE) and normalized standard deviation ratio (σ) of the model and reference data [52], which are commonly used to measure the agreement between simulated and observed fields [37, 53]. 

where N is the time series of ACA grid points (here N = 33), and M n,i and O n,i are the values of the extreme precipitation index i in period n of the model and GPCC, respectively. 

(2) Interannual Variability Skill Score (IVS) In terms of temporal simulation ability, the variation trend and interannual variability skill score (IVS) are mainly used to evaluate the simulation ability of the model [54]. The IVS is calculated as follows:

where STD M and STD O represent the interannual standard deviation of simulated and observed fields, respectively. The smaller the IVS value, the smaller the standard deviation between the simulated and observed fields, and the better the model simulation ability. 

(3) Comprehensive rating index (MR) We used a comprehensive rating index MR to evaluate the spatio-temporal simulation ability and overall performance of the model [55]. MR was mainly used in two places. First, MR 1 was calculated for 33 GCMs based on the three metrics (CC, RMSE, σ) used in the Taylor diagram. In addition, MR 2 was calculated for each mode IVS. MR is calculated as follows:

where n denotes the number of metrics (n = 3 in MR 1, n = 1 in MR 2 ), m indicates the number of GCMs involved in evaluation (m = 33 in this study) and RANK i is the ranking of 33 GCMs, ranging from 1 to 33. The closer MR is to 1, the better the simulation ability is.
section_name: Methodology 2.3.1. Extreme Precipitation Indices

section_content: The results of two unweighted and two weighted multi-model ensemble methods are compared to find the optimal scheme. The weights of 33 GCMs were calculated separately for eight extreme precipitation indices in the historical period. The better the pattern performance in the weighted multi-model ensemble, the higher the weight. 

(1) Unweighted multi-model ensemble methods (MME and AMME)

In this paper, we choose two equal-weight multi-model ensemble methods, namely, simple arithmetic mean AMME and median MME. Equal weight means that all GCMs are assigned the same weight. For AMME, the model weight coefficients w 1 = w 2 = • • • w 33 = 1; median multi-model ensemble is to sort each extreme precipitation index of 33 GCMs in the order from smallest to largest, and the index value of the 17th model is taken as the result of the multi-model ensemble. 

(2) Performance-based multi-model ensemble (MM-PERF) Four metrics (CC, RMSE, RB and σ) were used to calculate the rank of the individual model at each extreme index [56], which was calculated as follows:

Considering that the optimal RMSE and RB expected results should be close to zero, ranks are equal to the inverse of the absolute value (Equations ( 5 ) and ( 6 )). The optimal ratio of the normalized standard deviation between the model and the observation is 1; the inverse is taken when the ratio is greater than 1. To avoid a negative correlation between the model and observation, the rank of all model correlation coefficients is given plus 1 (Equation ( 4 )). 

where N is the model/observation time steps, M n and O n denote the model simulation and observation value, respectively, M and O refer to the mean value of model simulation and observation, respectively, m denotes the number of GCMs involved in the evaluation and i and j represent CMIP6 models. 

The ranks of the four metrics are normalized so that the sum of the rank of each metric is 1. Weights of the model are obtained by taking the mean of the four metrics and then normalizing them so that the sum of the weights is 1. 

(3) Independence weighted mean method (IWM) Bishop and Abramowitz developed a new independence weighted mean method (IWM); the central idea of this method is to minimize the mean square difference (MSD) between simulated and observed values by finding linear combinations of an ensemble. For more details, see Bishop and Abramowitz [22] and Bai et al. [23]. 

In Equation ( 10 ), (1,. . ., j,. . ., J) is the time steps of indices, (1,. . ., k,. . ., K) represents the selected CMIP6 GCMs in this paper, µ j e and y i denote the extreme indices of the multimodel ensemble and the extreme indices of the observation in the j time steps, respectively, w is the weight of each model from multi-model ensemble and x i is the model coefficient in linear composition. Moreover, to ensure that ∑ K k=1 w k = 1, a Lagrange multiplier (λ) is introduced as a constraint condition:

The calculation of w in the formula is based on the following formula:

where 1 T = [1, 1,. . ., 1]; A is the sample-based estimate of the covariance of the biascorrected errors between all of the ensemble members.
section_name: Multi-Model Ensemble Methods

section_content: J j=1 x jy j x jy j T J -1 (13)
section_name: A = ∑

section_content: The following is a flowchart of the data processing process and methods (Figure 2 ). 

Remote Sens. 2023, 15, 2376 8 of 26
section_name: Methodological Flowchart

section_content: The following is a flowchart of the data processing process and methods (Figure 2 ).
section_name: Methodological Flowchart

section_name: Results

section_name: Spatial Evaluation

section_content: Figures 3 4 5 6 show the percent bias of the spatial field over ACA extreme precipitation indices simulated by 33 GCMs during 1982-2014. The results show that most models can capture the spatial distribution of extreme precipitation better in NK and CD. In SX, almost all models simulate extreme precipitation indices with different degrees of wet bias, especially in the Kunlun Mountains region. Among the eight extreme indices, CMIP6 has the best ability to simulate SDII, with the percentage of spatial bias basically controlled at ±30% (Supplementary Material Figure S4 ). We also found that models from the same institution show similar regional biases and magnitudes in space. For CDD (Figure 3 ), nearly half of the models could depict the spatial distribution of ACA, among which ACCESS-CM2, KIOST-ESM, MRI-ESM2-0 and TaiESM1 showed the best spatial performance, while AWI-ESM-1-1-LR, IITM-ESM, MPI-ESM-1-2-HAM, MPI-ESM1-2-HR, MPI-ESM1-2-LR and NESM3 showed great positive percentage biases in TM, Altai Mountains in NX and the western and the southwestern of CD. The models showed similar spatial biases for PRCPTOT (Figure 4 ), CWD and R95pTOT simulations, with high positive percentage biases in most areas (Supplementary Materials Figures S1 and S2 ). In contrast, most GCMs have high negative percentage biases for R10mm simulations, with dry bias, which substantially underestimates heavy rainfall in the region (Figure 5 ). Most models have a good simulation ability for Rx1day and are able to capture Rx1day, even For CDD (Figure 3 ), nearly half of the models could depict the spatial distribution of ACA, among which ACCESS-CM2, KIOST-ESM, MRI-ESM2-0 and TaiESM1 showed the best spatial performance, while AWI-ESM-1-1-LR, IITM-ESM, MPI-ESM-1-2-HAM, MPI-ESM1-2-HR, MPI-ESM1-2-LR and NESM3 showed great positive percentage biases in TM, Altai Mountains in NX and the western and the southwestern of CD. The models showed similar spatial biases for PRCPTOT (Figure 4 ), CWD and R95pTOT simulations, with high positive percentage biases in most areas (Supplementary Materials Figures S1 and S2 ). In contrast, most GCMs have high negative percentage biases for R10mm simulations, with dry bias, which substantially underestimates heavy rainfall in the region (Figure 5 ). Most models have a good simulation ability for Rx1day and are able to capture Rx1day, even in the complicated terrain of southern Xinjiang (Supplementary Material Figure S3 ). The higher positive bias for Rx5day is mainly in the Kunlun Mountains as well as in central Kazakhstan (Figure 6 ). 

In general, these 33 GCMs have a certain ability to simulate the overall spatial distribution characteristics of extreme precipitation in ACA, with various performances in terms of precipitation intensity, precipitation duration and precipitation frequency. Unfortunately, no individual model was able to accurately capture all the spatial characteristics of the eight extreme precipitation indices in the ACA and each sub-region. Except for R10mm, EC-Earth3, EC-Earth3-CC, EC-Earth3-Veg-LR, KIOST-ESM and NorESM2-MM are the most powerful in simulating the other seven extreme precipitation indices.
section_name: Spatial Bias Analysis

section_content: To evaluate the ability of GCMs to simulate the spatial pattern of extreme precipitation indices, the Taylor diagram is further presented (Figure 7 ). It is found that CMIP6 models have the strongest ability to simulate the mean field of CDD and SDII; the CC of most models and the reference data are above 0. 5, and the maximum RMSE is less than 1. 3. For CWD, the consistency between model simulation and observation is poor on the whole, with the maximum CC being less than 0. 6 (mostly between 0. 1 and 0. 4), and the σ and RMSE of most models are above 3. Except for ACCESS-ESM1-5, the remaining models have slight differences in terms of simulation ability, so they are relatively concentrated in the Taylor diagram. Both IPSL-CM6A-LR and IPSL-CM6A-LR-INCA show higher CC values with respect to the reference data in terms of extreme precipitation intensity (PRCPTOT, Rx1day, Rx5day), R95pTOT and extreme precipitation frequency (R10mm) in ACA, but they also have larger σ and RMSE values at the same time. The ability to simulate the R10mm spatial field varies, with the CC between models and the reference data ranging from 0. 1 to 0. 7, the RMSE being basically below 1. 5 and σ ranging from 0. 3 to 2. 6. It is not hard to find that most GCMs have a similar spatial field simulation ability for Rx1day and Rx5day, with σ close to 1 and a smaller RMSE. 

The MR 1 of 33 GCMs was calculated using three indices (CC, RMSE and σ) in the Taylor diagram. Figure 8 shows the ranking of the model's ability to simulate the spatial field of eight extreme indices. The darker the red, the higher the value of MR 1 and the better the model performs. The results show that each model has a different ability to simulate the spatial fields of extreme precipitation indices in different regions. From the perspective of extreme precipitation indices, GCMs are more sensitive to the simulation of CDD and CWD in ACA and its sub-regions, with great uncertainty. 

In terms of sub-regions, the model performs more stable extreme precipitation indices in NK. EC-Earth3-CC and GFDL-ESM4 have the best overall simulation ability performance, while NorCPM1 and MPI-ESM1-2-LR have the worst performance. In SX, TM and ACA, there are large differences in the spatial simulation abilities of the models for different extreme precipitation indices, such as BCC-ESM1 and EC-Earth3-Veg. Based on MR 1, the models with the best performance in simulating the spatial mean field of extreme precipitation indices at CD, NK, NX, SX, TM and ACA are GFDL-ESM4, EC-Earth3-CC, CESM2-WACCM, BCC-ESM1, EC-Earth3-AerChem and EC-Earth3-AerChem, respectively. in the Taylor diagram. Both IPSL-CM6A-LR and IPSL-CM6A-LR-INCA show higher CC values with respect to the reference data in terms of extreme precipitation intensity (PRCPTOT, Rx1day, Rx5day), R95pTOT and extreme precipitation frequency (R10mm) in ACA, but they also have larger σ and RMSE values at the same time. The ability to simulate the R10mm spatial field varies, with the CC between models and the reference data ranging from 0. 1 to 0. 7, the RMSE being basically below 1. 5 and σ ranging from 0. 3 to 2. 6. It is not hard to find that most GCMs have a similar spatial field simulation ability for Rx1day and Rx5day, with σ close to 1 and a smaller RMSE.
section_name: Taylor Diagram-Based Spatial Simulation Ability Analysis

section_content: The ability of the model to simulate the trend change in the extreme precipitation indices is an important aspect in testing the performance of the model. Figure 9 displays the annual variation trend of extreme precipitation indices simulated by 33 GCMs and the multi-model ensemble average. Table 3 shows the correlation coefficients of the regional average between 33 CMIP6 models and GPCC. 3 shows the correlation coefficients of the regional average between 33 CMIP6 models and GPCC. Compared with the reference data, most CMIP6 models are able to detect the increasing trends of CWD, PRCPTOT, R10mm, R95pTOT, Rx1day, Rx5day and SDII, but for CDD, more than half of the simulated values of GCMs are contrary to the trends observed by GPCC. In terms of the magnitude of trend deviation, the vast majority of models un- Compared with the reference data, most CMIP6 models are able to detect the increasing trends of CWD, PRCPTOT, R10mm, R95pTOT, Rx1day, Rx5day and SDII, but for CDD, more than half of the simulated values of GCMs are contrary to the trends observed by GPCC. In terms of the magnitude of trend deviation, the vast majority of models underestimate the growth trends of CDD (0. 9 days/10a), R10mm (0. 2 days/10a), Rx1day (0. 3 mm/10a) and SDII (0. 07 mm/day/10a). About one-third of the models overestimate the trends of CWD (0. 008 days/10a), PRCPTOT (2. 5 mm/10a), R95pTOT (3. 2 mm/10a) and Rx5day (0. 4 days/10a). Compared with the spatial pattern, models have a relatively poor ability to simulate the trends of extreme precipitation indices. Taking CDD as an example, CESM2-FV2, E3SM-1-0 and AWI-ESM-1-1-LR all have high correlation coefficients with the reference data (0. 43, 0. 31 and 0. 37, respectively), but AWI-ESM-1-1-LR largely overestimates the variation trend of CDD, and CESM2-FV2 underestimates the trend of CDD. In contrast, E3SM-1-0 has a higher correlation coefficient and is closest to the variation trend of GPCC. Therefore, we believe that E3SM-1-0 has the best ability to simulate the trend change of CDD. For the remaining seven extreme precipitation indices, the models with the best simulation abilities are ACCESS-CM2, INM-CM4-8, EC-Earth3-CC, INM-CM4-8, MPI-ESM1-2-HR, MRI-ESM2-0 and IITM-ESM.
section_name: Temporal Evaluation 3.2.1. Trend Evaluation

section_content: According to Equation (2) in Section 2. 3. 2, we calculate the interannual variability skill score of the model spatial field relative to the GPCC observation field (Figure 10 ). Similar to the simulation results in the Taylor diagram, the CMIP6 model has the worst simulation of CWD in ACA (IVS between 0 and 20. 4), followed by R10mm (IVS between 0 and 12. 3). On the contrary, PRCPTOT, R95pTOT and Rx5day have better simulations of interannual variability. Furthermore, it is difficult to capture the interannual variability of the SDII observation field in ACA, and the difference between the models is not significant. 

precipitation indices, the models with the best simulation abilities are ACCESS-CM2, INM-CM4-8, EC-Earth3-CC, INM-CM4-8, MPI-ESM1-2-HR, MRI-ESM2-0 and IITM-ESM.
section_name: Evaluation of Interannual Variability

section_content: According to Equation (2) in Section 2. 3. 2, we calculate the interannual variability skill score of the model spatial field relative to the GPCC observation field (Figure 10 ). Similar to the simulation results in the Taylor diagram, the CMIP6 model has the worst simulation of CWD in ACA (IVS between 0 and 20. 4), followed by R10mm (IVS between 0 and 12. 3). On the contrary, PRCPTOT, R95pTOT and Rx5day have better simulations of interannual variability. Furthermore, it is difficult to capture the interannual variability of the SDII observation field in ACA, and the difference between the models is not significant. From the point of view of each sub-region, the interannual variability skill scores of models are relatively small in NK and CD, with the IVS of extreme precipitation indices ranging from 0 to 13 and from 0 to 8. 5, respectively, except for IITM-ESM. GCMs have the worst ability to simulate the interannual variability of extreme precipitation indices in SX, especially in CWD and SDII. The IVS values of CWD simulated by the model are relatively low in NX but have limited ability to capture interannual variability for R10mm and Rx1day, with significant differences in the simulation ability among models. GCMs perform better overall in reproducing interannual variability for PRCPTOT and R95pTOT in TM, with IVS values not exceeding 2. 5. 

In general, CMIP6 models are basically able to reproduce the interannual variability of PRCPTOT, R95pTOT and Rx5day in ACA and sub-regions, but the interannual variability of R10mm, SDII and CWD need to be improved. The best-performing models in capturing the interannual variability of extreme precipitation indices in CD, NK, NX, SX, From the point of view of each sub-region, the interannual variability skill scores of models are relatively small in NK and CD, with the IVS of extreme precipitation indices ranging from 0 to 13 and from 0 to 8. 5, respectively, except for IITM-ESM. GCMs have the worst ability to simulate the interannual variability of extreme precipitation indices in SX, especially in CWD and SDII. The IVS values of CWD simulated by the model are relatively low in NX but have limited ability to capture interannual variability for R10mm and Rx1day, with significant differences in the simulation ability among models. GCMs perform better overall in reproducing interannual variability for PRCPTOT and R95pTOT in TM, with IVS values not exceeding 2. 5. 

In general, CMIP6 models are basically able to reproduce the interannual variability of PRCPTOT, R95pTOT and Rx5day in ACA and sub-regions, but the interannual variability of R10mm, SDII and CWD need to be improved. The best-performing models in capturing the interannual variability of extreme precipitation indices in CD, NK, NX, SX, TM and ACA are: MPI-ESM1-2-HR, ACCESS-CM2, GFDL-CM4, NorESM2-MM, MPI-ESM1-2-HR and FGOALS-f3-L, respectively.
section_name: Evaluation of Interannual Variability

section_content: Given the inconsistency of spatial and temporal model rankings, in order to more intuitively determine the spatio-temporal simulation ability of the models, this study calculated MR based on Taylor and IVS for a comprehensive evaluation of the models (Figure 11 ). The study calculated MR based on Taylor and IVS to more intuitively determine the spatio-temporal simulation ability of the models. Taylor diagram-based MR 1 and IVS based MR 2 are significantly correlated in NK, NX and ACA; the correlation coefficients are 0. 57, 0. 69 and 0. 42, respectively (all passed the 95% significance test), while the spatial pattern and interannual variability of extreme precipitation indices simulated by the models in CD, SX and TM were not consistent, with correlation coefficients of 0. 32, 0. 21 and 0. 26, respectively. The points on the plane were divided into four quadrants. Models in the upper right quadrant have better spatial and temporal simulation ability, while models in the lower left quadrant have the worst simulation ability. The MR 1 and MR 2 correlations of NK and TM are not high, but the best-performing models are more concentrated. 

Considering the simulation capability in terms of spatial and temporal variability, ACCESS-CM2, ACCESS-ESM1-5, GFDL-CM4 and GFDL-ESM4 perform best in CD, among which ACCESS-CM2 and GFDL-CM4 also perform better in NK and are more suitable for extreme precipitation simulation in flat areas. In addition, models of the EC-Earth-Consortium institution (e. g., EC-Earth3-Veg, etc. ) all perform better in NK. GFDL-CM4, GFDL-ESM4 and CESM2-WACCM are more suitable for Xinjiang (SX and NX). For mountainous regions (TM), the models simulated extreme precipitation indices that are less consistent than those in other regions in terms of the spatial pattern and interannual variability, and the best-simulated models are NESM3, AWI-ESM-1-1-LR, MPI-ESM1-2-LR, MPI-ESM-1-2-HAM and E3SM-1-0. Across the ACA, the best-performing models are GFDL-CM4, GFDL-ESM4, EC-Earth3-Veg and EC-Earth3-AerChem.
section_name: Overall Model Performance

section_content: Compared with individual models, multi-model ensemble methods have significantly better simulation ability for extreme precipitation spatial fields in ACA. Except for CWD, the CCs of the four multi-model ensemble simulations with reference data are basically above 0. 6, and the RMSE is significantly reduced, while σ is closer to 1 (Figure 7 ). Figure 12 shows the spatial percentage biases of eight extreme precipitation indices simulated by four multi-model ensemble methods versus those observed by GPCC. 

Due to the inherent limitations of GCMs, even with different multi-model ensemble methods, the improvement of CWD, PRCPTOT, R10mm and R95pTOT is still very limited, and there are still very large wet biases in the southeastern TM and the southern and central SX. Two unweighted schemes have a similar simulation performance on CDD, CWD, PRCPTOT and R95pTOT. AMME is better than MME for R10mm, but MME is better than AMME in simulating extreme precipitation intensity SDII. Among the weighted multimodel ensemble methods, the spatial bias of the extreme precipitation indices simulated by MM-PERF is slightly smaller than that of unweighted AMME in ACA. With the exception of CDD, which has large biases in TM and NX, IWM-simulated extreme precipitation indices have an overall smaller spatial percentage bias than the other three multi-model ensembles. 

In conclusion, four multi-model ensemble methods have a certain ability to simulate extreme precipitation over ACA, and the performance of weighted multi-model ensemble is better than that of equal-weight multi-model ensemble. IWM has a relatively better spatial simulation ability for extreme precipitation indices. 0. 26, respectively. The points on the plane were divided into four quadrants. Models in the upper right quadrant have better spatial and temporal simulation ability, while models in the lower left quadrant have the worst simulation ability. The MR1 and MR2 correlations of NK and TM are not high, but the best-performing models are more concentrated. Multi-model ensemble can effectively reduce the systematic bias and variability within individual models. By comparing four multi-model ensemble methods, simulations of different extreme precipitation indices are improved to varying degrees (Figure 13 ). Compared to individual model simulations (Supplementary Material Figure S5 ), the bias of extreme precipitation indices across ACA was reduced by at least 50%, which is evident in both PRCPTOT and R95pTOT simulations. 

Compared with individual models, multi-model ensemble methods have significantly better simulation ability for extreme precipitation spatial fields in ACA. Except for CWD, the CCs of the four multi-model ensemble simulations with reference data are basically above 0. 6, and the RMSE is significantly reduced, while σ is closer to 1 (Figure 7 ). Figure 12 shows the spatial percentage biases of eight extreme precipitation indices simulated by four multi-model ensemble methods versus those observed by GPCC.
section_name: Multi-Modal Ensemble Performance Evaluation

section_content: By comparing the spatio-temporal simulation performance between CMIP6 models, it is found that most models have a large wet bias for TM and SX, which is consistent with the findings of Dong et al. [9, 50], while the simulation ability is better in plain regions. Liu et al. found that it was difficult for the CMIP6 model to simulate the overall spatial pattern of extreme precipitation in early spring and late winter, the inter-model variability of regional bias in summer was large and the areas of inconsistency that emerge always contain most of Xinjiang or its surrounding mountain regions [57]. This may be due to the From the perspective of correlation, the weighted multi-model ensemble simulation effect has improved significantly, and the correlation coefficient can reach up to 0. 7. In terms of the multimodal ensemble method performance, weighted multi-model ensemble effects are better than equal-weighted ensemble effects, and the temporal variation trend of IWM is closest to the reference data (Figure 8 ). The uncertainty of simulated precipitation duration (CDD, CWD), PRCPTOT and R95PTOT in MME was slightly less than that in AMME but less than that in R10mm, Rx1day, Rx1day and SDII. Although the MM-PERF method can also reduce the uncertainty of extreme precipitation simulation, the improvement effect is not as good as that of IWM. The results of the partitioning show that there are still great uncertainties in the simulation of extreme precipitation in SX and TM.
section_name: Discussion

section_content: By comparing the spatio-temporal simulation performance between CMIP6 models, it is found that most models have a large wet bias for TM and SX, which is consistent with the findings of Dong et al. [9, 50], while the simulation ability is better in plain regions. Liu et al. found that it was difficult for the CMIP6 model to simulate the overall spatial pattern of extreme precipitation in early spring and late winter, the inter-model variability of regional bias in summer was large and the areas of inconsistency that emerge always contain most of Xinjiang or its surrounding mountain regions [57]. This may be due to the complex conditions in the Tianshan Mountains and Southern Xinjiang, and the coarse resolutions of CMIP6 models make it difficult to capture the precipitation variability in this region. For extreme precipitation indices, most models overestimate the mean and variability of total precipitation PRCPTOT and CWD but underestimate the mean and variability of CDD and SDII, and a similar situation was found in East Africa and contiguous regions of the United States [58, 59]. Luo and Guo [60] pointed out that improved resolution has obvious advantages in simulating extreme precipitation, especially in some areas with complex topography. However, this explanation does not seem to hold up in some tropical regions, such as East Africa [58]. In this study, some models with a lower resolution are better than those with a higher resolution (such as IITM-ESM) in capturing the spatial distribution of extreme precipitation. As a result, the resolution of the model cannot fully explain the performance of the model, but the setting of the parameterization scheme of the model itself, the physical process and other factors should also be taken into account [61, 62]. 

Extreme precipitation is generally related to the convection and cloud microphysics [63]. It is difficult to attribute the difference in the simulation ability of each model in different regions. To reduce the uncertainty of future climate projection, scholars have been working on an ideal model ensemble method, but there is still a dispute over the weight of climate models in the world [56, 64]. Moreover, a weighted multi-model ensemble approach is very sensitive to the indicators used, and the selection of indexes can directly lead to the quality of simulation results. In this paper, four statistical indicators (CC, RMSE, RB and σ) are mainly used in the MM-PERF in the weighted multi-model ensemble method, and the IWM scheme outperforms the MM-PERF overall. However, if other indicators or methods are used, the results may be different. For example, Zhang et al. used a weighted multi-model ensemble to project the changes in extreme precipitation in Central Asia under different temperature rises and found that there was a robust consistency between the changes in precipitation intensity and the cumulative one, but the noise generated by the uncertainty between the models covered the signals of CDD changes [39]. In addition, the spread of summer precipitation extremes over northern Central Asia is large [40]. The simulation ability of models should be understood in the context of specific variables and spatio-temporality [65, 66]. 

Currently, the improvement of the model simulation ability is mainly focused on reducing uncertainty in climate change prediction by establishing some relationship between observations and simulations [67, 68]. However, this method is vulnerable to the influence of the quality of observed data, physical processes and even the subjective consciousness of researchers [64]. In addition, a simple assessment of the model's performance in the current climate is not sufficient to judge the model's performance [69]. Hence, the in-depth understanding and elimination of the CMIP6 models' bias is the key to improving the accuracy of regional climate projection, and it is also the focus of future research efforts. Zhang et al. used downscaling and two methods of local intensity scaling (LOCI) and quantile mapping (QM) and greatly reduced the deviation in extreme precipitation simulation [70]. Other multi-model ensemble methods or bias correction of the data, such as machine learning and downscaling, could be tried in the future to reduce model uncertainty in simulations and projections.
section_name: Discussion

section_content: Funding: This work was supported by the National Natural Science Foundation of China ( 42067062 ). The funder is Changchun Xu, and the funding number is 42067062.

section_name: Data Availability Statement:

section_content: The CMIP6 dataset that supports the findings of this study is openly available at https://esgf-node. llnl. gov/search/cmip6/ (accessed on 15 October 2021). GPCC data are available from https://opendata. dwd. de/climate_environment/GPCC/html/fulldata-daily_v2020_ doi_download. html (accessed on 5 November 2021). The vector data of administrative regionalization in Arid Central Asia are provided by the National Cryosphere Desert Data Center (http://www. ncdc. ac. cn) (accessed on 18 October 2021). The extreme index calculations defined by ETCCDI can be obtained at http://etccdi. pacificclimate. org/indices. shtml (accessed on 20 September 2021).
section_name: Data Availability Statement:

section_content: The authors declare no conflict of interest.
section_name: Conflicts of Interest:

section_name: Conclusions

section_content: In this paper, based on the GPPC V2020 daily precipitation data products and 33 GCMs datasets, eight extreme precipitation indices defined by ETCCDI were selected to quantitatively and comprehensively evaluate the spatial and temporal simulation ability of the models for extreme precipitation in ACA and each sub-region. By comparing and analyzing the simulation performance of four multi-model ensemble averaging methods and individual models, the models and multi-model ensemble averaging methods with relatively reliable simulation effects are selected. The main conclusions are as follows:

(1) In terms of spatial simulation, most models can simulate the spatial characteristics of extreme precipitation indices well. However, the simulation ability of the models varies for different indices. The models have a stronger simulation ability for CDD and SDII, and the correlation coefficients between the climate mean fields of extreme precipitation indices simulated by most of the models and the extreme precipitation indices observed by the reference data exceed 0. 6. Relatively speaking, the models have poor simulation ability for CWD and R10mm, the consistency of spatial distribution characteristics with the reference data is weak and the simulation ability varies greatly between models. 

(2) As for temporal variability simulation ability, the linear variation trends of the extreme precipitation indices simulated by the models are mostly consistent with the reference data in ACA, but the simulation ability is weak in terms of trend intensity. As for the interannual variability skill scores of each sub-region, the low IVS scores (0-8. 5 and 0-13, respectively) of the models in CD and NK indicated a relatively good simulation. However, it performs the worst in terms of reproducing the interannual variability of extreme precipitation in the SX. GCMs have limited ability to capture the interannual variability of R10mm, SDII and CWD in the ACA and each sub-region. 

(3) Considering the ability of the models to simulate the spatial pattern and interannual variability of the extreme precipitation indices, ACCESS-CM2 and GFDL-CM4 performed better in the relatively flat CD and NK. NESM3, AWI-ESM-1-1-LR, MPI-ESM1-2-LR, MPI-ESM-1-2-HAM and E3SM-1-0 have relatively good simulation capability in the TM region. For the whole ACA, the best-performing models are GFDL-CM4, GFDL-ESM4, EC-Earth3-Veg and EC-Earth3-AerChem. 

(4) Comparing the effects of four different multi-model ensemble averaging methods, it is found that the multi-model ensemble averaging methods have a significant improvement in the spatial and temporal simulation ability of extreme precipitation indices and can effectively reduce the internal bias and variability of individual models (especially in CD and NK). However, there are still large uncertainties in the extreme precipitation simulation in SX and TM. The non-equally weighted model ensemble averaging method is superior to the unweighted multi-model ensemble methods, among which IWM has the most obvious enhancement of the simulation ability of extreme precipitation indices and the best performance.
section_name: Conclusions

section_content: The following supporting information can be downloaded at: https: //www. mdpi. com/article/10. 3390/rs15092376/s1. Figure S1. Percentage bias of CWD simulated by 33 GCMs relative to the CWD observed by GPCC in the ACA and sub-regions from 1982 to 2014 (unit: %). Figure S2. Percentage bias of R95pTOT simulated by 33 GCMs relative to the R95pTOT observed by GPCC in the ACA and sub-regions from 1982-2014 (unit: %). Figure S3. Percentage bias of Rx1day simulated by 33 GCMs relative to the Rx1day observed by GPCC in the ACA and sub-regions from 1982 to 2014 (unit: %). Figure S4. Percentage bias of SDII simulated by 33 GCMs relative to the SDII observed by GPCC in the ACA and sub-regions from 1982 to 2014 (unit: %). Figure S5. Box-and-whisker plots of regional mean bias (model minus GPCC) of eight extreme precipitation indices (CDD, CWD, PRCPTOT, R10mm, R95pTOT, Rx1day, Rx5day, SDII) simulated by 33 GCMs relative to observed by GPCC from 1982 to 2014 in ACA. 

Author Contributions: Conceptualization, C. X. ; Data curation, C. X. and X. L. ; Formal analysis, X. L. ; Methodology, X. L., L. S., L. C. and N. S. ; Software, X. L. and F. L. ; Visualization, X. L. ; Writing-original draft, X. L. All authors have read and agreed to the published version of the manuscript.
section_name: Supplementary Materials: