Investigating statistical bias correction with temporal subsample of the upper Ping River Basin, Thailand

This study aims to investigate different statistical bias correction techniques to improve the output of a regional climate model (RCM) of daily rainfall for the upper Ping River Basin in Northern Thailand. Three subsamples are used for each bias correction method, which are (1) using full calibrated 30year-period data, (2) seasonal subsampling, and (3) monthly subsampling. The bias correction techniques are classified into three groups, which are (1) distribution-derived transformation, (2) parametric transformation, and (3) nonparametric transformation. Eleven bias correction techniques with three different subsamples are used to derive transfer function parameters to adjust model bias error. Generally, appropriate bias correction methods with optimal subsampling are locally dependent and need to be defined specifically for a study area. The study results show that monthly subsampling would be well established by capturing the monthly mean variation after correcting the model’s daily rainfall. The results also give the best-fitted parameter set of the different subsamples. However, applying the full calibrated data and the seasonal subsamples cannot substantially improve internal variability. Thus, the effect of internal climate variability of the study region is greater than the choice of bias correction methods. Of the bias correction approaches, nonparametric transformation performed best in correcting daily rainfall bias error in this study area as evaluated by statistics and frequency distributions. Therefore, using a combination of methods between the nonparametric transformation and monthly subsampling offered the best accuracy and robustness. However, the nonparametric transformation was quite sensitive to the calibration time period.


