Impacts of bias nonstationarity of climate model outputs on hydrological simulations

The assessment of climate change impacts on the hydrological cycle has been widely investigated during recent years (Graham et al. 2007; Li et al. 2015; Chen et al. 2016a; Marhaento et al. 2017; Shen et al. 2018; Lu & Qin 2020; Ragettli et al. 2020). Global climate models (GCMs) can provide climate variables (e.g. precipitation and temperature) for the future period used as inputs of hydrological models for hydrological impact studies. However, the coarse resolution of GCMs does not meet the need of high resolution for the hydrological models (Maraun et al. 2010). In parallel, GCMs are imperfect representations of reality with systematic biases found between climate model simulations and observations. To resolve these problems, a number of downscaling methods have been developed during the last two decades. Especially, bias correction becomes a standard procedure when using regional climate model (RCM) outputs for hydrological impact studies (van Pelt et al. 2009; Teutschbein & Seibert 2012; Olsson et al. 2015; Zhuan et al. 2019). During recent years, bias correction methods are also commonly used for GCM outputs. Since the resolution of climate model outputs is lower than that of observations, bias correction methods also act as downscaling methods. The usually used bias correction methods range from simple mean-based scaling to sophisticated distribution-based mapping and multivariate or/and multisite correction.

The commonly used bias correction methods are usually based on an assumption that biases of climate models outputs are stationary. In other words, these bias correction methods assume climate model outputs present the same biases in magnitude and direction between historical and future periods. However, this assumption may not be valid, as pointed out in a few recent studies (e.g. Buser et al. 2009; Ehret et al. 2012; Gutierrez et al. 2013; Maurer et al. 2013; Velázquez et al. 2015; Dixon et al. 2016). For example, Christensen et al. (2008) evaluated bias nonstationarity on simulated monthly precipitation and temperature from an ensemble of 13 RCMs over Europe. They found that biases of simulated precipitation and temperature vary as a second-order function of observations, suggesting that bias nonstationarity exists in different climatic regimes. Maraun (2012) verified the bias stationary assumption in RCMs for European seasonal mean temperature and precipitation sums using a pseudo-reality approach, which considers one climate model as pseudo observation to compare with other climate model simulations. The results showed that biases are relatively stable in general, but bias nonstationarity was identified in some regions where changes in potentially relevant physical variables are significant. This study was conducted in the climate model world, the transferability to real world needs to be further investigated. To test the bias nonstationarity in real-world climate, Chen et al. (2015) compared the biases between climate model simulations and corresponding observations over two historical periods. The results showed that biases of climate models simulated precipitation are nonstationary even for two close historical periods, while temperature biases are relatively stationary. This study attributed the bias nonstationarity of precipitation to natural climate variability (multi-decadal climate variability). A similar study was also carried out by Wang et al. (2017), who tested bias nonstationarity of precipitation in the eastern United States determined by a skill score, which compares the errors of a downscaling method over the validation period with the errors of observations between calibration and validation periods. The results show that precipitation biases are nonstationary at most of the stations, especially for the annual extreme precipitation. Taking into account multi-decadal climate variability, Nahar et al. (2017) investigated the bias stationarity of six GCMs simulated monthly and seasonal precipitation and temperature over multi-decadal time periods from 1900 to 1999 in Australia. One hundred samples of the observations and GCM simulations generated using bootstrapping was used to calculate uncertainty in biases with a 95% confidence, which represents a possible range of bias if bias stationarity exists. When the actual bias goes beyond this uncertainty, the bias was assumed to be nonstationary. This study showed that biases of precipitation and temperature are not stationary for some regions in Australia (e.g. east coast of Australia), because of the natural climate variability.

