Projection of future precipitation change using CMIP6 multimodel ensemble based on fusion of multiple machine learning algorithms: A case in Hanjiang River Basin, China

Projecting precipitation changes is essential for researchers to understand climate change impacts on hydrological cycle. This study projected future precipitation over the Hanjiang River Basin (HRB) based on the multimodel ensemble (ME) of six global climate models from Phase 6 of the Coupled Model Intercomparison Project (CMIP6). An ME method using the fusion of four machine learning (ML) algorithms (random forest [RF], K‐nearest neighbors [KNN], extra tree [ET], and gradient boosting decision tree [GBDT]) was proposed in this study. The future precipitation changes were investigated during 2023–2042 (Near‐term), 2043–2062 (Mid‐term), and 2081–2100 (Long‐term) periods, with reference to the base period 1995–2014, under three integrated scenarios (SSP1‐2.6, SSP2‐4.5, and SSP5‐8.5) of the Shared Socioeconomic Pathways (SSPs) and the representative concentration pathways (RCPs). The results show that: (1) the proposed ME method performs better than the ME mean and individual ML algorithms, with a correlation coefficient value reaching 0.88 and Taylor skill score reaching 0.764. (2) The precipitation under SSP5‐8.5 has the largest upward trend with the annual precipitation variation range of −9.27% to 112.84% from 2023 to 2100, followed by SSP2‐4.5 with −30.48% to 44.67%, and the smallest under SSP1‐2.6 with −37.19% to 37.78%, which show a significant trend of humidification over the HRB in the future. (3) The precipitation changes over the HRB are projected to increase over time, with the largest in the Long‐term, followed by Mid‐term, and the smallest in the Near‐term. (4) The northeastern parts of the HRB are projected to experience a large precipitation in the future, and the southeastern parts are smaller. (5) Uncertainties in the projected precipitation over the HRB still exist, which can be reduced by ME. The findings obtained in this study have important implications for hydrological policymakers to make adaptive strategies to reduce the risks of climate change.