GRAPHICAL ABSTRACT INTRODUCTION
One of the main issues in studying the impacts of climate change at a local scale is obtaining reliable output data from a climate model. Outputs from global climate models (GCMs) and/or regional climate models (RCMs) have limited capacity to capture catchment-scale climate variations (Teutschbein & Seibert ; Fang et al. ). RCMs are preferable to GCMs for providing more reliable results, specifically at a watershed scale. Maraun et al. () pointed out that the use of RCMs to downscale GCM outputs can add value to research because RCMs better resolve small-scale variability and regional processes over the intended local scale. However, RCM outputs still contain systematically inherited random error from GCMs (Maraun ) and/or RCMs. For example, some physical small-scale processes such as cloud formation, developing convective rainfall, and/or orographic rainfall are still not fully understood (Boé et  They found that gamma-based quantiles offer the best accuracy and robustness compared with linear (scaling correcting factor) and nonlinear (power transformation factor) bias correction techniques, while empirical QM could yield highly accurate results but was very sensitive to the calibration time period. Similar to a study by Fang et al.
(), they performed a comparative study of bias correction methods in downscaling RCM daily meteorological outputs to station scale over an arid area in China. However, this study showed that power transformation and QM performed equally well in correcting the precipitation frequency indices. Sharma () applied three different bias correction methods to correct monthly precipitation data retrieved from selected GCMs. Scaling bias correction and empirical-gamma and gamma-gamma transformation were tested and applied to the daily rainfall over nine years (1991)(1992)(1993)(1994)(1995)(1996)(1997)(1998)(1999)  Bias correction methods estimate parameters statistically from the observation data, the so-called transfer function or correcting factor, to adjust and correct the climate model bias error. However, there are no common bias correction techniques that are suitable for application to an RCM and to a specific area. Different regional climate characteristics and different models represent different bias errors, as evidenced by the many studies mentioned previously. A specific bias correction method might be appropriate for a certain region and also for a specific climate (GCM and/or RCM) model output (Chen et al. ).
As a result, it is of significance to find a proper bias correc- Moreover, suitable bias correction methods with an optimal subsample are locally dependent, and need to be defined for the specific catchment scale of the study area. In this study, basin average of daily rainfall product from RCM is tested instead of testing for all grid cells. After known appropriate bias correction methods with optimal subsamples, all the grid cell data can be manipulated to correct bias correction for other studies in detail and for a specific purpose.
Daily rainfall data are difficult to downscale due to their characteristics of non-normal distribution and discontinuity in both space and time (Cannon ). Eventually, it is sensitive in producing watershed runoff, and hence assessing climate change impacts on watershed hydrology (Chen et al. ). The present study will therefore address only the correction of daily rainfall bias error of RCM outputs. coverage. The study area is the drainage area of a major multipurpose reservoir, Bhumibol Dam, which serves to generate hydropower, supply irrigation water, and prevent floods in the downstream area, and is used for navigation, environmental conservation, and recreation purposes. More importantly, the dam also provides water to the industrial zone in the central region and domestic water resources to the capital city, Bangkok. The reservoir has a live storage capacity of 9.7 billion cubic meters (bm 3 ) and a total capacity of approximately 13.5 bm 3 . The average annual inflow is approximately 6.6 bm 3 , and the hydroelectric generation Historically, the spatial distribution of precipitation in the basin exhibits a pattern of higher rainfall intensity in mountainous areas than in the low-land area. Based on a completely observed data set during 1971-2000 (48 rain gauges), the average annual rainfall throughout the upper Ping Basin is approximately 1,106 mm. The average monthly rainfall variation is presented in Table 1. The table shows that 80% of rainfall occurs during the rainy season (mid-May to mid-Oct), whereas 20% of rainfall occurs during the dry season. This shows two distinct rainfall amounts during the wet and dry periods over a year.

Climate models and ground observation data
The climate model products were originally supplied by two selected GCMs, which were ECHAM5 and CCSM3. These  (1) Extract the RCM daily rainfall output from each climate model.
(2) Correct ground observation daily rainfall data and fill the gap as described above.
(3) Calculate basin average rainfall using the Thiessen polygon method for climate model outputs (P m ) and observation data (P o ).
(4) Prepare different time windows for both climate data and observation data as (5) Estimate the adjustment factor or transfer function parameters of each statistical bias correction method using each subsample to correct the climate model data.
(6) Apply the transfer function (TF) to the climate model data (P m ), and then obtain the corrected model bias error (P 0 m ) (Equation (2)). (7) Evaluate each bias correction approach with different subsamples.  The general expression when calculating bias is where P 0 m is the corrected model bias error. The details of each bias correction method for each group are summarized here.

Distribution-derived transformations
Distribution-derived transformations are a theoretical mixed distribution for estimating and assuming the observed and model data distribution parameters, conventionally known as QM. Bernoulli function with proper distribution was used to estimate the probability of rainfall events between dry days (no rainfall, or less than the prescriptive threshold value) and wet days. Common distributions for rainfall occurrence (wet days) are gamma, Weibull, log-normal, and exponential distributions. Therefore, four different mixed distributions, which were (1) Bernoulli-gamma (BG), (2) Bernoulli-Weibull (BW), (3) Bernoulli-log-normal (BL), and (4) Bernoulli-exponential (BE) were used to estimate the cumulative distribution function (CDF) parameters for both the observed CDF (P o ) and the model CDF (P m ) data.
The probability distribution function (pdf) of the Bernoulli function can be defined as the probability of a dry day (π) when rainfall data are less than 0.1 mm (threshold value) and the probability of a wet day (1 À π) when rainfall data are greater than 0.1 mm. The general formula is For rainfall occurrence (wet day), the Gamma, Weibull, log-normal, and exponential probability distributions (pdf written as f(P)) are commonly applied. Therefore, the mixed distribution between Bernoulli function with this four distributions CDFs, F(P), can be defined from where f (P) is a common rainfall occurrence probability distribution, as mentioned.
Theoretically, bias correction is intended to adjust the model data CDF (F m ) to have the same CDF as the observed data (F o ). This assumption can be stated as where α is the shape parameter and β the scale parameter which are derived from the observed data and models' data in order to find their fit with the theoretical distribution.
Then, the corrected model bias can be estimated as In this way, the inverse CDF of the observed data is applied as a transfer function to correct the model bias data such that the adjusting factor is obtained from the relationship between the two datasets' distribution parameters (i.e., α, β).

Parametric transformation
The parametric transformation used in this study consists of four different methods, which are (1) linear transformation (abbreviated as PTL) (2) power transformation (abbreviation as PTP) (3) scale transformation (abbreviated as PTS) (4) exponential tendency to an asymptote transformation (abbreviated as PTE) Equations (7) and (9) are simple scaling bias correction methods that are generally used in adjusting climate variables. Equations (8) and (10)  The general form of EQ is presented in Equation (2).   and CCSM3 (1970CCSM3 ( -1999, of the 30-year model hindcast run. Then, a 30-year calibrated period is used. Typically, the study area has three seasons, which are hot summer, rainy season, and winter season. However, from a hydrological point of view, the climate of the study area can be classified into a wet and dry season. This is because monthly rainfall variation is a clear-cut indicator for differentiating the trends of the wet and dry periods. Therefore, three different subsamples of data were chosen for testing in this study, as shown in Table 2. The total number of data (n) used to estimate the TF parameters are also displayed in Table 2.
The number of subsamples equally yields the number of TFs for each bias correction method and for each 30- year-period calibrated data set. The observation data were corrected from ground-based stations in the study area, as mentioned. A real averaging for the observed data and the raw outputs of RCMs using the Thiessen polygon approach are calculated in order to obtain the basin average daily rainfall. Three different subsamples were

PERFORMANCE EVALUATION
Performance metrics for assessing corrected bias daily rainfall products require some standard statistics such as mean, maximum, minimum, standard deviation, variance, root mean square error (RMSE), and mean absolute error (MAE). RMSE and MAE are used to measure the average magnitude of the corrected daily rainfall error as follows: where P i,obs is the observed daily rainfall, P i,corr the cor-

RESULTS
After taking the basin average daily rainfall for the two climate models (ECHAM5-MM5 and CCSM3-MM5) and the observation data, each climate product with the observed data is compared before bias correction in terms of statistical performance, monthly mean, and cumulative distribution curve to observe the differences between datasets. Next, the two models' daily rainfall after bias . Therefore, the two climate models' daily rainfall after bias corrections could also be compared to evaluate whether these models still yield a capability similar to that before the bias correction methods were applied.
Statistical performance before bias correction produced daily rainfall output with statistical values closer to the observed data than did the ECHAM5-MM5.

CDFs and monthly mean comparison before bias corrections
The CDFs of observed daily rainfall and the two climate models' raw datasets are plotted, as shown in Figure 3(a)).
The CDF magnitude from ECHAM5-MM5 and CCSM3-MM5 are much higher than the observed CDF, especially for the maximum cumulative distribution. This is because the two climate models produce higher maximum daily rainfall, as shown before. Thus, ECHAM5-MM5 gives a larger maximum CDF rainfall value (Max rain ¼ 243.81 mm) than CCSM3-MM5 (Max rain ¼ 94.46 mm) or the observed data (Max rain ¼ 54.90 mm).
The monthly mean variation of the three datasets is shown in Figure 3 It is important to point out that generally, the study area monthly average observed rainfall has two major changing points, in May and September.   Figure 5 shows that the distribution-derived transformation yielded a slightly better improvement in the CCSM3-MM5 CDFs. To conclude, the nonparametric transformation bias corrections give CDFs identical to the observed CDF for both climate models.
When considering the subsample aspect (by row), using the monthly subsample to adjust the two models' output yield almost matched the CDFs to the observations. However, the subsample performance was difficult to see and assess from the CDF plots.

Statistical performance after bias corrections
To confirm the bias correction performance with the opti-   transformation, Bernoulli-Gamma, was likely the best in the group for the corrected ECHAM5-MM5 model bias.

ECHAM5-MM5 performance after bias corrections
Like the bar plot in Figure 6, it can be clearly seen that the nonparametric transformation bias correction in group (3) (the right-most bar) shows the best statistics, average, maximum, standard deviation, and variance of all the methods.
Thus, for ECHAM5-MM5, the nonparametric bias correction group performs best in adjusting the model bias error, followed by the parametric bias correction method (middle column). From the same plot, it is obvious that the distribution-derived bias correction method yields inconsistency in bias adjustment, specifically for the BE bias correction, and gives the poorest performance. It may be that the model and observed data are not likely to fit with the BE distribution. For any derived distribution, the more unlikely the fit, the greater the residual error. After investigating the details, it is found that the BE and the BL were the two distributions most sensitive to fractional number. This is because

CCSM3-MM5 performance after bias corrections
The descriptive statistics for CCSM3-MM5's corrected daily rainfall are shown in Table 5. In the table, shaded values refer to the best statistical performance of each method in the group, and the underlined values with bold text are the closest statistical values to the observation data of all 11 bias correction adjustments.
The second and third columns in Table 5 present the raw output of CCSM3-MM5 daily rainfall and observed daily rainfall, respectively. It is clearly noted that CCSM3-MM5 produces higher daily rainfall than the observed data. However, after the bias corrections, the model results were almost equal to the observed data for all quantitative statistics. From the best bias-corrected results (underlined with shaded bold text) using a whole data subsample, the corrected model data of mean daily rainfall is 3.06, which is close to observed data value at approximately 2.98 mm; maximum daily rainfall is 54.9, which is equal to the observed data; the standard deviation yields 4.89, which is nearly the observation at 4.85 mm; the model variance is 23.9 mm, which is slightly higher than the observed data variance at 23.5 mm. Lastly, the model residual error after bias correction of RMSE improved from 8.2 to 5.61. The In each bias correction group in Table 5, the best practice for adjusting the CCSM3-MM5 model bias error (for all subsamples) using nonparametric transformation is EQL (maximum underlined with bold text values). For the parametric transformations, the exponential tendency to an asymptote transformation (PTE) performs the best of all methods in the group. Lastly, the distribution-derived transformation, BW, is likely the best in the group.
The bar chart in Figure 7 also displays all the descriptive statistics from Table 5 for CCSM3-MM5 bias correction.
One can see that all nonparametric transformation methods which are BG, BW, BL, and BE) produces quite inconsistent results between each method in the group. Specifically, BE tends to overestimate in correcting daily rainfall. As discussed before, the BE distribution is more sensible for adjusting model bias error.
To summarize, the statistical performance evaluations for the corrected ECHAM5-MM5 and CCSM3-MM5 daily rainfall show that nonparametric transformations are the best for improving the bias error of both climate models.
All three nonparametric transformations (EQ, REQ, and SS) showed equal performance. The second-best was the parametric transformation bias correction. The worst performance was by the distribution-derived transformation, as BE was the most sensitive to bias correction and finer subsamples (monthly), leading to some outliers in the maximum daily rainfall, as seen in the ECHAM5-MM5 and CCSM3-MM5 boxplots (last row of Figures 8 and 9). Thus, the use of BE methods is quite sensitive to correct climate model bias error.

Monthly mean comparison after bias corrections
To discover the monthly variation, daily rainfall data were aggregated by month and averaged out for a 30-year period, and then plotted as monthly rainfall. We attempted to find an optimum subsample and investigated the performance of different derivations of transformation function parameters. Thus, the monthly rainfall variation graphs disclosed better variation of both climate models after bias correction, as shown in Figures 10 and 11.   The CCSM3-MM5 monthly variation after bias correction is shown in Figure 11. Generally, bias-corrected Comparing between subsamples, this study showed that applying a monthly subsample yielded better performance for bias correction. Monthly subsamples are well established by capturing the internal monthly variation by giving the bestfitted parameter set. Therefore, this study could overcome the assumption that the reduced data when the sample size is split could reduce the robustness of bias correction, as claimed by Reiter et al. (). It may be that the sample size (n) when using a monthly subsample is still large enough to derive the bias function parameters. Moreover, the decreasing in the amount of data could be compensated for by using a monthly subsample that could remedy model bias error. However, applying the whole calibrated data set and seasonal subsamples yielded the least effectiveness, and also gave quite similar results for each climate model bias adjustment. This was because the study area had distinct wet/dry periods, which caused those two subsamples to be similar. The study area had two high rainfall intensities occurring in May and September. Thus, the effect of the region internal climate variability was also a major factor when implementing bias correction with an optimal subsample.
By comparing the ECMAM5-MM5 and the CCMS3-MM5 after bias correction, the study showed that the corrected ECHAM5-MM5 daily rainfall was greatly improved by more than CCSM3-MM5. The main reason was that the ECHAM5-MM5 produce better matching in monthly rainfall variation than CCSM3-MM5, although the CCSM3-MM5 raw output gave all descriptive statistical values closer to observation than ECHAM5-MM5. This evidence shows that the choice of RCMs/GCMs is a concern, because model uncertainties lead to a different result. Eventually, the output of RCMs/GCMs give different performance biascorrected results that affect the use of model data in studying the impact of climate change.
Thus, the effect of the study region internal climate variability is greater than the choice of bias correction method. The need for an optimal subsample and bias correction method is important when adjusting the model bias error. In this study, the monthly subsample is the most appropriate for the upper Ping River Basin. Of the bias correction approaches, the nonparametric transformation performs best in correcting daily rainfall bias error as evaluated by statistics and frequency distribution. In summary, using a combination of methods between the nonparametric transformation and monthly subsample offered the best accuracy and robustness for bias correction of the climate model daily rainfall in this study area. However, nonparametric transformation is quite sensitive to the calibration time period.