For a historical period with no anthropogenic forcing, the bias nonstationarity is attributed to natural climate variability. Natural climate variability refers to variations in the mean state and other statistics of the climate, due to natural internal processes or natural external forcing. It involves a wide range of time scales, from one day to the next, as well as from one year or multi-decade to the next. For the application of a bias correction method in climate change impact studies, the decadal and multiple decadal time scales are widely and implicitly used. At the decadal or multi-decadal time scales, the natural climate variability includes some natural modes of decadal or multi-decadal climate variability, such as El Niño/Southern Oscillation (ENSO), Atlantic Multidecadal Oscillation (AMO), Pacific Decadal Oscillation (PDO), and Interdecadal Pacific Oscillation (IPO). The phase of these modes has been found to be linked with changes in precipitation and temperature. For example, the PDO warm (cold) phase periods link to the decrease (increase) in precipitation in the majority of China (Ouyang et al. 2014). Precipitation decreases during IPO warm phase periods, while it increases during IPO cool phase periods in east Australia (Power et al. 1999). However, since the unpredictability and complexity of natural dynamic processes make it difficult for climate models to accurately capture the characteristics, GCMs from Coupled Model Intercomparison Project phase 5 (CMIP5) have been found to poorly represent the phase in natural modes of decadal or multi-decadal climate variability (Polade et al. 2013; Ruiz-Barradas et al. 2013; Bellenger et al. 2014; Fuentes-Franco et al. 2016). This poor representation for GCMs could lead to biases between simulations and observations when sampling in a finite decadal or multi-decadal time, and further lead to changes in bias among different decadal or multi-decadal periods (Nahar et al. 2017).

In climate change impact studies, natural climate variability used to be estimated based on multi-decadal climate model simulations in the absence of anthropogenic-induced climate change (Hulme et al. 1999; Arnell 2003; Chen et al. 2016b). In these methods, the differences among several separated multi-decadal periods obtained from a long time series of unforced simulations were supposed to represent the multi-decadal natural climate variability. Similarly, to evaluate bias nonstationarity of climate model outputs, the estimate of bias nonstationarity is depicted as the range of the biases among multi-decadal historical periods, while the timescale of each period may have an effect on the results of estimation. The different temporal scales of sampling can lead to changes in the magnitude and direction of bias, due to different phases between GCM outputs and observations. Therefore, the impacts of different temporal scales of natural climate variability on bias nonstationarity of climate model outputs need to be investigated, especially for their propagation in hydrology.

Even though few studies (Maraun 2012; Chen et al. 2015; Velázquez et al. 2015; Nahar et al. 2017) have investigated the bias nonstationarity of climate model outputs, the role of bias nonstationarity in future climate change was not investigated, especially in hydrological impact studies. Accordingly, the objectives of this study are to investigate: (1) the bias nonstationarity of climate model outputs in the context of natural climate variability, (2) the role of bias nonstationarity in future climate change, and (3) the propagation of bias nonstationarity in hydrological impact studies. The first objective is achieved by separating a long historical climate series (100 years) to decadal and multi-decadal time periods to calculate the difference in bias among different periods. Precipitation and temperature, which are primary variables for hydrological simulations, are investigated. The second objective is achieved by comparing the range of biases among multiple periods to climate change signals between future and historical periods. The last objective is achieved by analyzing streamflow time series simulated by a hydrological model using the above two approaches.

The rest of the paper is organized as follows: the following section introduces the study area and hydro-climate data followed by the methodology. The main results are then presented, and the discussion and conclusions are presented in the final section.

The case study was conducted at the Hanjiang River basin (Figure 1), located in south-central China, which is the largest tributary of the Yangtze River basin. The Hanjiang River flows through Shaanxi and Hubei Provinces with a drainage area of 159,000 km² and a length of 1,567 km. The watershed above the Danjiangkou Reservoir, with a sub-basin area of 95,200 km² and a length of 918 km, was used in this study. The Danjiangkou reservoir is the water source of the Middle Route of the South-to-North Water Transfer Project in China. The watershed has a subtropical monsoon climate. The mean annual precipitation is about 840 mm, of which 70–80% of the total amount falls in the wet season from May to September. The average maximum and minimum temperatures are 24–29 and 0–3 °C, respectively. The daily mean discharge of the Hanjiang River is approximately 1,150 m³/s.

This study used both observed and GCM-simulated precipitation and temperature for the Hanjiang watershed. The inflow to the Danjiangkou reservoir calculated using a water mass balance method was also used for hydrological model calibration and validation. The streamflow time series covering the 1961–2000 period was provided by the Bureau of Hydrology of the Changjiang Water Resources Commission.