| INTRODUCTION
Global climate change has become one of the most serious problems all over the world (Cook et al., 2020;Z. Li & Fang, 2016).By affecting meteorological and hydrological elements such as precipitation, evaporation, and runoff, climate change alters the water cycle of the region, directly or indirectly.This results in a more complex evolution of regional extreme hydrological events, which pose threats to public health, food security, flood control, water supply, and environmental security at the regional scale (Meaurio et al., 2017;Schilling et al., 2020;Zheng et al., 2018).Precipitation is a key component of the climate system and a necessary condition for the formation and transformation of runoff.Precipitation changes caused by climate change alter the original rainfall-runoff relationship in the hydrological cycle (L.Li et al., 2017;D. Wang et al., 2021).Therefore, understanding and projecting future precipitation changes are vital for mitigating the negative effects of hydrological fluctuations, reducing flooding risks, managing water resources, and addressing hydroclimatic variabilities (Z.Li, Li, et al., 2020;Xiong et al., 2017).
Global climate models (GCMs) are one of the most important tools to investigate climate change, which have been widely used especially since the operation of Coupled Model Intercomparison Project (CMIP) in 1995 (Eyring et al., 2016;Kim et al., 2020).The CMIP has developed rapidly since its implementation, and has already experienced five stages of development: Phase 1 of CMIP (CMIP1), CMIP2 (1997), CMIP3 (2004), CMIP4 (2005), andCMIP5 (2013).Recently, CMIP6 is currently in progress (Eyring et al., 2016).A major difference between CMIP6 and previous CMIPs is the future scenarios used to project climate change.Furthermore, compared with previous CMIP, GCMs participating in CMIP6 have higher spatial resolutions, extended historical run time, and more complicated physical processes (Eyring et al., 2016;O'Neill et al., 2013).With the gradual release of CMIP6 GCM data, numerous researchers have explored the future climate changes over different regions of the world (Cook et al., 2020;Hamed et al., 2022;Nashwan & Shahid, 2022;Shuaifeng & Xiaodong, 2022;Yue et al., 2021).
Although significant improvements have been observed in the latest CMIP6, due to incomplete understanding of the climate change system, simplified assumptions about the model structure, uncertainty of model parameters and unrealistic estimates of natural activities, etc., GCMs still have too large systematic biases and uncertainties to project future climate change at a regional scale (L.Chen & Frauenfeld, 2014;Eyring et al., 2016;Reichler & Kim, 2008;D. Wang, Liu, et al., 2022;S. Zhang & Chen, 2021;Zhao et al., 2022).To address this problem, dynamical and statistical downscaling methods were developed.Compared with the dynamical downscaling method, the statistical downscaling method has the advantages of being relatively computationally efficient, simple operation, and low cost.The statistical downscaling method is also called a bias correction technique.Many statistically downscaled methods have been developed and improved, like quantile mapping method (QM), Linear Scaling (LS), delta, transfer cumulative probability distribution method (CDF-t), equidistant cumulative distribution functions matching method (EDCDFm), quantile delta mapping (QDM), etc. (Guo et al., 2018;Patel et al., 2021;Tong et al., 2020).It must be noted that bias-corrected simulations of GCMs still contain some biases, leading to strong uncertainties in GCM projections.In addition, the spatial and temporal variability of precipitation makes accurate projection of precipitation a great challenge, and its simulation accuracy is often lower than that of temperature data (Yue et al., 2021).Therefore, effectively improving the GCMs precipitation simulation accuracy and reducing uncertainty have become the focus of many scholars.
In recent years, using an ensemble of multiple GCMs has become an effective way to reduce the biases and uncertainties of the bias-corrected simulations from individual GCMs, which can improve the reliability in projecting future climate change (Daron et al., 2018;Jose et al., 2022).Several previous studies have introduced multimodel ensemble (ME) methods, which can be categorized into multimodel ensemble mean (MME) and weighted ensemble method (WEM) (Ahmed et al., 2020;Jose et al., 2022).In WEM, each GCM was given different weight to calculate ME according to its ability to simulate historical climates, while the MME is the average of values simulated from all selected GCMs (Oh & Suh, 2016).Many studies have shown that both the above approaches perform better than individual GCM (Oh & Suh, 2016;D. Wang et al., 2021;D. Wang, Liu, et al., 2022;Yue et al., 2021).With the vigorous development of AI algorithms, machine learning (ML) algorithms have begun to be applied to develop MEs and have been proven more effective (Ahmed et al., 2020;Jose et al., 2022;R. Xu et al., 2020;Yu et al., 2022).For example, Crawford et al. (2019) evaluated the abilities of various MEs developed by random forest (RF), support vector regression (SVR), artificial neural networks (ANN), multiple linear regression (MLR), and weighted k-nearest neighbors (KNN) to simulate past temperature and precipitation in the Gulf Basin region of North America, and the results reported that RF, SVR, ANN, and MLR algorithms significantly improve the performances of precipitation simulation, as compared with the individual GCMs and the MME.Ahmed et al. (2020) evaluated the performance of MEs developed using ANN, KNN, SVR, and relevance vector machine (RVM), and found that MEs based on KNN and RVM show better skills.Jose et al. (2022) compared the MME method with various MLs like MLR, SVM, extra tree (ET), RF, and long short-term memory (LSTM) to develop MEs of precipitation, maximum temperature, and minimum temperature over a tropical river basin in India.The results showed that all ML approaches performed better than the MME approach, and the LSTM shows best in the case of precipitation.Acharya et al. (2013) compared the ME results developed by extreme learning machine (ELM), MME, and multiple linear regressions, and found that ELM expeditiously captures the northeast monsoon rainfall over south peninsular India compared with the other ME approaches.B. Wang et al. (2018) considered Bayesian model averaging (BMA), support vector machine (SVM), MME, and RF for ME, and found that RF and SVM are preferred ME approaches.As indicated above, ML algorithms are more effective in developing MEs compared with individual models and other approaches.However, it is still in its infancy to apply ML algorithms to the development of GCMs MEs for climate impact assessment.Additionally, each of the individual models has its own strengths and weaknesses, and has some limitations when it comes to capturing complex relationships in the datasets.Therefore, it is still necessary to tap the potential of various ML algorithms in developing MEs with GCMs simulation and increase the universalities of ML algorithms in different regions of the world.
The Hanjiang River Basin (HRB), the largest river basin over the Yangtze River Basin in China, has suffered considerable economic losses and human deaths due to extreme events, which are closely connected to frequent climate change (X.Zhang et al., 2022).In this study, the main objective is to project the changes of the future precipitation, which is one of the main hydrological elements over the HRB in China based on six GCMs from CMIP6 under future scenarios.To achieve this aim, an ME technique based on the fusion of four ML algorithms (RF, KNN, ET, and gradient boosting decision tree [GBDT]) using voting regressor (VR) algorithm is proposed to develop ME of GCMs simulations.To the best of our knowledge, the implementation of VR to fuse multiple ML algorithms in developing GCMs ME fields is never reported previously.The structure of this article is as follows: Section 2 describes 'Study area and data' by introducing the study area and datasets considered.Section 3 is about 'Methodology', presenting the bias correction method, ML algorithms, GCMs performance metrics, and projected changes analysis method.Section 4 is about 'Results', describing the results of this article's findings, including the comparison of various ME approaches and future precipitation projection over HRB.Section 5 is about 'Discussion', where the findings of this study and the scope for future work were discussed.Section 6 is about 'Conclusions', providing the summary of this article.The framework of this study is shown in Figure 1.

| Study area
The HRB located between 30 80 0 -34 110 0 N and 106 120 0 -114 140 0 E covers approximately 1.59 Â 10 5 km 2 (Figure 2).The Hanjiang River, with a length of 1577 km, is the largest tributary of the Yangtze River, China, flowing through six provinces including Shanxi, Hubei, Henan, Sichuan, Gansu, and Chongqing (Qin et al., 2019;D. Wang, Liu, et al., 2022).The HRB has three distinct levels from west to east, resulting in complex topography.Mountainous and semi-mountainous regions dominate the upper reaches of the HRB, while plain areas cover the middle and lower reaches, ranging in elevation from 0 to 3577 m (Figure 2).There is a subtropical monsoon climate in the HRB, which has an average annual temperature of 15-17 degrees Celsius in the climate transition zone between the north and south of China (Hao et al., 2019).Approximately 70%-80% of the precipitation in the HRB falls between May and October, with the mean annual precipitation range between 700 and 1800 mm.Thus, approximately 75% of the annual runoff is caused by precipitation between May and October (Hao et al., 2019;Qin et al., 2019;D. Wang, Liu, et al., 2022).There are several major water development projects in the HRB, such as Danjiangkou Reservoir, Hanjiang to Weihe River Project, and the Southto-North Water Diversion Project.With the impact of climate change and water diversion projects in recent years, extreme hydrological events have increased in probability and magnitude, resulting in more conflicts in the HRB (Hao et al., 2019;Y. Zhou et al., 2015;Y. P. Zhu et al., 2008).To predict the water availability and provide disaster prevention measures in the future, it is essential to accurately project the future precipitation change in the HRB.

| Observation dataset
The observational dataset used in this study comes from the CN05.1 dataset, which is a gridded daily dataset with a spatial resolution of 0.25 Â 0.25 (J.Wu & Gao, 2013).The CN05.1 dataset released by the National Meteorological Center of China Meteorological Administration is developed by J. Wu and Gao (2013) incorporating ANUS-PLIN and angular distance weights methods based on observations from 2416 national meteorological stations between 1961 and 2021 (D.Wang, Liu, et al., 2022).In addition to precipitation data, the CN05.1 dataset also includes average temperature, maximum and minimum temperatures.With high spatial resolution and complete time series, the CN05.1 dataset covers the entire Chinese mainland (Y.Luo et al., 2020;J. Wu & Gao, 2013).In China, this is the most accurate dataset of gridded nearsurface meteorological fields, and it represents changes in meteorological elements accurately (Y.Luo et al., 2020;D. Wang & Wang, 2017;B. Zhou et al., 2016).A strict quality control process has been implemented for this dataset, including inspections of missing values and extreme value tests.It is noteworthy that the CN05.1 dataset has been extensively used as a reference for evaluating and verifying GCM outputs (N.Luo et al., 2021;D. Wang, Liu, et al., 2022;Y.-Y. Zhu & Yang, 2020).Thus, the CN05.1 dataset is optimal for verifying the accuracy of GCM precipitation simulations.Comparison and verification of the GCMs precipitation simulation outputs were conducted using daily precipitation data from the CN05.1 dataset dating from 1970 to 2015 in this study.

| CMIP6 model datasets
In our previous study, we evaluated the performance of CMIP6 and CMIP5 models for precipitation simulation over the HRB and found that the CMIP6-MME exhibits superior performance to the CMIP5-MME (D. Wang, Liu, et al., 2022).Therefore, the GCMs from CMIP6 were selected to project future precipitation over the HRB in this study.Since the availability of data at the beginning of our study, the daily precipitation from six CMIP6 models was obtained.The names, horizontal resolutions, respective modeling centers, countries, and references of the GCMs used in this study are listed in Table 1.The historical precipitation simulations coving the period of 1970-2014 and future projections during the period of 2015-2100 in the CMIP6 archive (https://esgf-node.llnl.gov/projects/cmip6/) were used in this study.For consistency, only the first realization (r1i1p1f1) from each model was used.More detailed information about these GCMs can be found at https://esgf-node.llnl.gov/projects/input4mips/.The future scenario in CMIP6 is a combination of different representative concentration pathways (RCPs) and shared socioeconomic pathways (SSPs) (H.Wang, Wang, et al., 2022).Three integrated scenarios, SSP1-2.6 (SSP1 + RCP2.6),SSP2-4.5 (SSP2 + RCP4.5), and SSP5-8.5 (SSP5 + RCP8.5), were used in this study.SSP1-2.6 is a low-emission scenario with a stable radiative forcing of 2.6 W/m 2 in 2100, representing sustainable development.SSP2-4.5 is a medium-emission scenario with a stable radiative forcing of 4.5 W/m 2 in 2100, representing moderate development.SSP5-8.5 is a high-emission scenario with an anthropogenic radiative forcing of 8.5 W/m 2 in 2100, representing conventional development (Eyring et al., 2016;O'Neill et al., 2013;Riahi et al., 2017).A bilinear interpolation method was utilized to make the precipitation from the selected models comparable to the observations above by placing them on a regular grid of 0.25 Â 0.25 .It has been widely accepted that the bilinear interpolation method is the most efficient one for resampling grid datasets (Liu et al., 2018;Yue et al., 2021;H. Zhu et al., 2020).Bilinear interpolation may affect the results due to its effect on data quality (Rajulapati et al., 2021).Without a standard spatial resolution, however, it would be impossible to assess and compare the precipitation simulation capabilities of the GCMs against the data that are available for observation.
F I G U R E 2 Location of the study region and the distribution of the analyzed grid points.

| Bias correction
In the bias correction method, the relationship between the climatic variables simulated by GCMs or RCMs and the corresponding observed climatic variables is employed to correct the biases of simulated climatic variables.The method defines the difference in statistical parameters between observation data and historical simulated data as model bias, which must be eliminated when simulating climate variables in the future.The bias correction method was originally developed for the climate variables simulated by RCMs.With its development, it has been gradually applied to downscaling and bias correction of climate variables simulated by GCMs in recent years (Enayati et al., 2021;Patel et al., 2021;Supharatid et al., 2021;Yue et al., 2021).Among the various bias correction methods, the QM technique based on probability distribution is one of the most widely used methods (Cannon, 2017;Enayati et al., 2021;Han et al., 2018;Maraun, 2013;Murdock et al., 2015).In the QM method, the frequency distributions of the observational climate variables and the GCM simulation variables are assumed to be consistent.In this method, a transfer function is established by fitting the cumulative distribution functions (CDFs) of the simulated and observed data in the historical period to correct the simulated data in the future period (H.Li et al., 2010).The emphasis of this approach is the establishment of the transfer function.There are three types of approaches based on their different transfer functions: distributionderived transformations, parametric transformations, and nonparametric transformations.The nonparametric transformations were recommended for most bias correction applications (Gudmundsson et al., 2012).A nonparametric transformation method called empirical quantile mapping (EQM) technique was applied to remove the systematic precipitation biases of GCM simulations over the HRB in this study, and it is widely used to correct the biases in precipitation simulated by GCM (Mishra et al., 2020;Yue et al., 2021).In the EQM method, instead of assuming specific or parametric distributions, the quantiles of observed and simulated empirical CDFs are matched for each month at each grid.The EQM method can be mathematically expressed as Equation ( 1) (Patel et al., 2021;Yue et al., 2021).
where x m-p,pre represents the pre-corrected GCM data in a projected period, and x m-p,cor is the corresponding postcorrected result; F m-h and F o-h denote the empirical CDFs of the GCM data and observed data in a historical period, respectively, and

| Multi-model ensemble
After bias correction, six methods were used to create ME of precipitation simulated by all the GCMs.These methods were MME and four ML algorithms, including RF, ET, KNN, and GBDT.Additionally, an ME method using a fusion of the above four ML algorithms was developed to obtain ME of GCMs precipitation simulations.All of these approaches attempt to reduce the uncertainties of GCMs simulations and improve the accuracy of future precipitation simulations.All these ML algorithms were implemented in Python 3.7 using the scikit-learn library.A separate calculation was used for each grid point.The following provides a brief overview of each ML algorithm used in this study.
The RF algorithm proposed by Breiman (2001) is a nonparametric ML algorithm with supervision based on decision trees, which is widely used for classification and regression.Based on the characteristics of nonparametric statistical regression and randomness, the RF model can construct multiple independent decision trees using limited samples, thus improving the accuracy of the final ensemble model and preventing overfitting (Breiman, 2001;L. Xu et al., 2019).The principle of RF is shown in Figure 3a.In the first step, the bootstrap sampling method is used to obtain n sub-training sets with the same sample size from the original training set.In the next step, a random forest is formed from n decision trees, which are generated from n sub-training sets by randomly selecting a certain number of features.Finally, the mean output of each decision tree is taken as the final result for n test samples (L.Xu et al., 2019).The RF algorithm utilizes a random feature subspace and a bagging algorithm to enhance generalization.More details can be found in Breiman (2001) and Niu et al. (2022).The ET algorithm proposed by Geurts et al. (2006) is also a nonparametric ML algorithm with supervision based on decision trees, similar to RF.The ET is an abbreviation for extremely randomized trees.As a variation of RF, this algorithm adds a level of randomness to the splitting of trees (Geurts et al., 2006).In comparison with RF, ET has two major differences (Figure 3b): (1) rather than using the bootstrap method for sampling, ET trains each decision tree with the whole original training data; (2) the RF always chooses the best possible cut point based on the principles of Gini coefficient, mean square error, etc., but ET is more aggressive to select random cut points (Jose et al., 2022).As ET selects the cut points of features randomly rather than determining the optimal position, the decision trees generated by ET will generally have a larger scale than that generated by RF, which results in the variance of the ET being further reduced relative to RF, but the bias is further increased.It has been reported that ET generalizes more effectively than RF, thereby preventing overfitting.Since the calculation process of ET is similar to that of RF, it will not be repeated here.Detailed information about ET can be found in Geurts et al. (2006) and R. Xu et al. (2020).
The KNN algorithm proposed by Cover and Hart (1967) is a common algorithm in the ML field, which is also widely used for classification and regression problems.It is a nonparametric ML algorithm based on measuring the distance between different characteristics (Manocha & Girolami, 2007).The KNN uses distance measurements to find K neighbors closest to the sample data and uses these neighbors for new data predictions (Figure 3c).When the KNN is used in regression prediction, the result is typically given by the average of K-nearest values, which is the characteristic of the prediction.It is important to choose the number of neighbors (K) carefully, as it will have an impact on the results of KNN.It has been shown in previous research that when K is greater than 5, the performance of the KNN model does not improve.Therefore, the value of K in this study was set to 5 as a result (Ahmed et al., 2020;S. Zhang et al., 2018).The formula and other information on the KNN algorithm can be found in Cover and Hart (1967) and S. Zhang et al. (2018).
The GBDT algorithm proposed by Friedman ( 2001) is an iterative algorithm based on decision trees, which has been widely used in data analysis and prediction (Ke et al., 2017).This algorithm takes the decision trees as the base learner, and gets the results by continuous iteration.Each iteration produces a result in one decision tree, and the next decision tree is based on the residual error of the previous one (Figure 3d).When all decision trees have been trained, the final predictions can be generated (Ke et al., 2017;Rao et al., 2019).GBDT follows the direction of a negative gradient to make the algorithm converge globally.The advantages of GBDT can be summarized into three points: (1) just the number, depth, and learning rate of the decision trees are the parameters that need to be considered, not complex hyperparameters; (2) the GBDT algorithm can easily produce high-accuracy predictions without having to test many parameters; (3) in this algorithm, strong loss functions are allowed, which are robust to abnormal values, and the learning efficiency is excellent (Niu et al., 2022).In recent years, the GBDT algorithm has been widely used in the field of the power industry and communication technology.The technique was introduced to develop the ME of GCMs simulations in this study.The detailed rationale of this ML algorithm can be found in Friedman (2001) and Rao et al. (2019).
For establishing a powerful ML model to develop the ME of GCMs simulations, an ensemble ML algorithm named voting regressor (VR) was utilized to fuse RF, ET, KNN, and GBDT algorithms.To the best of our knowledge, this is the first time an ME method is used to develop VR model by fusing multiple ML algorithms.VR proposed by Winkler and Makridakis (1983) is an ensemble meta-estimator that comprises several machine learning models and averages their individual predictions across all models to form a final prediction (Rajesh et al., 2022).This method is useful for a set of wellperforming models to compensate for their individual weaknesses to build a single model that can better generalize, and it has been widely used in hydrology research (Arsenault et al., 2015;Rajesh et al., 2022;Troin et al., 2022).Different individual ML models are trained on the same dataset and generate independent predictions with different learning capabilities (Yulisa et al., 2022).However, by voting on predictions from individual ML models, VR can aggregate their predictions and produce a less error-prone final prediction (An & Jiang, 2010;Geron, 2019;Yulisa et al., 2022).The advantage of the VR algorithm is that it can outperform most individual ML models because it is compatible with many individual algorithms.The principle of the VR algorithm can be seen in Figure 3e.Firstly, the RF, ET, KNN, and GBDT models individually were trained in terms of tuning parameter and hyper-parameter, to improve their performance.Wellperformed trained models were selected to input into the VR model.The contribution of each feature variable to prediction was determined using the GBDT algorithm since VR has no attribute for feature importance (Yulisa et al., 2022).There are four different types of VR algorithms based on the voting mechanism, namely absolute majority voting, relative majority voting, weighted voting, and average probability voting (Wijaya & Afianti, 2021;Z.-H. Zhou, 2012).The final output is predicted based on the average probability of voting in this study.

| Performance evaluation
In this study, the Taylor diagram, widely used to assess climate impact (Rivera & Arnould, 2020;D. Wang, Liu, et al., 2022;Yue et al., 2021), was used to assess and compare the performance of bias-corrected results and different ME development methods.The Taylor diagram is a polar coordinate diagram that shows the SD, the rootmean-square difference (RMSD), and the correlation coefficient (R), providing a comprehensive and intuitive summary of how well GCMs and observations match (Taylor, 2001).The mathematical representation of the Taylor diagram can be found in D. Wang, Liu, et al. (2022).Normalization of the SD and RMSD of GCMs data was performed using the SD of observed data in this study.When the normalized SD is close to 1, the R is equal to 1, and the normalized RMSD is close to 0, the model simulation is optimal.Furthermore, Taylor skill score (TSS) along with scatter plots were used to investigate the agreement between the precipitation simulations based on different ME algorithms and observations.According to Taylor (2001), the TSS is a measure of SD and R value in Taylor diagrams and is widely used to rank GCM simulations of precipitation performance (L.Chen & Frauenfeld, 2014;Rivera & Arnould, 2020).The mathematical representation of the TSS can be found in Ahmed et al. (2020).When the TSS is close to 1, the model simulation is optimal.
In this study, the observation and simulation precipitation datasets were divided into a calibration period and a validation period.The first 35 years  of the historical data were used for calibrating the bias correction results and the MEs developed with multiple methods.The rest 10 years (2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014) of the historical data were used for validation.The performance of monthly precipitation for the bias correction results and MEs in the validation period of 2005-2014 was evaluated in this study.

| Projected precipitation changes analysis
Based on historical and future precipitation datasets from the ME of CMIP6 GCMs bias correction results, projected future precipitation changes under three different scenarios over HRB are investigated.Percentage changes were used to express the projected precipitation changes (Equation ( 2)).The period 1995-2014 was used as the base period in the IPCC Sixth Assessment Report (AR6) on climate change projection (IPCC, 2021; T. Zhou et al., 2021).To provide a more comprehensive relative and statistical analysis consistent with the IPCC AR6 report, the period of 1995-2014 was also considered as the base period in this paper.Under the three scenarios, the projection periods of 2023-2042, 2043-2062, and 2081-2100 were deemed 'Near-term', 'Mid-term', and 'Long-term', respectively.
where PC represents the percentage changes in GCM outputs in the projection periods compared with the base period; P GCM is the precipitation of GCM outputs in the projection periods; P o is the precipitation in the base period.

| Performance evaluation of the biascorrected results and MEs
The performance of the bias-corrected results was first evaluated using the Taylor diagram and Taylor skill score in this study.Monthly precipitation of the individual GCMs before and after bias correction was evaluated against the observations in the validation period of 2005-2014.Figure 4 shows the Taylor diagram for precipitation of the individual GCMs before bias correction, after bias correction, and the MEs.As it depicts the R values between the corrected models and the observations are larger than that between the uncorrected models and the observations.The R values of all uncorrected models are below 0.6 except BCC-CSM2-MR, while that of all corrected models are above 0.65.The normalized SDs range vastly from 0.86 to 1.29 before bias correction, while that of the bias-corrected precipitation mainly ranges from 0.95 to 1.12, which reveals that the SDs of precipitation are closer to the observed SD after bias correction.As for the RMSD (normalized), the range of the uncorrected models is from 0.79 to 1.16, while that of the corrected models is from 0.73 to 0.87, which can be seen that the values are getting smaller after bias correction.The TSS for precipitation of the individual GCMs before bias correction, after bias correction, and the MEs results can be seen in Figure 5.The TSS range of all uncorrected models is from 0.294 to 0.425, while that of all corrected models is from 0.477 to 0.577, which indicates that the TSS of precipitation is closer to 1 after bias correction.Furthermore, Figures 4 and 5 describe that the performance of BCC-CSM2-MR is the best whether before or after bias correction compared with other models used in this study.In summary, the performances of these models become better after bias correction.
From Figures 4 and 5, it also can be seen that the performance is relatively limited in adjusting precipitation using the bias correction method in the HRB, which is attributed to the spatial-temporal variability and uncertainty of precipitation.However, using the ME of several GCMs can effectively reduce the biases and uncertainties of GCMs, which are present in Figures 4 and 5.It shows that the performances of all the MEs developed using different algorithms are better than individual models for precipitation.In detail, the TSSs of all the MEs are above 0.6, which are closer to 1 than that of individual models.The R values between the different MEs and observations are all above 0.80, which are larger than that of individual models.The RMSD (normalized) of different MEs in precipitation are all smaller than that of individual models, while the normalized SDs of ME are close to that of individual models.In conclusion, the above results indicate that the performances of MEs developed using different algorithms are superior to that of individual models.Compared with the individual model, the MEs show certain stabilities and can reduce the uncertainties of precipitation simulations.For accurately projecting the future precipitation changes in the HRB based on the SSP-RCP scenarios from CMIP6, the performances of different MEs are compared in the following section.

| Performance comparisons of MEs based on different algorithms
The MEs of six CMIP6 models for precipitation developed using MME, RF, ET, GBDT, and fusion of four ML algorithms were compared for selecting the best ME method to project the future precipitation changes in the HRB.The Taylor diagram on monthly precipitation during the validation period of 2005-2014 given in Figure 4 indicates that the ML approaches have improved the performance of MEs when compared with the MME approach.In detail, the ML approaches can improve the R values between the MEs and observation from 0.80 to 0.88, reduce the RMSD (normalized) values from 0.59 to 0.49, and bring the SD (normalized) closer to that of the observation from 0.85 to 0.94 compared with the MME approach.As to the performance inter-comparisons of MEs based on different ML algorithms, the fusion of RF, ET, KNN, and GBDT algorithms shows the best performances for R and RMSD, and as to the SD, ET is the best algorithm.However, Figure 5 shows that the performances of MEs developed with all ML algorithms are very close, making it difficult to distinguish which is best.
To further understand the agreement between precipitation simulated by the MEs developed with ML algorithms and observed precipitation, a comparison was conducted using scatter plots along with the Taylor skill score, as shown in Figure 6.The visual inspection of Figure 6 clearly shows that MEs developed using different ML algorithms performed much better than the MME.In other words, MEs based on ML algorithms can simulate the pattern of observed precipitation better and shows higher TSS values.The performance of all the ME approaches is more or less the same with TSS values in the range of 0.716-0.764.The ME based on the fusion of four ML algorithms refers to the highest TSS value of 0.764 during the validation period of 2005-2014.In addition, Figure 6 also shows that MEs based on all ML algorithms behave poorly in simulating extremes.The poor representation of extreme precipitation values by MEs based on ML algorithms is also reported in many previous studies (Ahmed et al., 2020;Jose et al., 2022).These results above demonstrate that ME developed using the fusion based on RF, ET, KNN, and GBDT ML algorithms matches better with the observed data than MEs developed using the MME and individual ML algorithms.
Figure 7 presents the spatial distribution of the mean annual precipitation for observation, ME developed using the fusion of four ML algorithms, and MME during the validation period of 2005-2014 over the HRB.It is apparent that both MEs results can better capture the spatial pattern of observed precipitation.Precipitation simulated by MME in the southwest region of the HRB was underestimated, while this bias in ME developed using the fusion of four ML algorithms was significantly reduced.Therefore, the mean annual precipitation of ME developed using the fusion of four ML algorithms is basically consistent with the observed value in spatial distribution, and the bias between the two is small.Additionally, the mean monthly, seasonal, and annual precipitation of observation, MME, and ME developed using the fusion of four ML algorithms were compared during the validation period of 2005-2014 over the HRB (Figure 8).As Figure 8a shows that the average monthly precipitation of ME developed using the fusion of four ML algorithms in all months is closer to the observed data except for in March, April, and October compared with that of the MME.From the shaded area shown in Figure 8a, it can be found that the uncertainties in July and August are relatively larger than in other months.
From Figure 8b, it can be seen that the mean annual precipitation of observation is 990 mm, the ME developed using the fusion of four ML algorithms is 983 mm, and the MME is 938 mm during the validation period of 2005-2014, which reveals that the precipitation simulated by ME developed using the fusion of four ML algorithms is more correctly than MME at the annual scale.At a seasonal scale, precipitation from ME developed using a fusion of four ML algorithms in spring, summer, autumn, and winter is closer to observed data than MME for all seasons, which indicates that the precipitation simulated by ME developed using the fusion of four ML algorithms is more correct than MME at the seasonal scale.These results indicate that the precipitation of ME developed using the fusion of four ML algorithms performs better in temporal variation characteristics compared with MME.
In summary, all the results above indicated that the biases in CMIP6 GCMs were substantially reduced by the EQM method and ME developed with the fusion of four ML algorithms, and the results can be used for the relevant climate change analysis and climate impact assessment.Therefore, the ME developed using the fusion of four ML algorithms was employed here to explore the future changes in precipitation over the HRB during 2023-2100.
Figure 10 shows the intra-annual distribution of projected precipitation changes over the HRB in the Nearterm, Mid-term, and Long-term under SSP1-2.6,SSP2-4.5, and SSP5-8.5 scenarios relative to the baseline period of 1995-2014.In spring, the precipitation over the HRB shows an increasing trend relative to the baseline period of 1995-2014, with a maximum increase of 98.04% in the Long-term under SSP5-8.5.In summer, the precipitation over the HRB shows an increasing trend, except for a decreasing trend in the Near-term (SSP1-2.6 and SSP2-4.5).In autumn, the precipitation over the HRB shows an increasing trend, except for a decreasing trend in the Near-term under SSP2-4.5, with a decrease of 11.88%.In winter, the precipitation over the HRB shows a decreasing trend, except for a slightly increasing trend in the Long-term (SSP5-8.5 and SSP2-4.5)and Mid-term (SSP1-2.6).As to the monthly precipitation change, the variation law of precipitation at the monthly scale is consistent with that at the seasonal scale in different scenarios and different periods.The projected precipitation shows an increasing trend from January to June and October under all the SSP-RCP scenarios relative to the baseline period of 1995-2014, with a maximum increase of 121.53% in the Long-term under SSP5-8.5.In general, the projected precipitation changes at seasonal and monthly scales show an increase over time with the largest in the Long-term, followed by the Mid-term, with the smallest change in the near term.

| Spatial distribution of projected precipitation changes
The spatial distribution of the projected changes in the mean annual precipitation over the HRB during the three different future periods (2023-2042, 2043-2062, and 2081-2100) relative to the base period of 1995-2014, under the three scenarios (SSP1-2.6,SSP2-4.5, and SSP5-8.5)were explored in this study, as shown in Figure 11.As clearly shown in the figure, the spatial distribution patterns of annual precipitation over the HRB during the three different future periods under all the SSP-RCP scenarios are relatively consistent, with larger precipitation in the northeastern and smaller in the southeastern.The precipitation in the Long-term under all the SSP-RCP scenarios shows an increasing trend relative to the baseline period of 1995-2014.In detail, the ranges of the projected precipitation changes in the Near-term under respectively.In the Mid-term, the ranges of the projected precipitation changes under respectively.In the Long-term, the ranges of the projected precipitation changes under SSP1-2.6,SSP2-4.5, and SSP5-8.5 are 4.45%-24.89%,8.70%-34.86%,and 30.41%-90.59%, respectively.From the above, it can be found that the precipitation in the Near-term and Mid-term is the largest under SSP5-8.5,followed by SSP1-2.6, and the least under SSP2-4.5, while the precipitation in the Long-term is the largest under SSP5-8.5,followed by SSP2-4.5, and the smallest under SSP1-2.6.In addition, it can also be found that under the three SSP-RCP scenarios, the precipitation in the Long-term is the largest, followed by Mid-term, with the smallest in the Near-term.

| DISCUSSION
There still exist large systematic biases in the GCMs from CMIP6 (Figure 4).It is impossible to study the impact of climate change on hydrological processes, especially in regional areas (S.Zhang & Chen, 2021).For improving the accuracies of the CMIP6 GCMs outputs, the EQM method was used to correct the biases.The evaluation of the bias-corrected results shows that the performances of the six GCMs have been greatly improved after bias correction compared with before bias correction (Figure 4).The bias correction results of precipitation in this study are comparable to numerous previous studies (Ahmed et al., 2020;Anjum et al., 2019;R. Chen et al., 2022;Supharatid et al., 2021).However, it also can be seen that the performance is relatively limited in adjusting precipitation using the EQM method in the HRB due to the spatial-temporal variability and uncertainty of precipitation, as shown in Figure 4.This result has also been reported in previous studies (Almazroui et al., 2020;R. Chen et al., 2022;Yue et al., 2021).In order to further reduce the biases and uncertainties of GCMs outputs, the MEs are often used to project future climate change.For example, Yue et al. (2021) investigated future climate changes using the MME of the bias-corrected datasets, because they found that the MME performs better than individual models for each climatic variable.R. Chen et al. (2022) found that the MME of the five best models was in stronger agreement with observations than the single model, and used the MME under four SSP-RCP scenarios to project the seasonal precipitation changes in 2015-2099 over the Tibetan Plateau.
In recent years, with the rapid development of AI algorithms, ML algorithms have begun to be used to develop MEs and have been proven to be more effective compared with MME (Ahmed et al., 2020;Jose et al., 2022).There still exists a certain bias in the precipitation from GCMs after the bias correction using the EQM method.Similarly, after bias correction, to further reduce the biases and uncertainties of GCMs, the same approach was adopted in this study.Different from the previous studies, this study developed an ME approach using VR algorithm to fuse RF, ET, KNN, and GBDT algorithms, and compared with the ME using the MME, RF, ET, KNN, and GBDT methods.The results show that the ME method developed with the fusion of RF, ET, KNN, and GBDT algorithms proposed in this study performs better than other methods (Figures 5 and 6).In addition, Figures 7 and 8 present the comparison results of the observation, MME, and ME developed with the fusion of four ML algorithms.Results indicate that the ME based on the fusion of four ML algorithms performs better than MME.The VR algorithm used in this study to fuse multiple ML algorithms for developing ME can aggregate the predictions from several individual ML models by voting and produces a final prediction that is less prone to error (An & Jiang, 2010;Geron, 2019;Yulisa et al., 2022).The VR algorithm can provide higher accuracy and outperform most individual ML models because it is compatible with many individual algorithms.Therefore, this method can be extended to other basins for ME development.
Based on the ME developed with the fusion of four ML algorithms, the differences in the projected precipitation between the data in the baseline period (1995)(1996)(1997)(1998)(1999)(2000)(2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014) and the three future periods (2023-2042, 2043-2062, and 2081-2100) under SSP1-2.6,SSP2-4.5, and SSP5-8.5 scenarios were explored in this study .Results indicate that projected precipitation changes show an increase significantly over the HRB during the period 2023-2100, with most increase under SSP5-8.5,followed by SSP2-4.5, and the least increase under SSP1-2.6 by the end of the 21st century.The results indicate that the future climate of the HRB will show a significant trend of humidification.With the continuous humidification of the climate, the frequency of abnormal precipitation events in the HRB will increase in the future, mainly in the form of a sudden increase in precipitation.It can also be found that the projected precipitation increases with time over the HRB, with the largest in the Longterm, followed by Mid-term, and the smallest in the Near-term.Similar results can also be found in the previous studies in the Yangtze River Basin, to which the HRB belongs (X.Li et al., 2022;Yue et al., 2021).On the one hand, with the precipitation increase in the future, the water resources of underground, farmland, and reservoirs can be effectively supplemented over the HRB.In particular, there are many water diversion projects in the HRB, and the precipitation increase can supplement the continuously increasing volume of transferred water.On the other hand, the precipitation increase may increase the risk of urban flooding and waterlogging over the HRB.The results displayed that the northeastern areas of the HRB are projected to experience a large precipitation increase in the future, which needs to be paid attention to for disaster prevention and mitigation in these areas.In addition, it is worth noting that the projected precipitation change in winter is negative relative to the base period, which indicates that the HRB may experience droughts in winter in the future.
Uncertainties in GCM simulations and future climate change projections are often a limiting factor, which is unavoidable in the research field of climate change, especially in regional areas (Yue et al., 2021).Typically, the sources of uncertainties in future precipitation projection mainly include the selection of GCMs, natural variability within the selected GCMs, and the SSP-RCP scenarios.Figure 9 shows that the uncertainties of precipitation variable increase over time under the three SSP-RCP scenarios.Additionally, the uncertainty ranges of the average ME for precipitation variable are smaller than that of the six biascorrected model projections.Similar results can also be found in the previous studies over the Yangtze River Basin (X.Li et al., 2022;Yue et al., 2021).However, it should be noted that uncertainties still exist in the projected precipitation changes, and efforts are still needed to further reduce the uncertainties.Therefore, future climate change projections should be interpreted with caution.These findings revealed the projections and uncertainties of CMIP6 precipitation, which will be helpful for a better understanding of the future evolution of precipitation in the HRB.
In recent years, global climate change has impacted the occurrence of extreme precipitation and shifts in precipitation patterns, which caused severe flooding disasters and other extreme weather (D.Wang, Liu, et al., 2022).The intensity-duration-frequency (IDF) curves can represent the relationship between extreme rainfall intensity, duration, and return period, which are often used for flood risk analysis and various urban infrastructure designs (Crévolin et al., 2023).IDF curves are usually constructed based on the historical characteristics of rainfall and assuming that it would remain constant in the future.However, characteristics of extreme rainfall would change under global climate change.Therefore, there is a need to better understand future changes in extreme precipitation values and to update the methodology for constructing IDF curves.Under changing environment, the covariate-based nonstationary method and the GCM-based method are widely applied to construct future IDF curves (Crévolin et al., 2023;Rootzén & Katz, 2013;Yan et al., 2021;Yan et al., 2022).Among them, the GCM-based method can fully consider the future physical processes in the atmosphere, ocean, cryosphere, and land surface (Crévolin et al., 2023).However, the GCM-based methods highly rely on the projections of local high-resolution extreme precipitation, which inevitably contains different degrees of uncertainties, especially for regional scales (Yan et al., 2021).Therefore, more efforts need to be made to update the future IDF curves constructed for analyzing flooding risk under the changing environment.
In addition, there are some limitations in this study.For example, only six GCMs from CMIP6 were used in this study.More GCMs from the CMIP6 should be considered in future studies.The proposed method for developing ME of GCMs simulations can be strengthened and extended from different aspects.For example, the mechanisms of different ML algorithms in the current combinatorial algorithms were not considered in this study.A certain ML algorithm may play a negative impact on the ME results.Therefore, ablation experiments need to be supplemented to quantify the contribution of each ML algorithm in the current combinatorial algorithms in the future, in order to find the best combination.Moreover, only the RF, ET, KNN, and GBDT algorithms were fused using VR algorithms for ME of GCM simulations in this study.However, there are numerous ML algorithms.More ML algorithms can be fused to develop ME.Due to the strong ability of ML algorithms to mine data features, more ME methods based on ML algorithms need to be developed to provide more reliable climate projections.Recently, some ME methods based on other advanced ML algorithms such as LSTM (Jose et al., 2022) and ANN (Ahmed et al., 2020) have been developed and successfully applied in the field of GCMs ensemble.

| CONCLUSIONS
In this study, we explored the future precipitation changes over the Hanjiang River Basin (HRB) in the 21st century, relative to the baseline period of 1995-2014.The major conclusions can be summarized as follows: 1.The comparison of multimodel ensembles (MEs) developed using six methods shows that the ME proposed in this study developed by fusing random forest, extra tree, K-nearest neighbors, and gradient boosting decision tree machine learning (ML) algorithms performed better than the mean and individual ML approaches.The Taylor skill score value of the ME based on the fusion of four ML algorithms refers to the highest with 0.764 during the validation period of 2005-2014, compared with other methods.Therefore, the fusion of four ML algorithms is used for the creation of ME to project the future precipitation changes over the HRB in this study.2. Projections display a significant upward trend of the annual precipitation over the HRB during 2023-2100 under the three scenarios.The increase in the projected precipitation over the HRB is the largest under SSP5-8.5 (À9.27% to 112.84%), followed by SSP2-4.5 (À30.48% to 44.67%), and the smallest under SSP1-2.6 (À37.19% to 37.78%).In addition, the increase in the projected precipitation over the HRB is the largest in the Long-term, followed by Mid-term, with the smallest in the Near-term.3. Uncertainties still exist in the projected precipitation over the HRB.Although uncertainties were reduced to some extent using the ME based on the ML algorithms in this study, efforts should also be made to further reduce the uncertainties.More research studies on reducing uncertainties in global climate models are critical to better understand climate change in the future.

F
I G U R E 4 Taylor diagram for precipitation of the individual global climate models before bias correction, after bias correction, and the multimodel ensembles at the monthly scale during the validation period of 2005-2014.The 'model_cor' refers to the model after bias correction.For example, the 'BCC-CSM2-MR_cor' refers to the BCC-CSM2-MR after bias correction.ET, extra tree; GB, gradient boosting; KN, K-nearest neighbors; RF, random forest; RMSD, root-mean-square difference.F I G U R E 5 Taylor skill score for precipitation of the individual global climate models before bias correction, after bias correction, and the multimodel ensembles at the monthly scale during the validation period of 2005-2014.ET, extra tree; GB, gradient boosting; KN, K-nearest neighbors; RF, random forest; TSS, Taylor skill score.

F
I G U R E 6 Scatter plot of observation and multimodel ensembles simulated monthly precipitation during the validation period of 2005-2014.ET, extra tree; GB, gradient boosting; KNN, K-nearest neighbors; RF, random forest.

F
I G U R E 7 Spatial distribution of the mean annual precipitation for (a) observation, (b) multimodel ensemble (ME) developed using the fusion of multiple algorithms (ME_Fusion), and (c) MME (ME_Mean) during the validation period of 2005-2014 over the Hanjiang River Basin.F I G U R E 8 The mean (a) monthly, (b) seasonal, and annual precipitation of observation, MME (ME_Mean), and multimodel ensemble (ME) developed using the fusion of multiple algorithms (ME_Fusion) during the validation period of 2005-2014 over the Hanjiang River Basin.The shaded area in (a) represents the uncertainties (the value range of the models) of the six CMIP6 models.

F
I G U R E 9 Annual mean precipitation changes (PC) (%) over the Hanjiang River Basin in SSP1-2.6,SSP2-4.5, and SSP5-8.5 from the multimodel ensemble (developed with the fusion of random forest, extra tree, K-nearest neighbors, and gradient boosting decision tree machine learning algorithms) relative to the base period of 1995-2014.The shaded areas in (a) represent the uncertainties (the value range of the models) of the six CMIP6 models in SSP1-2.6,SSP2-4.5, and SSP5-8.5.

F
I G U R E 1 0 Projected percentage changes of precipitation (PC) over the Hanjiang River Basin at different time scales in Near-term, Mid-term, and Long-term under SSP1-2.6,SSP2-4.5, and SSP5-8.5 from the multimodel ensemble (developed with the fusion of random forest, extra tree, K-nearest neighbors, and gradient boosting decision tree machine learning algorithms) relative to the base period of 1995-2014.F I G U R E 1 1 Spatial distribution of projected percentage changes (PC) for the mean annual precipitation over the Hanjiang River Basin from the multimodel ensemble (developed with the fusion of random forest, extra tree, K-nearest neighbors, and gradient boosting decision tre machine learning algorithms) in the (a-c) Near-term, (d-f) Mid-term, and (g-i) Long-term, relative to the base period of 1995-2014.
List of the CMIP6 models analyzed in this study.
T A B L E 1