Abstract

This study aims to explore possible distributional changes in annual daily maximum rainfalls (ADMRs) over South Korea using a Bayesian multiple non-crossing quantile regression model. The distributional changes in the ADMRs are grouped into nine categories, focusing on changes in the location and scale parameters of the probability distribution. We identified seven categories for a distributional change in the selected stations. Most of the stations (28 of 50) are classified as Category III, which is characterized by an upward trend with an increase in variance in the distribution. Moreover, stations with a downward trend with a decrease in the variance pattern (Category VII) are mainly distributed on the southern Korean coast. On the other hand, Category I stations are mostly located in eastern Korea and primarily show a statistically significant upward trend with a decrease in variance. Moreover, this study explored changes in design rainfall estimates for different categories in terms of distributional changes. For Categories I, II, III, and VI, a noticeable increase in design rainfall was observed, while Categories IV, V, and VII showed no evidence of association with risk of increased extreme rainfall.

HIGHLIGHTS

  • A Bayesian multiple non-crossing quantile regression (BMN-QR) model is developed.

  • A BMN-QR model provides a reliable way of detecting distribution changes.

  • We identified seven categories for a distributional change in the selected stations.

  • Most of the stations have an upward trend with an increase in variance in the distribution.

  • A BMN-QR nonstationary frequency model can identify dynamic distributional changes.

Graphical Abstract

Graphical Abstract
Graphical Abstract

INTRODUCTION