The observed precipitation and temperature were obtained from the gridded Climatic Research Unit (CRU) Time-series (TS) data, produced by CRU at the University of East Anglia, UK. This dataset provides monthly gridded data with a high resolution of 0.5 × 0.5°. The latest version CRU TS version 4.01 (Harris et al. 2017) covering the period 1901–2000 was used in this study.

The GCM-simulated daily precipitation and maximum and minimum temperatures were extracted from the database of CMIP5. In order to capture the uncertainty related to climate models, 17 GCMs were used in this study (Table 1). The time period covers from 1901 to 2100. All GCMs' simulations during the 1901–2005 period were driven by historical climate forcing (including natural solar, volcanic variations and anthropogenic radiative forcing), while those during the 2006–2100 period were generated under the representative concentration pathway (RCP) 4.5. RCP4.5 corresponds to the medium anthropogenic forcing for the 21st century. The other greenhouse gas emission scenarios (RCP8.5) are also presented and discussed in the Discussion section. Since this study only used monthly and annual data, the daily precipitation and temperature were summed (for precipitation) or averaged (for temperature) to monthly and annual values. In order to run a lumped hydrological model, all gridded precipitation and temperature within and surrounding the Hanjiang watershed were averaged to a single time series using the Thiessen polygon method for both observations and GCM simulations.

A climate time series for the past–present period may contain both natural climate variability and anthropogenic climate change signal. In order to investigate the impacts of natural climate variability on bias correcting climate model outputs, the climate change signal needs to be first removed. A simple linear detrending method of Zhuan et al. (2018) was used to remove the anthropogenic climate change signal for observed and GCM-simulated monthly precipitation and temperature for the 1901–2000 period. This detrending method first detects the trend and breakpoint using the Mann–Kendall test (Mann 1945; Kendall 1975). If a trend exists, it is removed for the period after the breakpoint using a linear method. To preserve the seasonality, the detrending method was applied for each month, respectively.

After detrending, the variation of observed and GCM-simulated precipitation and temperature within the 1901–2000 period is considered to be only attributed to the natural climate variability. Since natural climate variability is inherently complex and manifests itself over various temporal and spatial scales, this study only investigates its impacts on bias correcting climate model outputs over decadal and multi-decadal temporal and watershed spatial scales. In order to investigate the impacts of natural climate variability at the temporal scale, the 100-year precipitation and temperature time series are divided into several independent decadal and multi-decadal periods. The time window consists of 10, 20, 33 and 50 years for the 100-year period. For example, when the time window is 20 years, the 100-year period is divided into five 20-year periods (1901–1920, 1921–1940, 1941–1960, 1961–1980 and 1981–2000). The first period of each time window (e.g. 1901–1910 for 10-year window) is used as the baseline period.

The differences in bias of precipitation (DBP) and temperature (DBT) between the baseline period and all other historical periods are then calculated for each GCM. The use of the difference in bias is beneficial to compare with the results from different GCMs. The differences in bias represent the possible changes of GCM bias due to natural climate variability at decadal and multi-decadal scales. The range of difference in biases is used to judge the bias stationarity of climate model outputs. The impacts of natural climate variability on bias nonstationarity is investigated by using multiple time windows. The main steps are illustrated in Figure 2 for a 20-year time window.

To identify the role of bias nonstationarity in future climate change, the signal-to-noise ratio (SNR) is used as a criterion to quantify the extent of bias nonstationarity relative to future climate change. The SNR is defined as the ratio of climate change signal (CCS) of future periods and the range of bias (RB) among multiple historical periods.

Then, SNR is calculated as the absolute ratio of CCS and RB for both monthly and annual precipitation and temperature. If the CCS is relatively smaller than the RB, e.g. the SNR being smaller than 1, the RB is considered to have large impacts on future climate changes and their hydrological impacts and vice versa.

The observed monthly streamflow was used to model calibration (1961–1980) and validation (1981–2000). Model calibration was done automatically using shuffled complex evolution-University of Arizona (SCE-UA) algorithm (Duan et al. 1992). The optimal combination of parameters was chosen on the basis of the Nash–Sutcliffe efficiency (NSE) coefficient (Nash & Sutcliffe 1970). The values of NSE are 0.79 for both calibration and validation periods for the Hanjiang watershed, indicating a good performance.

