Changes in extreme rainfall and its implications for design rainfall using a Bayesian quantile regression approach

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.

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.

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⁠.

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.

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.

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.

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.

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.

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.

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.

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.

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.

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.

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.

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 ([email protected]).

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.

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