It has been reported that precipitation changes associated with a climate change can differ spatially and temporally for a region due to the geographic distribution of rain, which is largely driven by climate and topography (Schmidli et al. 2006; Pall et al. 2007; Giorgi et al. 2011; Trenberth 2011; So et al. 2015; Donat et al. 2017). Especially, changes in extreme rainfall events such as floods and droughts have been the main causes of natural disasters in many parts of the world (Keellings & Hernández Ayala 2019; Kim et al. 2011; King et al. 2017; Roderick et al. 2019). It has been acknowledged that hydrologic circulation is likely to be accelerated with the increase of water vapor under climate change such that extreme rainfall events are more likely to occur more frequently (Trenberth et al. 2003; Trenberth 2011; Zhang et al. 2013; Fischer & Knutti 2014, 2016; O'Gorman 2015; Roderick et al. 2019). The Intergovernmental Panel of Climate Change (IPCC) reported that occurrences of frequent extreme rainfall events have been observed in many regions and have been recognized as an important global issue (Gray 2007). Overall, the societal infrastructure will become more vulnerable to extreme rainfall conditions in a changing climate (Kunkel et al. 1999; Easterling et al. 2000; Jahn 2015; Nissen & Ulbrich 2017).

Nonstationarity in precipitation extremes has been a major concern for decades in the currently adopted stationary modeling procedure. Changes in rainfall distribution affect current operations and future design criteria for water resource systems; thus, it might be necessary to explore rainfall trends. In this context, statistical modeling of extreme rainfall trends (or nonstationarity) has been extensively conducted. In particular, an increasing number of studies have been carried out to explore rainfall trends (Matti et al. 2009; Alexander et al. 2011; Clarke et al. 2011; Dravitzki & McGregor 2011; Jung et al. 2011; Río et al. 2011; Rana et al. 2012; Westra et al. 2013). To be more specific, various approaches such as the Mann–Kendall (MK) test (Mann 1945; Kendall 1975; Gilbert 1987), Kwiatkowski–Phillips–Schmidt–Shin (KPSS) test (Kwiatkowski et al. 1992), and Augmented Dickey–Fuller (ADF) test (Dickey & Fuller 1981) have been extensively used to confirm stationarity. However, these methods have certain limitations as they are both sensitive to outliers and serial correlations (i.e., persistence) within time series that can lead to inflated estimates of temporal trends. Moreover, a traditional approach for the trend analysis is not appropriate for extreme rainfall events (Muhlbauer et al. 2009; Timofeev & Sterin 2010) as noticeable distribution changes in the extremes and the hydro-meteorological variables are often undetectable in traditional approaches (Shiau & Huang 2015).

Many studies have applied trend analysis to hydrologic variables such as precipitation, temperature, and discharge in South Korea (Chang & Kwon 2007; Kwon et al. 2008b; Lee et al. 2012), showing that long-term trends in annual total rainfall are increasing in the northern part of South Korea and decreasing in the southwestern part of South Korea. Jung et al. (2011) employed the ordinary linear regression and the MK test to identify the linear trends in rainfall for South Korea and used a point rainfall frequency analysis based on the generalized extreme value (GEV) distribution to explore changes in the magnitude of 2-day maximum rainfall. They concluded that, although the existing studies revealed significant trends in rainfall over South Korea, several of the methods were insufficient to detect substantial distribution changes in the extreme since the existing approaches aim to identify possible monotonic trends in any given time series.

For these reasons, quantile regression was developed to better identify distinct changes in time-series data at any percentile values (Koenker & Bassett Jr 1978). A traditional ordinary regression method is designed to estimate the conditional expectation (mean) of the response variable given certain values of the independent variables. In contrast, a quantile regression method allows flexibility to estimate either the conditional median or any quantiles of the response variable (Sankarasubramanian & Lall 2003). Quantile regression has been applied to explore changes in annual precipitation over time in Zimbabwe (Chamaillé-Jammes et al. 2007; Mazvimavi 2010), in annual precipitation over the Midwestern United States (Villarini et al. 2011), and in monthly precipitation over the Southeastern United States (Wang et al. 2014). Moreover, quantile regression has been used to estimate predictive uncertainty in hydrological models (Weerts et al. 2011; Roscoe et al. 2012; Dogulu et al. 2015) and nonstationary quantiles in floods from a Montana basin (Sankarasubramanian & Lall 2003; Wasko & Sharma 2014; Fan & Chen 2016). However, most of the studies have been primarily based on specific distribution quantiles (e.g., 0.9 and 0.1). Such limited results in terms of precipitation changes do not reveal changes in the entire distribution (Friederichs & Hense 2007; Roth et al. 2015; Malik et al. 2016). More recently, a technique was applied to explore possible distributional changes of rainfall characteristics using quantile regression (Shiau & Huang 2015). Although quantile regression has been successfully applied in hydrological applications, the systematic estimation of parameter uncertainty has not been properly addressed. In this context, this study applies a Bayesian parameter estimation approach to a quantile regression model to better explore the distributional changes in extreme rainfall for given years. More specifically, Bayesian inference for the quantile regression is used to categorize the distributional changes in extreme rainfall into certain classes and then further utilized to provide a nonstationary estimation of the design rainfall corresponding to the classes for a better understanding of changes in design rainfall estimates in a changing climate. However, estimated conditional quantile functions often fail to be a monotonically increasing function of quantile, leading to quantile crossing problems. To alleviate crossing issues in quantile functions, this study develops a Bayesian multiple non-crossing quantile regression (BMN-QR) approach within a fully Bayesian framework, which is applied to multiple rainfall stations in South Korea.

The theoretical background for the BMN-QR and the rainfall data used in this study is discussed in the following section. The BMN-QR is then applied to the annual daily maximum rainfalls (ADMRs) over the period of 1973–2015, and a comprehensive discussion of the results is provided in the following section. Finally, a summary and conclusion from the study are offered.

DATA AND METHODOLOGY

Study area and data

South Korea tends to have a humid continental and subtropical climate. The country's terrain is mostly mountainous, with lowlands located in the western and southeastern regions. Due to the mountainous terrain, the country experiences uneven rainfall distribution. Annual precipitation ranges from 1,000 to 1,800 mm in the southern part of South Korea and from 1,100 to 1,400 mm in the central part. More than half of the annual precipitation happens during the monsoon season when a stationary front forms across the Korean Peninsula for about a month in summer (Kwon et al. 2008b). Winter precipitation is less than 10% of the total annual precipitation. The typhoon and the movement of air masses containing a relatively high concentration of water vapor from the Asian continent and the Pacific Ocean over the Korean Peninsula have a significant influence on extreme rainfall in the Korean Peninsula (Oh et al. 2011; Son et al. 2017; Lima et al. 2018). Out of approximately 30 typhoons generated annually in the northwest Pacific Ocean, two or three influence the Korean Peninsula from June to October (Kwon et al. 2008b). In general, there is a relatively large seasonal variation of water resources in South Korea.

This study obtained daily precipitation data of 93 stations, subjected to a certain level of quality assurance by the Korea Meteorological Administration (KMA). As seen in Figure 1, rainfall stations uniformly cover South Korea, including Jeju Island. We considered length and the proportion of missing values for quality assessment criteria. The periods of record vary by station, but for this study, we select stations with less than 5% missing value and with continued observed data of 43 years from 1973 to 2015.

Figure 1

Spatial distribution of rainfall observation stations. Blue circles indicate the stations selected for the study, and digital elevation map (DEM) is in yellow. Please refer to the online version of this paper to see this figure in colour: http://dx.doi.org/10.2166/nh.2020.003.

Figure 1

Spatial distribution of rainfall observation stations. Blue circles indicate the stations selected for the study, and digital elevation map (DEM) is in yellow. Please refer to the online version of this paper to see this figure in colour: http://dx.doi.org/10.2166/nh.2020.003.

Bayesian multiple non-crossing quantile regression

For the analysis of distributional changes in ADMRs, we propose the following quantile regression model. Let denote the ADMRs for year i and p denote the quantile level , then the generic structure of the model is given by: 
formula
(1)
where is a vector of regressors, is a vector of regression coefficients for the conditional quantile , is an unspecified quantile-specific error term where no distribution is specified other than the constraint and estimated by minimizing the asymmetrically weighted absolute deviations through linear programming. Here, the prime () denotes the transpose of the matrix x. As discussed in Koenker & Bassett (1978), no distribution is specified for . The only existing condition is , with the estimates obtained by minimizing the error according to: 
formula
(2)
where is the loss function given by the following equation (Friederichs & Hense 2007): 
formula
(3)
In general, the above estimates for the quantile regression model are solved separately through linear programming for each of the q quantiles, . In this setting, the estimated conditional quantile functions for a set of predictors often fail to be a monotonically increasing function of the quantile, which ultimately leads to a quantile crossing problem (Koenker & Regression 2005). In particular, the crossing intensifies in the upper quantile and in the presence of a larger number of predictors. In order to minimize the crossing problem, we develop a multiple non-crossing quantile regression model in which all quantile functions are simultaneously estimated under the non-crossing constraints within a fully Bayesian framework, which is subject to: 
formula
(4)
The use of Bayesian inference in generalized linear and additive models is standard due to its ability to estimate the posterior distribution of parameters even in more complex models, which allows us to explicitly account for parameter uncertainty in time-series analysis (Yu & Moyeed 2001). For these reasons, this study adopted a more flexible Bayesian modeling approach. A Bayesian inference commonly requires a likelihood; thus, a typical approach to Bayesian quantile regression is to assume an asymmetric Laplace distribution (ALD) for (Yu & Moyeed 2001), which produces minimization of the loss function (Equation (3)) equivalent to the maximization of a likelihood function of independently distributed ALD (Yue & Rue 2011). The ALD is characterized by a set of parameters such as location (), precision (), and skewness () by letting to ensure. The density is given by: 
formula
(5)
where denotes the indicator function.
According to features of ALD (Yue & Rue 2011), can be written as a linear combination of two random variables: 
formula
(6)
where W is an exponentially distributed random variable with the rate parameter and Z has a standard normal distribution N(0,1). Posterior estimates can be subsequently obtained using Bayesian updates conditioned on the exponential random variable W. Using the above stochastic representation with Equation (1), the proposed quantile regression model can be rewritten as: 
formula
(7)
where and are subject-specific values of W and Z, respectively, and is defined in Equation (1).

Bayesian inference for parameter estimation

To enable the inclusion and improved consideration of the uncertainties in the proposed quantile regression model parameters, we adopted a Bayesian inference framework. Following Bayes' rule, the joint posterior for the set of parameters given the observed data and the prior distribution of the parameter can be written as: 
formula
(8)
where is the likelihood function.

The prior distribution of the parameter should reflect our knowledge of the data before observing it. For the parameters of the proposed model, we do not have any prior knowledge of their values other than that regression parameters can be any real number, and the ALD precision parameter should be greater than zero. Hence, we assume a conjugate normal prior to regression coefficients and a gamma prior to ALD precision parameter.

The posterior inference for the desired quantiles, simultaneously proceeds via data augmentation by introducing observation specific latent weights,, as specified in Equation (7). Furthermore, let, where for brevity and include the regression coefficients and regression functions, as seen in Equation (1). Therefore, given and, are independent and identically distributed random variables with normal distribution: 
formula
(9)
Therefore, the joint posterior for all model parameters is given by: 
formula
(10)

Analytical integration of the joint posterior distribution, Equation (10), for all model parameters is not tractable, and this study adopts the widely used Markov Chain Monte Carlo (MCMC) method approach, specifically a Gibbs sampler. The Gibbs sampler investigates the parameter space by fixing all the parameters except one of the parameters. In other words, the conditional distribution of the parameter relies on the fixed values of the remaining parameters (Kwon et al. 2008a). After setting initial values to the model parameters, the algorithm iterates the following steps. (1) Sample from its full conditional normal distribution. (2) Sample from its full conditional gamma distribution. (3) For all, update the inverse latent weights for its exponential distribution using previously sampled. (4) Steps 1–3 are repeated until we obtain a sufficient number of samples for the proposed model in year. Here, we performed 30,000 simulations with 20,000 burnouts for three chains. The convergence of the parameter was evaluated by the Gelman–Rubin statistic (Gelman et al. 2014) with all values less than 1.1, suggesting convergence after 30,000 iterations with a 1,000 cycle burn-in.

Categorical approach to distributional change detection based on BMN-QR

Many studies have demonstrated strong evidence of precipitation changes under climate change. Such changes in precipitation patterns can lead to misleading effects on water resource management and can alter future design criteria of water resource systems.

Rainfall has very skewed probability distributions, with the distribution becoming more skewed to the right for enhanced extreme rainfall events such that the use of existing approaches based on the normality assumption might not be valid. At the beginning of the 1990s, following the increasing interest in extreme rainfall analysis, the use of quantile regression has become increasingly popular to better characterize the dynamics behind extreme hydrological events (Chamaillé-Jammes et al. 2007; Mazvimavi 2010; Wang et al. 2014) and has demonstrated a significant improvement over existing approaches. Based on the previous study (Shiau & Huang 2015), we propose a Bayesian quantile regression-based distributional change detection approach. More specifically, this study uses the entire marginal posterior distribution of the trend for a particular quantile to examine changes in rainfall characteristics and categorize them into classes. Figure 2 shows an example of the basic concept, which exemplifies the detection of distribution change using the Bayesian quantile regression.

Figure 2

An illustration of the construction of an empirical PDF of the ADMRs through quantile regression for a wide range of quantiles.

Figure 2

An illustration of the construction of an empirical PDF of the ADMRs through quantile regression for a wide range of quantiles.

Since the sign of the slope parameter on extreme rainfall can vary with quantiles, different behaviors in probabilistic distributions can be produced from year to year. Detecting the distribution changes in rainfall over time is thus made possible by comparing the shapes of derived empirical probability density functions (PDFs), as illustrated in Figure 2. To analyze the shapes of the empirical PDFs, a simplified comparison is performed to demonstrate changes in the location (central) and scale (dispersion) of the PDFs. Changes in the PDF scale focus on whether the shape of the PDF sharpens (i.e., less dispersion or extent of peakedness) or not, while location changes emphasize shifts of the PDF. Thus, there are nine possible categories in terms of distributional changes depending on both changes of scale and the location of PDFs, which are summarized in Table 1 and illustrated in Figure 3. We implemented the BMN-QR over a wide range of quantiles in order to drive the empirical PDFs, which are further utilized to categorize the distribution changes, as shown in Figure 3.

Table 1

Summary of the nine categories associated with distributional changes in the ADMRs

CategoryTrend typeChange in PDFs
ScaleLocation
Upward-convergent lines Yes (peakness) Yes (rightward) 
II Upward parallel lines No Yes (rightward) 
III Upward-divergent lines Yes (flatness) Yes (rightward) 
IV Horizontally convergent lines Yes (peakness) No 
Horizontally parallel lines No No 
VI Horizontally divergent lines Yes (flatness) No 
VII Downward-convergent lines Yes (peakness) Yes (leftward) 
VIII Downward parallel lines No Yes (leftward) 
IX Downward divergent lines Yes (flatness) Yes (leftward) 
CategoryTrend typeChange in PDFs
ScaleLocation
Upward-convergent lines Yes (peakness) Yes (rightward) 
II Upward parallel lines No Yes (rightward) 
III Upward-divergent lines Yes (flatness) Yes (rightward) 
IV Horizontally convergent lines Yes (peakness) No 
Horizontally parallel lines No No 
VI Horizontally divergent lines Yes (flatness) No 
VII Downward-convergent lines Yes (peakness) Yes (leftward) 
VIII Downward parallel lines No Yes (leftward) 
IX Downward divergent lines Yes (flatness) Yes (leftward) 
Figure 3

A graphical illustration of the nine categories associated with changes of PDFs for the ADMRs. The blue line represents the PDF in the first year, while the red line indicates the PDF in the last year. Please refer to the online version of this paper to see this figure in colour: http://dx.doi.org/10.2166/nh.2020.003.

Figure 3

A graphical illustration of the nine categories associated with changes of PDFs for the ADMRs. The blue line represents the PDF in the first year, while the red line indicates the PDF in the last year. Please refer to the online version of this paper to see this figure in colour: http://dx.doi.org/10.2166/nh.2020.003.

As seen in Figure 3, a graphical representation of the distributional changes in the variables of interest provides qualitative information on how the changes in the scale and location of the PDFs differ; thus, this study considers goodness-of-fit tests on the estimated empirical PDFs. More specifically, we propose a two-tier diagnostic approach to determine the distribution changes informed by the BMN-QR model. First, the BMN-QR is applied to the ADMRs over a wide range of quantiles. Second, the empirical PDFs for the first year (1973) and final year (2015) are derived from the entire marginal posterior distribution of the trend for a set of quantiles. Third, two-sample t-test and two-sample F-test are used to determine whether the means and variances of the two populations differ, respectively, in order to assign a single-specified category out of nine categories. Finally, two-sample KS and two-sample AD tests are further utilized to test the hypothesis that two samples could have come from the same underlying distribution. While the KS test has more power against deviations in the middle of the distribution, the AD test has better power against heavier tails (Razali & Wah 2011). A flowchart of the two-tier approach to determine the distribution change in the ADMRs is illustrated in Figure 4.

Figure 4

A flowchart of the two-tier diagnostic approach to determine the distribution changes informed by the BMN-QR model. In tier 1, two-sample t-test and two-sample F-test are used to determine whether the means and variances of the two populations differ. In tier 2, two-sample Kolmogorov–Smirnov (KS) and two-sample Anderson-Darling (AD) tests are used to test the significance of the difference in the distribution.

Figure 4

A flowchart of the two-tier diagnostic approach to determine the distribution changes informed by the BMN-QR model. In tier 1, two-sample t-test and two-sample F-test are used to determine whether the means and variances of the two populations differ. In tier 2, two-sample Kolmogorov–Smirnov (KS) and two-sample Anderson-Darling (AD) tests are used to test the significance of the difference in the distribution.

RESULTS AND DISCUSSION

Stationarity test

For comparison purposes, we adopted well-known and widely used statistical methods such as the MK test (Mann 1945; Kendall 1975; Gilbert 1987) for determining the existence of linear trends in the rainfall extremes and the KPSS test (Kwiatkowski et al. 1992) and ADF test (Dickey & Fuller 1981) for the existence of stationarity in the rainfall extremes. The results of the stationarity tests are presented in Table 2. Based on the MK test, three of 50 stations, such as Ulleungdo, Uljin, and Jecheon, showed a statistically significant linear trend during 1973–2015 at a 5% significance level. Six stations similarly showed statistically significant nonstationarity at a 5% significance level based on the KPSS test. On the other hand, 16 stations were found to be nonstationary by the ADF test. Using these tests, the presence of nonstationarity over stations might not be recognized or consistent, as it is affected by the stationarity test and outliers. Therefore, this study further explored the evidence of nonstationarity in the ADMRs by the proposed Bayesian quantile regression in the following section.

Table 2

Summary of trend and stationarity test results at the 5% significance level

St. No.Trend test
Stationary test
St. No.Trend test
Stationary test
MK
KPSS
ADF
MK
KPSS
ADF
Trendp-valueNonstationaryp-valueNonstationaryp-valueTrendp-valueNonstationaryp-valueNonstationaryp-value
No 0.363 No 0.100 No 0.058 26 No 0.090 No 0.100 No 0.061 
No 0.221 No 0.100 No 0.074 27 No 0.198 No 0.085 No 0.076 
No 0.420 No 0.100 Yes 0.003 28 Yes 0.013 Yes 0.034 No 0.150 
No 0.082 No 0.100 No 0.069 29 No 0.875 No 0.100 Yes 0.037 
Yes 0.003 No 0.100 No 0.131 30 No 0.516 Yes 0.048 No 0.113 
No 0.325 No 0.061 No 0.052 31 No 0.778 No 0.100 Yes 0.028 
No 0.464 Yes 0.027 No 0.107 32 No 0.834 No 0.100 Yes 0.021 
No 0.127 Yes 0.022 No 0.137 33 No 0.630 No 0.100 No 0.182 
Yes 0.016 No 0.100 No 0.077 34 No 0.842 No 0.100 No 0.130 
10 No 0.054 No 0.100 No 0.059 35 No 0.158 No 0.100 No 0.093 
11 No 0.859 Yes 0.028 No 0.079 36 No 0.062 No 0.100 Yes 0.045 
12 No 0.660 No 0.062 No 0.075 37 No 0.209 No 0.100 Yes 0.027 
13 No 0.615 No 0.100 Yes 0.005 38 No 0.402 No 0.100 Yes 0.032 
14 No 0.060 No 0.100 No 0.091 39 No 0.683 No 0.100 Yes 0.024 
15 No 0.834 No 0.100 Yes 0.012 40 No 0.158 Yes 0.013 No 0.155 
16 No 0.762 No 0.100 Yes 0.040 41 No 0.079 No 0.100 No 0.250 
17 No 0.933 No 0.100 Yes 0.023 42 No 0.983 No 0.100 Yes 0.027 
18 No 0.699 No 0.100 No 0.058 43 No 0.867 No 0.100 No 0.142 
19 No 0.645 No 0.100 Yes 0.043 44 No 0.124 No 0.094 No 0.133 
20 No 0.439 No 0.088 No 0.086 45 No 0.586 No 0.100 No 0.244 
21 No 0.451 No 0.100 No 0.061 46 No 0.402 No 0.100 No 0.058 
22 No 0.738 No 0.100 Yes 0.037 47 No 0.408 No 0.100 No 0.062 
23 No 0.834 No 0.100 No 0.120 48 No 0.088 No 0.100 No 0.152 
24 No 0.079 No 0.100 No 0.051 49 No 0.992 No 0.100 Yes 0.046 
25 No 0.103 No 0.100 No 0.050 50 No 0.445 No 0.100 No 0.086 
St. No.Trend test
Stationary test
St. No.Trend test
Stationary test
MK
KPSS
ADF
MK
KPSS
ADF
Trendp-valueNonstationaryp-valueNonstationaryp-valueTrendp-valueNonstationaryp-valueNonstationaryp-value
No 0.363 No 0.100 No 0.058 26 No 0.090 No 0.100 No 0.061 
No 0.221 No 0.100 No 0.074 27 No 0.198 No 0.085 No 0.076 
No 0.420 No 0.100 Yes 0.003 28 Yes 0.013 Yes 0.034 No 0.150 
No 0.082 No 0.100 No 0.069 29 No 0.875 No 0.100 Yes 0.037 
Yes 0.003 No 0.100 No 0.131 30 No 0.516 Yes 0.048 No 0.113 
No 0.325 No 0.061 No 0.052 31 No 0.778 No 0.100 Yes 0.028 
No 0.464 Yes 0.027 No 0.107 32 No 0.834 No 0.100 Yes 0.021 
No 0.127 Yes 0.022 No 0.137 33 No 0.630 No 0.100 No 0.182 
Yes 0.016 No 0.100 No 0.077 34 No 0.842 No 0.100 No 0.130 
10 No 0.054 No 0.100 No 0.059 35 No 0.158 No 0.100 No 0.093 
11 No 0.859 Yes 0.028 No 0.079 36 No 0.062 No 0.100 Yes 0.045 
12 No 0.660 No 0.062 No 0.075 37 No 0.209 No 0.100 Yes 0.027 
13 No 0.615 No 0.100 Yes 0.005 38 No 0.402 No 0.100 Yes 0.032 
14 No 0.060 No 0.100 No 0.091 39 No 0.683 No 0.100 Yes 0.024 
15 No 0.834 No 0.100 Yes 0.012 40 No 0.158 Yes 0.013 No 0.155 
16 No 0.762 No 0.100 Yes 0.040 41 No 0.079 No 0.100 No 0.250 
17 No 0.933 No 0.100 Yes 0.023 42 No 0.983 No 0.100 Yes 0.027 
18 No 0.699 No 0.100 No 0.058 43 No 0.867 No 0.100 No 0.142 
19 No 0.645 No 0.100 Yes 0.043 44 No 0.124 No 0.094 No 0.133 
20 No 0.439 No 0.088 No 0.086 45 No 0.586 No 0.100 No 0.244 
21 No 0.451 No 0.100 No 0.061 46 No 0.402 No 0.100 No 0.058 
22 No 0.738 No 0.100 Yes 0.037 47 No 0.408 No 0.100 No 0.062 
23 No 0.834 No 0.100 No 0.120 48 No 0.088 No 0.100 No 0.152 
24 No 0.079 No 0.100 No 0.051 49 No 0.992 No 0.100 Yes 0.046 
25 No 0.103 No 0.100 No 0.050 50 No 0.445 No 0.100 No 0.086 

Distributional changes in the ADMRs and their relations to design rainfall

Here, we aimed to identify the distributional changes in the desired quantiles using the entire posterior distributions inferred from the time-varying mean and variance. First, a comparison between the BMN-QR and Bayesian Quantile Regression (BQR) is provided to demonstrate that the BMN-QR is able to avoid the quantile crossing issue with the Chupungyeong and Geoje stations. As described in Figure 5(a) and 5(b), quantile crossing is seen for the Geoje station in the BQR for the upper quantile, while no quantile crossing issues are seen with the proposed BMN-QR approach. Similarly, as shown in Figure 5(c) and 5(d), for the Chupungyeong station, the quantile crossing issue is problematic in the BQR for the lower quantile, which is also solved by the BMN-QR approach. Therefore, we applied the BMN-QR to the ADMRs for a wide range of quantiles (e.g., 0.05, 0.15, 0.25, 0.35, 0.45, 0.55, 0.65, 0.75, 0.85, and 0.95) for all stations. The associated trends, along with Bayesian credible intervals, are illustrated in Figure 6.

Figure 5

A comparison between the BQR (a, c) and BMN-QR (b, d) on crossing quantiles at Geoje (a, b) and Chupungyeong (c, d) stations. The dotted lines indicate trends at different quantiles, while shaded areas represent the associated Bayesian credible intervals.

Figure 5

A comparison between the BQR (a, c) and BMN-QR (b, d) on crossing quantiles at Geoje (a, b) and Chupungyeong (c, d) stations. The dotted lines indicate trends at different quantiles, while shaded areas represent the associated Bayesian credible intervals.

Figure 6

The observed ADMRs and their trends along with the credible interval over different quantiles. The x-axis indicates the time in years, and the y-axis represents the ADMRs in millimeters. The coloured dash line presents quantiles from 0.5 to 0.95 with 0.1 increments. Please refer to the online version of this paper to see this figure in colour: http://dx.doi.org/10.2166/nh.2020.003.

Figure 6

The observed ADMRs and their trends along with the credible interval over different quantiles. The x-axis indicates the time in years, and the y-axis represents the ADMRs in millimeters. The coloured dash line presents quantiles from 0.5 to 0.95 with 0.1 increments. Please refer to the online version of this paper to see this figure in colour: http://dx.doi.org/10.2166/nh.2020.003.

The quantile crossing problem violates the basic principle of probability distributions. Moreover, the presence of crossing can lead to bias in time-varying quantile distributions. The quantile crossing issue is mainly due to the insufficient sample size, so that the crossing problems can be eventually alleviated as the sample size increases. In this study, the quantile crossing problem occurs at 40 stations out of 50 stations. To demonstrate quantile crossing at different regions, we selected two stations, Chupungyeong and Geoje, with the quantile crossing at the beginning and end of the period, respectively. Our results indicated that the overall shape of the distribution with quantile crossing was not significantly different from that of the non-crossing. However, there was a difference in the tail of the distribution, which in turn may lead to a large bias in the design rainfall. For instance, design rainfalls estimated from samples with the quantile crossing at the beginning of the period (i.e., Chupungyeong) are consistently smaller (or larger) than those obtained from the non-crossing at the first year (or last year), as summarized in Supplementary material, Table S1. On the other hand, the design rainfall quantities taken from samples with the quantile crossing at the end of the period (i.e., Geoje) are subsequently larger (or smaller) than those obtained from the non-crossing at the first year (or last year).

The two-tier diagnostic approach is introduced to detect the changes in ADMRs over South Korea, as illustrated in Figure 4. From the tier-1 criteria using the two-sample t-test and the F-test, we identified seven categories (i.e., Categories I, II, III, IV, V, VI, and VII) of a distributional change in the selected stations, as summarized in Table 3. Categories I and III showed an increasing trend in the mean value of AMDRs. In contrast, Categories VI and VII showed no change or a decreasing pattern in the mean value of AMDRs with associated variance. To ensure the identified category associated with distribution changes, the tier-2 criteria, two-sample KS, and AD tests were applied to assess the significance of the difference in the estimated distributions at the first and last years of the records for all stations (Table 4). The two-sample KS test showed that 25 out of 50 stations have significantly different distributions at a 5% significance level, while the two-sample AD test detected 32 stations with AD statistics over the critical value of 2.492, indicating that they have significantly different distributions.

Table 3

Categorization of distributional changes in the ADMRs using ‘Tier-1’ approach and summary of the significance tests for the parameters of the distributions at the 5% significance level

St. No.Station nameCategoryTwo-sample t-test
Two-sample F-test
ST. No.Station nameCategoryTwo-sample t-test
Two-sample F-test
Hp-valueHp-valueHp-valueHp-value
Sokcho 1 0.000 1 0.000 26 Inje 1 0.000 1 0.000 
Chuncheon III 1 0.000 1 0.048 27 Hongcheon III 1 0.000 1 0.000 
Gangneung III 1 0.000 1 0.000 28 Jecheon II 1 0.000 0.154 
Wonju III 1 0.000 1 0.000 29 Boeun VI 0.627 1 0.000 
Ulleungdo 1 0.000 1 0.000 30 Cheonan III 1 0.000 1 0.000 
Suwon III 1 0.000 1 0.000 31 Boryeong III 1 0.000 1 0.001 
Chungju III 1 0.000 1 0.000 32 Buyeo VII 1 0.000 1 0.000 
Seosan III 1 0.000 1 0.000 33 Geumsan III 1 0.000 1 0.000 
Uljin 1 0.000 1 0.000 34 Buan III 1 0.007 1 0.000 
10 Cheongju 1 0.000 1 0.022 35 Imsil III 1 0.000 1 0.000 
11 Daejeon III 1 0.000 1 0.000 36 Jeongeup III 1 0.000 1 0.000 
12 Chupungyeong III 1 0.000 1 0.000 37 Jangheung VII 1 0.003 1 0.000 
13 Pohang III 1 0.000 1 0.000 38 Haenam VII 1 0.069 1 0.000 
14 Daegu 1 0.000 1 0.000 39 Goheung 1 0.000 1 0.049 
15 Ulsan IV 0.339 1 0.000 40 Yeongju III 1 0.000 1 0.000 
16 Gwangju 0.140 0.152 41 Mungyeong III 1 0.000 1 0.000 
17 Busan IV 0.138 1 0.000 42 Yeongdeok VI 1 0.012 1 0.000 
18 Tongyeong VII 1 0.000 1 0.000 43 Uiseong III 1 0.000 1 0.000 
19 Mokpo VII 1 0.000 1 0.000 44 Gumi III 1 0.000 1 0.000 
20 Yeosu III 1 0.000 1 0.000 45 Yeongcheon 1 0.000 1 0.000 
21 Wando VII 1 0.000 1 0.000 46 Hapcheon III 1 0.000 1 0.000 
22 Jeju III 1 0.000 1 0.000 47 Miryang VII 1 0.000 1 0.000 
23 Ganghwa III 1 0.000 1 0.000 48 Sancheong III 1 0.000 1 0.000 
24 Yangpyeong III 1 0.000 1 0.000 49 Geoje VII 1 0.000 1 0.000 
25 Incheon III 1 0.000 1 0.000 50 Namhae III 1 0.000 1 0.000 
St. No.Station nameCategoryTwo-sample t-test
Two-sample F-test
ST. No.Station nameCategoryTwo-sample t-test
Two-sample F-test
Hp-valueHp-valueHp-valueHp-value
Sokcho 1 0.000 1 0.000 26 Inje 1 0.000 1 0.000 
Chuncheon III 1 0.000 1 0.048 27 Hongcheon III 1 0.000 1 0.000 
Gangneung III 1 0.000 1 0.000 28 Jecheon II 1 0.000 0.154 
Wonju III 1 0.000 1 0.000 29 Boeun VI 0.627 1 0.000 
Ulleungdo 1 0.000 1 0.000 30 Cheonan III 1 0.000 1 0.000 
Suwon III 1 0.000 1 0.000 31 Boryeong III 1 0.000 1 0.001 
Chungju III 1 0.000 1 0.000 32 Buyeo VII 1 0.000 1 0.000 
Seosan III 1 0.000 1 0.000 33 Geumsan III 1 0.000 1 0.000 
Uljin 1 0.000 1 0.000 34 Buan III 1 0.007 1 0.000 
10 Cheongju 1 0.000 1 0.022 35 Imsil III 1 0.000 1 0.000 
11 Daejeon III 1 0.000 1 0.000 36 Jeongeup III 1 0.000 1 0.000 
12 Chupungyeong III 1 0.000 1 0.000 37 Jangheung VII 1 0.003 1 0.000 
13 Pohang III 1 0.000 1 0.000 38 Haenam VII 1 0.069 1 0.000 
14 Daegu 1 0.000 1 0.000 39 Goheung 1 0.000 1 0.049 
15 Ulsan IV 0.339 1 0.000 40 Yeongju III 1 0.000 1 0.000 
16 Gwangju 0.140 0.152 41 Mungyeong III 1 0.000 1 0.000 
17 Busan IV 0.138 1 0.000 42 Yeongdeok VI 1 0.012 1 0.000 
18 Tongyeong VII 1 0.000 1 0.000 43 Uiseong III 1 0.000 1 0.000 
19 Mokpo VII 1 0.000 1 0.000 44 Gumi III 1 0.000 1 0.000 
20 Yeosu III 1 0.000 1 0.000 45 Yeongcheon 1 0.000 1 0.000 
21 Wando VII 1 0.000 1 0.000 46 Hapcheon III 1 0.000 1 0.000 
22 Jeju III 1 0.000 1 0.000 47 Miryang VII 1 0.000 1 0.000 
23 Ganghwa III 1 0.000 1 0.000 48 Sancheong III 1 0.000 1 0.000 
24 Yangpyeong III 1 0.000 1 0.000 49 Geoje VII 1 0.000 1 0.000 
25 Incheon III 1 0.000 1 0.000 50 Namhae III 1 0.000 1 0.000 
Table 4

Summary of results of ‘Tier-2’ approach, the significance test for the difference in the distributions at the 5% significance level

ST. No.Station nameCategoryTwo-sample KS test
Two-sample AD test
ST. No.Station nameCategoryTwo-sample KS test
Two-sample AD test
Hp-valueHAD staticsHp-valueHAD statics
Sokcho 0.241 1.793 26 Inje 1 0.000 1 7.480 
Chuncheon III 1 0.004 1 5.880 27 Hongcheon III 1 0.009 1 4.592 
Gangneung III 0.095 1 3.096 28 Jecheon II 1 0.000 1 17.068 
Wonju III 1 0.000 1 9.410 29 Boeun VI 0.095 1.669 
Ulleungdo 1 0.000 1 21.647 30 Cheonan III 0.095 1 2.821 
Suwon III 1 0.001 1 6.470 31 Boryeong III 0.358 1.082 
Chungju III 0.155 1 2.832 32 Buyeo VII 0.358 2.162 
Seosan III 1 0.000 1 9.349 33 Geumsan III 1 0.017 1 3.633 
Uljin 1 0.000 1 15.074 34 Buan III 0.678 0.855 
10 Cheongju 1 0.000 1 8.982 35 Imsil III 1 0.017 1 6.512 
11 Daejeon III 0.095 2.255 36 Jeongeup III 1 0.000 1 10.186 
12 Chupungyeong III 0.241 2.112 37 Jangheung VII 1 0.000 1 5.302 
13 Pohang III 0.678 1.113 38 Haenam VII 1 0.017 1 4.641 
14 Daegu 1 0.000 1 10.066 39 Goheung 0.241 1.437 
15 Ulsan IV 0.678 0.659 40 Yeongju III 1 0.004 1 5.431 
16 Gwangju 0.954 0.334 41 Mungyeong III 1 0.000 1 10.623 
17 Busan IV 0.155 2.270 42 Yeongdeok VI 0.954 0.793 
18 Tongyeong VII 1 0.001 1 7.224 43 Uiseong III 0.508 1.458 
19 Mokpo VII 0.358 1.786 44 Gumi III 1 0.000 1 10.688 
20 Yeosu III 0.056 1 2.597 45 Yeongcheon 0.056 2.144 
21 Wando VII 1 0.032 1 4.305 46 Hapcheon III 0.056 1 4.148 
22 Jeju III 1 0.004 1 3.510 47 Miryang VII 0.095 1 2.693 
23 Ganghwa III 0.056 1 4.588 48 Sancheong III 1 0.001 1 8.966 
24 Yangpyeong III 1 0.000 1 10.071 49 Geoje VII 0.508 1.339 
25 Incheon III 1 0.000 1 9.888 50 Namhae III 0.056 2.446 
ST. No.Station nameCategoryTwo-sample KS test
Two-sample AD test
ST. No.Station nameCategoryTwo-sample KS test
Two-sample AD test
Hp-valueHAD staticsHp-valueHAD statics
Sokcho 0.241 1.793 26 Inje 1 0.000 1 7.480 
Chuncheon III 1 0.004 1 5.880 27 Hongcheon III 1 0.009 1 4.592 
Gangneung III 0.095 1 3.096 28 Jecheon II 1 0.000 1 17.068 
Wonju III 1 0.000 1 9.410 29 Boeun VI 0.095 1.669 
Ulleungdo 1 0.000 1 21.647 30 Cheonan III 0.095 1 2.821 
Suwon III 1 0.001 1 6.470 31 Boryeong III 0.358 1.082 
Chungju III 0.155 1 2.832 32 Buyeo VII 0.358 2.162 
Seosan III 1 0.000 1 9.349 33 Geumsan III 1 0.017 1 3.633 
Uljin 1 0.000 1 15.074 34 Buan III 0.678 0.855 
10 Cheongju 1 0.000 1 8.982 35 Imsil III 1 0.017 1 6.512 
11 Daejeon III 0.095 2.255 36 Jeongeup III 1 0.000 1 10.186 
12 Chupungyeong III 0.241 2.112 37 Jangheung VII 1 0.000 1 5.302 
13 Pohang III 0.678 1.113 38 Haenam VII 1 0.017 1 4.641 
14 Daegu 1 0.000 1 10.066 39 Goheung 0.241 1.437 
15 Ulsan IV 0.678 0.659 40 Yeongju III 1 0.004 1 5.431 
16 Gwangju 0.954 0.334 41 Mungyeong III 1 0.000 1 10.623 
17 Busan IV 0.155 2.270 42 Yeongdeok VI 0.954 0.793 
18 Tongyeong VII 1 0.001 1 7.224 43 Uiseong III 0.508 1.458 
19 Mokpo VII 0.358 1.786 44 Gumi III 1 0.000 1 10.688 
20 Yeosu III 0.056 1 2.597 45 Yeongcheon 0.056 2.144 
21 Wando VII 1 0.032 1 4.305 46 Hapcheon III 0.056 1 4.148 
22 Jeju III 1 0.004 1 3.510 47 Miryang VII 0.095 1 2.693 
23 Ganghwa III 0.056 1 4.588 48 Sancheong III 1 0.001 1 8.966 
24 Yangpyeong III 1 0.000 1 10.071 49 Geoje VII 0.508 1.339 
25 Incheon III 1 0.000 1 9.888 50 Namhae III 0.056 2.446 
Accurate design rainfall estimation is vital for the effective performance of water resource management. Design rainfall is defined as the rainfall intensity for a given duration and return period. The primary objective of frequency analysis is to estimate the exceedance probability (or non-exceedance probability) of rainfall intensity for each duration. For instance, a 100-year return period corresponds to an extreme rainfall event that an exceedance probability of 0.01 (or 1% chance) that the extreme event will exceed in one year. The proposed BMN-QR method allows the direct computation of conditional quantiles, which result in the estimation of time-varying quantiles in the nonstationary frequency analysis context. Moreover, a stationary frequency analysis based on the GEV distribution is performed to compare with design rainfalls obtained from the BMN-QR-based nonstationary frequency analysis. The stationary frequency analysis involves the fit of the GEV to the ADMRs, and the distribution is then used to provide an estimate of the exceedance probability of given rainfall intensity. The GEV probability distribution and its inverse function can be defined as follows: 
formula
(11)
 
formula
(12)
where y is the ADMRs and k, σ, and μ are the shape, scale, and location parameters, respectively.

Changes in design rainfall for a set of representative stations under the identified categories are illustrated in Figure 8. In that figure, we clearly identified the differences in design rainfall for the first year (1973) and the last year (2015). The design rainfalls based on the stationary rainfall frequency analysis using the GEV distribution are mostly comparable to the estimates of design rainfall for the first year.

It has been acknowledged that the trends obtained from quantile regression for intermediate quantiles are more reliable against outliers than those of the ordinary linear regression for trends in the mean (Koenker & Regression 2005; Donner et al. 2012). However, the advantage decreases as regression quantiles increase due to reduced sample size, leading to an increase of uncertainty in the estimated trend. In other words, the estimated conditional distributions on the upper quantiles may be susceptible to single peaks (or outliers). In other words, uncertainty in trends can be affected by the location and magnitude of the peak, as illustrated in Figure 6. Furthermore, we experiment with different locations and magnitudes of the peak to explore changes in the uncertainty, as shown in Supplementary material, Figure S1. The obtained results are within our expectations with respect to the association between locations and magnitudes of the peak. The results demonstrated that uncertainty in trends is substantially narrower when maxima of precipitation extremes occur at the beginning or end of the period. Moreover, uncertainty in trends becomes more substantial as the magnitude of single peak increases at the same location. To consider this issue, we consider the entire posterior distributions for classifying possible categories associated with distributional changes in the ADMRs.

Discussion of nonstationarity in the ADMRs

Assessing distributional changes in the ADMRs over the various quantiles contributes to identifying the trends involved in different aspects of water resource management. Moreover, information on the distributional changes for a specific quantile supports estimating design rainfalls due to the fact that appropriate design rainfalls required for various water users differ with regard to quantity and frequency. Specifically, understanding the lower tails of the rainfall distribution function is more beneficial for irrigated agriculture, hydro-power generation, and water supply users. In contrast, the information obtained from the middle part of the rainfall distribution is useful for urban drainage design practice.

As shown in Figure 6, a few stations showed divergent trends in the ADMRs over the quantiles, while similar trends were identified for all quantiles at the remaining stations. For example, Incheon station showed decreasing trends in the lower quantiles but an increasing trend in the upper quantiles. For Buan and Cheonan stations, the means of the ADMRs were more or less equal, but the degree of dispersion varied significantly over the quantiles of interest. In this context, understanding the trend for different quantile levels can provide an effective strategy for mitigating potential damage from hydrological extremes. Some stations showed substantial interannual variability in the ADMRs, while their mean and variance were largely stationary over time.

On the other hand, Chuncheon, Wonju, Seosan, and Incheon stations showed noticeable differences in the distribution of ADMRs over different quantiles. For the significance of trends, few stations showed relatively narrow uncertainty bounds for all quantiles (Suwon, Daegu, and Jecheon stations), meaning that the identified trends over entire quantiles are more statistically significant. In contrast, other stations (Cheongju and Yeongdeok) exhibited a broader uncertainty. This study further utilized the uncertainty to assess changes in design rainfalls for a given specific category of nonstationarity.

The return period is an estimate of the likelihood that an extreme event of any given magnitude occurs in any given year. It is a statistical measure typically based on historical data denoting the average recurrence interval over an extended period of time and is usually used for risk analysis. In this sense, estimating the return period of an extreme rainfall event using BMN-QR suggests the collective behavior of the ADMR quantiles in a changing climate. Furthermore, it captures the key attributes of the dynamic distributional changes in the ADMR. Seven stations were selected to represent the different distributional changes for each category in order to demonstrate the behavior of trend and design rainfall: Daegu (Category I), Jecheon (Category II), Suwon (Category III), Ulsan (Category IV), Gwangju (Category V), Boeun (Category VI), and Wando (Category VII).

For Daegu station (Category I), a mild increasing trend (or no trend) is observed in the upper quantiles, while a steep increasing trend is found in the lower and middle quantiles, as illustrated in Figure 7(a). Furthermore, the lower uncertainty of the increasing trend associated with the middle quantiles clearly appears to be related to significant changes in the posterior ADMRs. In addition to the mild increasing trends in the upper quantiles, trends in the ADMRs can be collectively seen as strong upward-convergent and are also suggestive of a possible climate transition. The marginal posterior distributions of quantiles vary over time, which is closely related to the changes in location and scale parameters in the distribution of ADMRs. In 1973, there was a 1-in-5 year return period of a less than 150 mm/day ADMR, while the occurrence of the same amount of rainfall in 2015 has a higher frequency and thus a lower return period. On the other hand, it can be seen that the design rainfalls become almost identical, with a slightly increased dispersion, after the 1-in-10 return period.

Figure 7

Distributional changes in the ADMRs compared with design rainfall estimates of frequency analysis. The left panel represents the time series of ADMR and its quantile with Bayesian credible bounds. The right panel represents nonstationary design rainfall estimates corresponding to different return periods with those of the stationary frequency approach. The stations are: (a) Daegu, (b) Jecheon, (c) Suwon, (d) Ulsan, (e) Gwangju, (f) Boeun, and (g) Wando, which are selected from different distributional change categories. In the right panel, ‘UB’ stands for uncertainty bound and ‘FA’ stands for frequency analysis.

Figure 7

Distributional changes in the ADMRs compared with design rainfall estimates of frequency analysis. The left panel represents the time series of ADMR and its quantile with Bayesian credible bounds. The right panel represents nonstationary design rainfall estimates corresponding to different return periods with those of the stationary frequency approach. The stations are: (a) Daegu, (b) Jecheon, (c) Suwon, (d) Ulsan, (e) Gwangju, (f) Boeun, and (g) Wando, which are selected from different distributional change categories. In the right panel, ‘UB’ stands for uncertainty bound and ‘FA’ stands for frequency analysis.

For Jecheon station (Category II), there is an overall increase over all quantiles in AMDRs (Figure 7(b)). The marginal posterior distribution of quantiles appears to be largely shifted upward in 2015 compared with 1973, resulting in a comparable increase in design rainfall at all return periods. Hence, it is quite evident that a stationary estimate of design rainfall for the Category II is placed in the middle of the nonstationary design rainfalls between the first and the last year.

For Suwon station (Category III), all quantiles showed a marked upward trend in ADMRs during 1973–2015 (Figure 7(c)). The presence of exceptionally high rainfall in the early 2000s is leading to an enhanced trend in the current decade. As expected, a nonstationary estimate of design rainfall for 2015 is noticeably higher than that of the stationary frequency model. Hence, design rainfall estimates for the Category III stations are expected to be significantly enhanced compared with those of the stationary model. Under this circumstance, formal risk analysis, including the nonstationarity issue, needs to be conducted in a changing climate. Subsequently, several strategies need to be explored and implemented for appropriate management of the hydrologic risk of the existing water-related structures during their assigned lifespans.

For Ulsan station (Category IV), there is a significant decrease in variance over time, whereas there is no significant change in the mean of the ADMRs (Figure 7(d)). Therefore, the design rainfall estimates for 2015 are likely to be similar to those obtained in the stationary or less than those for the first year.

For Gwangju station (Category V), the same design rainfalls between for the first year (1973) and the final year (2015) indicate that there are no changes in the scale and location of the PDFs, as illustrated in Figure 7(e). Moreover, the design rainfall estimates for Category V are mostly fallen within the 95% credible intervals for all return periods.

For Boeun station (Category VI), there is a substantial increase in variance over time, which could lead to an increased tendency for design rainfall estimates to respond to the asymmetric trend of changes in ADMRs (Figure 7(f)). Therefore, the design rainfall estimates for Category VI under the nonstationary assumption are likely to be significantly higher than those of the stationary frequency approach, but the degree of increase varies less markedly than in Category III. Again, similar to the findings for Category III, risk-informed management strategies are generally required to ensure reliable operations beyond the planning horizon.

For Wando station (Category VII), there was a substantial decreasing trend in the ADMRs on average. Their associated design rainfall estimates for 2015 are significantly smaller than that estimated for 1973, as illustrated in Figure 7(g). Moreover, a noticeable decreasing trend was seen in the upper quantile with a strong negative slope; hence, a comparable decrease in the design rainfall estimates over 20–100-year return periods was identified.

Finally, the spatial distribution of the identified categories associated with the distributional changes in the ADMRs is illustrated in Figure 8. The Category III stations are extensively distributed over the western and central parts of South Korea. Interestingly, the stations with a downward-convergent pattern (Category VII) are largely distributed on the southern Korean coast, but the changes in the quantiles are mostly not significant. On the other hand, Category I stations are primarily located in the eastern part of South Korea and mainly show a statistically significant upward-convergent trend.

Figure 8

Spatial distribution of the categories associated with the distributional changes in ADMRs.

Figure 8

Spatial distribution of the categories associated with the distributional changes in ADMRs.

SUMMARY AND CONCLUSION

This study proposed a BMN-QR model to assess nonstationarity by exploring distributional changes in the ADMRs over the period of 1973–2015 in South Korea. To be more specific, the distributional changes were categorized into a number of classes in which nonstationarity was largely defined by changes in location and scale parameters of the probability distribution. The BMN-QR model-based nonstationarity detection scheme was then utilized to understand the impact of distributional changes on design rainfall estimates. The key findings of this study are summarized as follows.

  1. The proposed two-tier diagnostic approach based on the BMN-QR model provided a reliable and efficient way to identify distributional changes by simultaneously fitting models over the quantiles to be consistent with the monotonic response desired under the non-crossing constraints within a fully Bayesian framework.

  2. We identified seven categories for a distributional change in the selected stations. Most of the stations (28 of 50) were classified as Category III, which is characterized by an upward-divergent pattern in the distribution. Several of the stations classified as Category VII were distributed in the southern part of South Korea, in which there was a downward-convergent pattern over the period. However, the differences between the distributions of ADMRs over the quantiles were found not to be statistically significant, meaning that the differences in most of the stations classified as Category VII could be attributed to a random sample error rather than a systematic regime change.

  3. It was confirmed that the BMN-QR-based nonstationary frequency analysis model provided a method for understanding the key attributes of the dynamic distributional changes in the ADMR in a changing climate. In order to focus more clearly on nonstationarity and its direct impact on design rainfall estimates, this study explored changes in design rainfall estimates for different nonstationarity categories. For Categories I, II, III, and VI, a noticeable increase in design rainfalls was observed. In contrast, Categories IV, V, and VII showed no evidence of association with the risk of increased extreme rainfalls.

  4. It should be noted that the categorization process was still adapted in an ad hoc manner according to the distribution change in the ADMRs, although the two-tier significance test was considered. In this regard, future work will focus on developing a more systematic approach that explicitly combines a likelihood function of the categorization process within a hierarchical Bayesian modeling framework. Moreover, the climate interpretation and diagnosis of the distribution changes in extreme rainfalls over South Korea will be further explored in future work.

  5. It should be noted that the categorization process was still adapted in an ad hoc manner according to the distribution change in the ADMRs although the two-tier significance test was considered. In this regard, future work will focus on developing a more systematic approach that explicitly combines a likelihood function of the categorization process within a hierarchical Bayesian modeling framework. Recent advances in data science can play a significant role in understanding key climatic dynamics, which can be involved in the classification of pairs of regions that show mutually interconnected climate patterns. Thus, the well-identified clusters may subsequently help us better identify and discover distinct patterns (or clusters) in climate, and insights from the classification can also assist in refining existing models to better reproduce hydrologic processes. Moreover, the climate interpretation and diagnosis of the distribution changes in extreme rainfalls over South Korea will be further explored in future work.

ACKNOWLEDGEMENTS

This work was supported by a Grant (127568) from the Water Management Research Program funded by the Ministry of the Environment of the Korean Government. Climate data are freely available from https://data.kma.go.kr. The data used in this study are also available upon request from the corresponding author via email (hkwon@sejong.ac.kr).

AUTHOR CONTRIBUTIONS

We jointly conceived the project and designed the study. S.U. and H.K. wrote the majority of the manuscript. B.K. and T.K. critically revised the manuscript. All authors reviewed the manuscript.

SUPPLEMENTARY MATERIAL

The Supplementary Material for this paper is available online at https://dx.doi.org/10.2166/nh.2020.003.

REFERENCES

REFERENCES
Alexander
L. V.
Wang
X. L.
Wan
H.
Trewin
B.
2011
Significant decline in storminess over southeast Australia since the late 19th century
.
Australian Meteorological and Oceanographic Journal
61
(
1
),
23
.
Chang
H.
Kwon
W.-T.
2007
Spatial variations of summer precipitation trends in South Korea, 1973–2005
.
Environmental Research Letters
2
(
4
),
045012
.
Clarke
H. G.
Smith
P. L.
Pitman
A. J.
2011
Regional signatures of future fire weather over eastern Australia from global climate models
.
International Journal of Wildland Fire
20
(
4
),
550
562
.
Dickey
D. A.
Fuller
W. A.
1981
Likelihood ratio statistics for autoregressive time series with a unit root
.
Econometrica: Journal of the Econometric Society
49
(
4
),
1057
1072
.
Dogulu
N.
López
P. L.
Solomatine
D. P.
Weerts
A.
Shrestha
D.
2015
Estimation of predictive hydrologic uncertainty using the quantile regression and UNEEC methods and their comparison on contrasting catchments
.
Hydrology and Earth System Sciences
19
(
7
),
3181
3201
.
Donat
M. G.
Lowry
A. L.
Alexander
L. V.
O'Gorman
P. A.
Maher
N.
2017
Addendum: more extreme precipitation in the world's dry and wet regions
.
Nature Climate Change
7
(
2
),
154
158
.
Donner
R. V.
Ehrcke
R.
Barbosa
S. M.
Wagner
J.
Donges
J. F.
Kurths
J.
2012
Spatial patterns of linear and nonparametric long-term trends in Baltic sea-level variability
.
Nonlinear Processes in Geophysics
19
(
1
),
95
111
.
Dravitzki
S.
McGregor
J.
2011
Extreme precipitation of the Waikato region, New Zealand
.
International Journal of Climatology
31
(
12
),
1803
1812
.
Easterling
D. R.
Meehl
G. A.
Parmesan
C.
Changnon
S. A.
Karl
T. R.
Mearns
L. O.
2000
Climate extremes: observations, modeling, and impacts
.
Science
289
(
5487
),
2068
2074
.
Fischer
E.
Knutti
R.
2014
Detection of spatially aggregated changes in temperature and precipitation extremes
.
Geophysical Research Letters
41
(
2
),
547
554
.
Fischer
E.
Knutti
R.
2016
Observed heavy precipitation increase confirms theory and early models
.
Nature Climate Change
6
(
11
),
986
991
.
Gelman
A.
Carlin
J. B.
Stern
H. S.
Rubin
D. B.
2014
Bayesian Data Analysis, 2
.
Chapman & Hall/CRC
,
Boca Raton, FL
,
USA
.
Gilbert
R. O.
1987
Statistical Methods for Environmental Pollution Monitoring
.
John Wiley & Sons
,
New Jersey
.
Giorgi
F.
Im
E. S.
Coppola
E.
Diffenbaugh
N. S.
Gao
X. J.
Mariotti
L.
Shi
Y.
2011
Higher hydroclimatic intensity with global warming
.
Journal of Climate
24
(
20
),
5309
5324
.
Jung
I. W.
Bae
D. H.
Kim
G.
2011
Recent trends of mean and extreme precipitation in Korea
.
International Journal of Climatology
31
(
3
),
359
370
.
Kendall
M.
1975
Multivariate Analysis
.
Charles Griffin & Co
,
London
.
King
A. D.
Karoly
D. J.
Henley
B. J.
2017
Australian climate extremes at 1.5 C and 2 C of global warming
.
Nature Climate Change
7
(
6
),
412
416
.
Koenker
R.
Bassett
G.
Jr
1978
Regression quantiles
.
Econometrica: Journal of the Econometric Society
41
(
1
),
33
50
.
Koenker
R.
Regression
Q.
2005
Quantile Regression
.
Cambridge University Press
,
New York
.
Kunkel
K. E.
Pielke
R. A.
Jr
Changnon
S. A.
1999
Temporal fluctuations in weather and climate extremes that cause economic and human health impacts: a review
.
Bulletin of the American Meteorological Society
80
(
6
),
1077
1098
.
Kwon
H. H.
Khalil
A. F.
Siegfried
T.
2008b
Analysis of extreme summer rainfall using climate teleconnections and typhoon characteristics in South Korea
.
JAWRA Journal of the American Water Resources Association
44
(
2
),
436
448
.
Malik
N.
Bookhagen
B.
Mucha
P. J.
2016
Spatiotemporal patterns and trends of Indian monsoonal rainfall extremes
.
Geophysical Research Letters
43
(
4
),
1710
1717
.
Mann
H. B.
1945
Nonparametric tests against trend
.
Econometrica: Journal of the Econometric Society
13
(
3
),
245
259
.
Matti
C.
Pauling
A.
Küttel
M.
Wanner
H.
2009
Winter precipitation trends for two selected European regions over the last 500 years and their possible dynamical background
.
Theoretical and Applied Climatology
95
(
1–2
),
9
26
.
Mazvimavi
D.
2010
Investigating changes over time of annual rainfall in Zimbabwe
.
Hydrology and Earth System Sciences
14
(
12
),
2671
2679
.
Muhlbauer
A.
Spichtinger
P.
Lohmann
U.
2009
Application and comparison of robust linear regression methods for trend estimation
.
Journal of Applied Meteorology and Climatology
48
(
9
),
1961
1970
.
O'Gorman
P. A.
2015
Precipitation extremes under climate change
.
Current Climate Change Reports
1
(
2
),
49
59
.
Rana
A.
Uvo
C. B.
Bengtsson
L.
Sarthi
P. P.
2012
Trend analysis for rainfall in Delhi and Mumbai, India
.
Climate Dynamics
38
(
1–2
),
45
56
.
Razali
N. M.
Wah
Y. B.
2011
Power comparisons of Shapiro-Wilk, Kolmogorov-Smirnov, Lilliefors and Anderson-Darling tests
.
Journal of Statistical Modeling and Analytics
2
(
1
),
21
33
.
Río
S. d.
Herrero
L.
Fraile
R.
Penas
A.
2011
Spatial distribution of recent rainfall trends in Spain (1961–2006)
.
International Journal of Climatology
31
(
5
),
656
667
.
Roderick
T. P.
Wasko
C.
Sharma
A.
2019
Atmospheric moisture measurements explain increases in tropical rainfall extremes
.
Geophysical Research Letters
46
(
3
),
1375
1382
.
Roscoe
K. L.
Weerts
A. H.
Schroevers
M.
2012
Estimation of the uncertainty in water level forecasts at ungauged river locations using quantile regression
.
International Journal of River Basin Management
10
(
4
),
383
394
.
Roth
M.
Buishand
T.
Jongbloed
G.
2015
Trends in moderate rainfall extremes: a regional monotone regression approach
.
Journal of Climate
28
(
22
),
8760
8769
.
Sankarasubramanian
A.
Lall
U.
2003
Flood quantiles in a changing climate: seasonal forecasts and causal relations
.
Water Resources Research
39
(
5
),
1
10
.
Schmidli
J.
Frei
C.
Vidale
P. L.
2006
Downscaling from GCM precipitation: a benchmark for dynamical and statistical downscaling methods
.
International Journal of Climatology
26
(
5
),
679
689
.
Shiau
J.-T.
Huang
W.-H.
2015
Detecting distributional changes of annual rainfall indices in Taiwan using quantile regression
.
Journal of Hydro-Environment Research
9
(
3
),
368
380
.
Timofeev
A.
Sterin
A.
2010
Using the quantile regression method to analyze changes in climate characteristics
.
Russian Meteorology and Hydrology
35
(
5
),
310
319
.
Trenberth
K. E.
2011
Changes in precipitation with climate change
.
Climate Research
47
(
1–2
),
123
138
.
Trenberth
K. E.
Dai
A.
Rasmussen
R. M.
Parsons
D. B.
2003
The changing character of precipitation
.
Bulletin of the American Meteorological Society
84
(
9
),
1205
1218
.
Villarini
G.
Smith
J. A.
Baeck
M. L.
Vitolo
R.
Stephenson
D. B.
Krajewski
W. F.
2011
On the frequency of heavy rainfall for the Midwest of the United States
.
Journal of Hydrology
400
(
1–2
),
103
120
.
Westra
S.
Alexander
L. V.
Zwiers
F. W.
2013
Global increasing trends in annual maximum daily precipitation
.
Journal of Climate
26
(
11
),
3904
3918
.
doi:10.1175/jcli-d-12-00502.1
.
Yu
K.
Moyeed
R. A.
2001
Bayesian quantile regression
.
Statistics & Probability Letters
54
(
4
),
437
447
.
Yue
Y. R.
Rue
H.
2011
Bayesian inference for additive mixed quantile regression models
.
Computational Statistics & Data Analysis
55
(
1
),
84
96
.
Zhang
X.
Wan
H.
Zwiers
F. W.
Hegerl
G. C.
Min
S. K.
2013
Attributing intensification of precipitation extremes to human influence
.
Geophysical Research Letters
40
(
19
),
5252
5257
.

Supplementary data