Using the calibrated hydrological model, the streamflow time series are generated for each time window using the detrended observed and GCM-simulated precipitation and temperature for the historical period. Meanwhile, the streamflow time series are also generated for each time window of the future periods. The propagation of precipitation and temperature bias nonstationarity in hydrology was investigated using the same methods as presented above under ‘Calculation of difference in bias’ and ‘Calculation of climate change signals relative to differences in bias’.

Figure 3 presents the difference in bias of mean annual precipitation between each historical period and the baseline period for four different time window lengths: 10, 20, 33 and 50 years. For each time window length, the dash lines represent differences in bias of all historical periods relative to the baseline period for all 17 GCMs, and the thick line is their median value. The x axis represents the midpoints of the time window. The results show that the difference in bias of mean annual precipitation varies with time over historical periods, indicating that the bias of precipitation is not stationary. The difference in bias reflects the impacts of natural climate variability on bias nonstationarity, because the trend of precipitation change was removed. The results in Figure 3 show that the natural climate variability has great impacts on bias nonstationarity for mean annual precipitation, especially for the short time window. However, this problem becomes less important as the time window gets longer, as indicated by the fact that the difference in bias decreases with the extension of the time window. For example, in terms of the median value over all GCMs, the difference in bias changes between –11.5 and 3.8% for the 10-year window, between 0.5 and 5.5% for the 20-year window, between –2.7 and –0.4% for the 33-year window, and it is –4.4% for the 50-year window. This indicates that the decadal variability has more impacts on bias nonstationarity than multi-decadal variability. Not only the median of difference in bias but also the uncertainty (the difference between maximum and minimum values over the whole period) related to GCMs decreases with the longer time window. The uncertainty of difference in bias is 56.4% for the 10-year window, 27% for the 20-year window, 17.5% for the 33-year window, and 15.8% for the 50-year window. Therefore, using a long reference period to filter out low frequency modes of variability may lessen the impacts of natural climate variability on bias correction methods. However, if just looking at a single climate model, the difference in bias can reach 10%, even for the 50-year time window.

Figure 4 presents the difference in bias of mean annual temperature between each historical period and the baseline period for four different time windows. The results show that the differences in temperature bias vary with time. In particular, the magnitude and uncertainty of differences in bias are larger for short than long time windows. For example, the median value of difference in temperature bias changes between –0.7 and –0.05 °C for the 10-year window, between –0.5 and –0.05 °C for the 20-year window, between –0.28 and –0.06 °C for the 33-year window, and it is –0.13 °C for the 50-year window. The negative values of the difference in bias indicate the bias over the first baseline period is larger than that over the following several historical periods. The large uncertainty related to climate models is observed, especially for the short time window. Since the temperature time series were detrended before calculating the difference in bias, the variation resulted from natural climate variability. Even though the bias varies with time due to natural climate variability, its impacts on future climate change need to be investigated by comparing future climate change signals.

In order to investigate the role of bias nonstationarity in future climate change, the SNR of mean annual precipitation and temperature is calculated for all four time windows. Figure 5 presents the SNR of mean annual precipitation and temperature over multiple time periods from 2001 to 2100 for 10-, 20-, 33- and 50-year windows. Each boxplot depicts the distribution of SNR from an ensemble of 17 GCMs. The larger SNR values suggest that the bias range is relatively less important compared with the climate change signal and vice versa. The unity of SNR means that the range of bias is equal to the climate change signal. In this case, the calculated climate change signal may just be the difference in bias between two periods. Figure 5(a) shows that SNR values of precipitation are consistently less than one for the 10-year window of all future periods. In this case, it is hard to judge the existence of climate change signal. The calculated climate change signal may be just a part of the range of bias caused by the natural climate variability, especially taking into account the fact that the SNR values are mostly constant over time. With the increase in the length of time window, the bias nonstationarity of mean annual precipitation becomes less important, as indicated by some GCMs presenting SNR values being greater than one. In particular, the uncertainty related to GCMs becomes larger with the increase in the length of time window. However, the median values of SNR are still less than one, even for the 33- and 50-year time window. The above results indicate that great attention should be paid to natural climate variability when investigating the change of precipitation, in the face of the increase in the greenhouse gas emission scenario.

However, different results are observed for annual mean temperature. Even though the difference in bias varies with time, the anthropogenic climate change is significantly greater than the bias nonstationarity for annual mean temperature, especially for the future periods (Figure 5(b)). In other words, even though the bias nonstationarity can affect the detection of the real anthropogenic climate change of annual mean temperature, the influence is limited. Comparing to future climate change signal, the temperature bias can be considered as stationary. However, for the first and second decades, the natural climate variability still plays an important role in anthropogenic climate change or the climate change signal may not exist, as indicated by the relatively small SNR. The SNR values increase with time, which is caused by the increase in the anthropogenic climate change, as the range of bias is supposed to be constant. Meanwhile, the dispersion of boxplot is gradually amplified with time, resulting from the increase in the GCM uncertainty in the far future. Additionally, the magnitude and dispersion of the SNR become larger with the extension of time window. This implies that the impacts of bias nonstationarity on future climate changes are important for the long time window.

The SNR of precipitation and temperature is also calculated at monthly scale. Figure 6 shows the portrait diagram of the median value of the SNR for mean monthly precipitation and temperature over multiple time periods from 2001 to 2100 with 10-, 20-, 33- and 50-year windows, respectively. The y-axis represents the 12 months from January to December. The median value of the SNR was calculated over 17 GCMs for each period. Generally, the bias of monthly precipitation is nonstationary, while the temperature is stationary compared to future climate change. In other words, the impacts of bias nonstationarity on future climate changes are large for monthly precipitation, but it is limited for monthly temperature. Similar to annual precipitation and temperature, the impacts of bias nonstationarity on future climate changes lessen for the longer time window, as indicated by the increase in the SNR. However, the monthly variability is observed for these impacts. In terms of the monthly precipitation (Figure 6(a)), the SNR in the wet season (May–September) is larger than in that the dry season (October–April), especially for the mid and far future, while the SNR is still smaller than one for 10-, 20- and 33-year windows. When the time window increases to 50 years, the SNR in the wet season becomes larger than 1. This indicates that the impacts of bias nonstationarity on future climate changes are less important in the wet season than in the dry season. Therefore, the impacts of bias nonstationarity in the dry season need to be paid more attention than in the wet season and annual scale. In terms of the mean monthly temperature (Figure 6(b)), the SNR from May to November is larger than other months for the mid and far future, and the maximum values of the SNR present between May and August over each decadal or multi-decadal period. In addition, the comparison of the SNR between annual and monthly temperature shows that the SNR at the monthly scale is smaller than that at the annual scale, especially for short time windows, while for the 50-year window, the SNR between June and October is comparable to or even larger than that at the annual scale. In the near future (2001–2030), the SNR in each month is comparable to or smaller than one for the short time window, because of the weak climate change signal.

The propagation of bias nonstationarity of climate models simulated precipitation and temperature in hydrology is investigated by running a hydrological model over the Hanjiang watershed. Figure 7 presents the difference in bias of mean annual streamflow between each historical period and the baseline period for 10-, 20-, 33- and 50-year time windows. The difference in bias between GCM-driven and observation-driven streamflow also varies with time, which reflects the impacts of natural climate variability on hydrological simulations. The results show that the impacts of natural climate variability are amplified from climate world to hydrological world. For example, the difference in bias for all GCMs varies between –30 and 20% for mean annual precipitation, while that varies between –60 and 60% for mean annual streamflow for the 10-year time window. This is because of the non-linear process from climate to hydrology. Similar to precipitation and temperature, the extent of the impacts of natural climate variability gradually mitigates as the time window gets longer. For example, the median value of difference in streamflow bias changes between –12.3 and 25.6% for the 10-year window, between 8.5 and 20.2% for the 20-year window, between –0.9 and 5.9% for the 33-year window, and it is –6.2% for the 50-year window. The uncertainty of difference in streamflow bias related to GCMs also reduces as the extension of the time window. The uncertainty is 121.6% for the 10-year window, 58.2% for the 20-year window, 35% for the 33-year window and 27.1% for the 50-year window.

The impacts of bias nonstationarity of climate model outputs on future streamflow changes are also evaluated using the SNR as a metric. Figure 8 shows the SNR of mean annual streamflow over multiple periods from 2001 to 2100 for 10-, 20-, 33- and 50-year windows. In the future, the magnitude of the streamflow change between the future and baseline period is much smaller than the range of bias, as indicated by the small SNR values, especially for the 10- and 20-year time windows. This indicates that the natural climate variability plays an important role in future streamflow changes, especially for the short time window. In other words, previous studies that attribute the streamflow change to anthropogenic climate change may not be correct when taking into account the impacts of climate change. Similar to precipitation and temperature, the importance of natural variability reduces with the extension of the time window, as indicated by the increase in the SNR. For example, the median values of the SNR change between 0.07 and 0.21 for the 10-year window, between 0.09 and 0.3 for the 20-year window, between 0.38 and 0.85 for the 33-year window, and between 0.61 and 1.01 for the 50-year window. The uncertainty related to GCMs becomes greater for the longer time window. With the exception of the second period of the 50-year time window, the median values of SNR are smaller than one, indicating that the natural climate variability is still more important than the anthropogenic climate change for future streamflow simulations. These results of streamflow are similar to those of mean annual precipitation, but the SNR of the former is smaller than the latter, due to the larger range of bias and lower streamflow change.

Figure 9 presents the median of SNR of mean monthly streamflow across 17 GCMs for 10-, 20-, 33- and 50-year windows. The results show that the bias nonstationarity has large impacts on the change of future monthly streamflow, even though the impacts mitigate with the extension of the time window. For example, the values of SNR are consistently smaller than one for 12 months and all future periods in 10- and 20-year windows. This is also the case for most months and future periods in the 33-year window. For the 50-year window, the values of SNR increase to larger than 1 in the wet season, indicating that the anthropogenic climate change may be more important than the bias nonstationarity resulting from natural climate variability for the wet season. In addition, the impact of bias nonstationarity on future climate change is amplified in hydrological simulations, as indicated by the fact that the SNR of streamflow is smaller than that of precipitation and temperature.

This study first investigated the bias nonstationarity of GCM-simulated precipitation and temperature in the context of natural climate variability for multiple time windows. The role of bias nonstationarity in future climate changes was estimated by comparing the range of bias in the historical periods to climate change signal in the future periods. Finally, the propagation of bias nonstationarity in hydrological simulations is investigated by running a hydrological model. The main findings are as follows: (1) the biases of GCM precipitation and temperature vary with time, due to natural climate variability, especially for the short time windows; (2) precipitation bias nonstationarity plays an important role in future precipitation changes at annual and monthly scales, while the importance of temperature bias nonstationarity in future temperature change is not significant; (3) the bias nonstationarity of climate model outputs is amplified when driving a hydrological model for hydrological simulations. In other words, the bias nonstationarity of precipitation and temperature has considerable impacts on future streamflow changes for this specific watershed. This implies that the usually predicted hydrological change in the future may be just the result of bias nonstationarity due to natural climate variability. In addition, the increase in the length of time window mitigates the impacts of bias nonstationarity on streamflow projections. Thus, it suggests using a long period (e.g. 50 years) for calibrating a bias correction method in hydrological climate change impact studies. The climatology of 20 or 30 years may not be enough for calibrating a bias correction method for impact studies.

In this study, the CRU dataset was used as the baseline to calculate the bias of GCM-simulated precipitation and temperature. The CRU data is gridded data obtained by interpolating gauged precipitation and temperature to regular grids. Due to limited gauges used for interpolation, the CRU dataset may also be biased. In order to investigate the uncertainty related to the baseline dataset, another dataset produced by the University of Delaware (version 4.01, available at www.esrl.noaa.gov/psd/data/gridded/data.UDel_AirT_Precip.html), with a spatial resolution of 0.5° from 1901 to 2000, the same as the CRU dataset, is also used to estimate the difference in bias between each historical period and the baseline period. Figure 10 shows the difference in bias for mean annual precipitation and temperature obtained from CRU and the University of Delaware (USA) dataset for 10-, 20-, 33- and 50-year windows. Similarly, the difference in bias varies with time for all time windows. Meanwhile, the range of difference in bias decreases with the extension of the time window. Generally, the range of difference in bias is similar for two datasets, even though the magnitude is not exactly the same.

The role of bias nonstationarity on future climate changes may rely on greenhouse gas emission scenarios, as different scenarios predict different future climate change signal. All the above results are based on the scenario of RCP4.5. In order to investigate the uncertainty related to the greenhouse gas emission scenario, the climate change signal of mean annual precipitation and temperature predicted by RCP4.5 and RCP8.5 are compared (Figure 11). The results showed that climate change signals predicted by RCP4.5 and RCP8.5 are comparable for precipitation, while the former is much smaller than the latter for temperature in terms of both median value and uncertainty related to GCMs for all time windows. Thus, the use of different greenhouse gas emission scenarios would not change the conclusion that the role of precipitation bias nonstationarity is important in future precipitation change. However, the impacts of temperature bias nonstationarity may depend on emission scenarios. For a high emission scenario (e.g. RCP8.5), the bias nonstationarity will become even less important, while for a low emission scenario (e.g. RCP4.5), it may become important for near future periods.

For the historical period, the difference in bias between two periods is caused by natural climate variability, because the anthropogenic climate change is pre-removed. Thus, to investigate the role of bias nonstationarity in future climate change, the range of bias is defined as ‘noise’, i.e. the difference between maximum and minimum values of biases across multiple periods with outliers deleted. An alternative approach using the standard deviation of bias over multiple periods can also be used to estimate natural climate variability. Table 2 presents a comparison of the median of range of bias and standard deviation of bias over historical periods. The median value is obtained from all 17 GCMs. The results show that the range of bias is much larger than the standard deviation of bias for both precipitation and temperature, especially over the short time window. When using the standard deviation of bias instead of the range of bias, the lower values of ‘noise’ can lead to larger SNR, indicating the impacts of bias nonstationarity would be less important. However, the impacts of precipitation bias are still great for most GCMs in the near and mid future, as indicated by the SNR being smaller than 1.

This study only investigated the bias nonstationarity of precipitation and temperature for the historical period. Thus, the bias nonstationarity only resulted from natural climate variability. However, for a future period, bias nonstationarity resulted from both natural climate variability and climate model sensitivity (Chen et al. 2015; Velázquez et al. 2015; Hui et al. 2019). The climate model sensitivity is more significant than natural climate variability, especially for precipitation bias nonstationarity (Hui et al. 2019). Climate models were proposed using different structures and parameterization schemes to represent the atmospheric physics, due to the limitation of inadequate knowledge of climate system, such as cloud and convective precipitation mechanisms, surface albedo feedback and land–atmosphere interactions. Thus different models simulate somewhat different responses to the same external forcing. This may lead to different biases between different climate models and observations. As the attribution of both natural climate variability and climate model sensitivity, the performance of a bias correction method may be more deteriorated in the future climate change than historical period. This will transfer to hydrological simulations, implying that hydrological simulations forced by bias-corrected precipitation and temperature may perform even worse when taking into account the impacts of climate model sensitivity. The impacts of climate model sensitivity on bias nonstationarity can be investigated in future studies.

Due to the limited length of observed data, 100 years data were divided into 10, five, three and two non-overlapping periods respectively corresponding to four different window lengths (10, 20, 33 and 50 years). For each time window, the biases over these finite periods were used to represent the range of biases (RB) caused by natural climate variability. However, small sample size, e.g. two non-overlapping periods for the 50-year window, may lead to large uncertainty of the estimated RB. To investigate the impacts of sample size on the estimation of bias variability, the use of multiple member GCMs may be a solution (Chen & Brissette 2019). In addition, this study only used the Two-Parameter Monthly Water Balance Model to simulate monthly runoff. The hydrological model uncertainty was not considered. All these can be avenues for future studies.

This work was partially supported by the National Natural Science Foundation of China (Grant No. 51779176 and 51539009) and the National Key R&D Program of China (2019YFC0408903). The authors would like to acknowledge the contribution of the World Climate Research Program Working Group on Coupled Modelling, which is responsible for CMIP. We wish to thank the climate modeling groups (listed in Table 1 of this paper) for producing and making available their model output.

Data cannot be made publicly available; readers should contact the corresponding author for details.