Integrating nonstationary behaviors of typhoon and non-typhoon extreme rainfall events in East Asia

Extreme rainfall events in East Asia can be derived from the two subcomponents of tropical cyclones (TC) and non-TC based rainfall (mostly summer monsoons). Critical natural hazards including floods and landslides occur repeatedly due to the heavy rainfall associated with the two subcomponents, and disaster losses are increasing because global warming has caused changes in the extreme rainfall characteristics of two subcomponents. Subsequently, the frequency and intensity of extreme rainfall have reportedly become nonstationary. The majority of literature on nonstationary frequency analyses do not account for the different behaviors (stationarity or nonstationarity) of annual maximum rainfall (AMR) from the two subcomponents (PMTC and PMNTC). To carry out a nonstationary frequency analysis considering the different behaviors of the PMTC and PMNTC series, this study proposes a novel approach of integrating the fitted PMTC and PMNTC series after modeling the nonstationarity of the PMTC and PMNTC series individually. The presented results conclude that the proposed approach provides more reliable estimates than existing nonstationary approaches by reflecting the different features of the PMTC and PMNTC series. We suggest that the proposed approach provides a reasonable design rainfall in constructing hydraulics to mitigate the different nonstationary effects of two TC and non-TC rainfall extremes.

over Japan and the Korean peninsula and concluded an increase of rainfall intensity using future climate simulations (A1B emission scenario) by Couple Model Intercomparison Project 3 (CMIP 3) MME.
To mitigate the impact of extreme TC and non-TC rainfall events, flood control systems (e.g., hydraulic structures) have been installed with design rainfall and flood estimates. Hydraulic structures are generally designed based on the AMR series with the stationary assumption [9][10][11] . However, the assumption of stationarity in frequency analysis is questionable, and new frequency analysis methods that allow for nonstationarity in a given distribution parameters are required [12][13][14][15] . To treat nonstationarity of extreme events (e.g., rainfall and flood), several frequency analyses have been published in the literature, in which the parameters of a given distribution may vary over time 13,[16][17][18][19][20][21][22][23][24][25][26][27][28] . Regardless, the AMR events obtained from different subcomponents, TC and non-TC rainfall, extreme rainfall events from TC and non-TC have not been studied separately in the majority of literature available on nonstationary frequency analyses 9,16,18,25,29 . Recent studies related with the nonstationary frequency analysis for South Korea does not take those differences into account 29,30 because the difference of two weather systems (TC and non-TC extreme rainfall events) has not been further studied and the current model capability does not allow combining different systems.
In cases in which one of the two subcomponents contains an increasing trend whereas the other shows a decreasing trend, the traditional nonstationary frequency analysis could lead to an unrealistic representation of the design rainfalls (or floods). More specifically, the combined weather system might present stationarity or ambiguous nonstationarity by cancelling each other out in the trends. In addition, the AMR series for the rainfall frequency analysis are routinely constructed from the largest rainfall totals corresponding to different durations in each year so that one of the subcomponents with an increasing trend might contribute to the AMR more in the future and vice versa. For these reasons, a traditional nonstationary approach may not be appropriate for such systematic nonstationary characteristics.
We propose a design rainfall estimation approach with a nonstationary frequency analysis that separately treats the regions (e.g., East Asia) that are influenced by clearly different subcomponents for a better understanding of future climate change. Specifically, the current study aims to develop a novel approach in (1) partitioning extreme rainfall data (PM TOT ) into TC (PM TC ) and non-TC (PM NTC ) rainfall with a recently developed technique by Son, et al. 31 [details in supplementary material]; (2) obtaining two AMR series from the separated TC and non-TC rainfall data; (3) modeling the PM TC and the PM NTC series with nonstationary GEV (generalized extreme value) model [details in Methods] separately; and (4) integrating the fitted two nonstationary GEV models to obtain a design rainfall estimate to adapt the nonstationarities in the PM TC and PM NTC series. Details of the overall procedure are provided in the 'Methods' section. The proposed approach is then applied to every rainfall station in South Korea and 10 rainfall stations in the Tokyo region that are usually in the path of TC in Japan.

Results
Among the East Asian countries, South Korea presents a serious seasonal rainfall deviation in that most of the annual average rainfall (approximately 66%, 854.4 out of 1,292 mm) occurs during the rainy season (June-September), which makes it vulnerable to flooding 1 . In addition, this region has experienced extreme rainfall events causing water-related disasters from the two subcomponents, TC and non-TC rainfall. Therefore, South Korea was selected as the study area. Additionally, we included the Tokyo region as another study area because numerous TCs in the WNP have passed through this region (Supplementary Table S3). Figure 1 shows the geographical distribution of a vertically integrated moisture flux (VIMF) [details in Methods] 32 related to a TC rainfall event that occurred in September 2003 ( Fig. 1(a)) and a non-TC rainfall event that occurred in June 2005 ( Fig. 1(b)) over East Asia. Moreover, TCs commonly generate over the WNP in conditions of warm temperatures (above 26.5 °C), high relative humidity, and low vertical wind shear 33,34 . After their formation, TCs accompanied by strong winds and high moisture content tend to move across the WNP Ocean and approach East Asia in an arched pattern along the edges of the North Pacific high air mass.
For example, the TC Maemi that occurred on September [10][11][12]2003, is presented in the top panels of Fig. 1. The TC originated over the WNP and headed northwestward toward Taiwan with intensifying magnitude on September 10, 2003, and hesitated near Taiwan (its center was located at 125°E, 25°N) as shown in Fig. 1(a-1). Afterward, the TC recurved northeastward toward the Korea Strait and the southern Sea of Japan on September 11 (see Fig. 1(a-2)). Eventually, the TC emerged on the Korean peninsula on September 12 (see Fig. 1(a-3)) and poured out tremendous amounts of rainfall, causing severe destruction of hydraulic structures and economic losses as well as 130 casualties.
A non-TC event (a typical case of monsoon rainfall, also called Changma in Korea) is a part of the quasi-stationary front of the East Asia summer monsoon season. For example, the non-TC event that occurred on June 27-29, 2005, is presented in the bottom panels of Fig. 1(b-1~b-3). The accumulation of moisture was maximized in the East China Sea and the WNP Ocean and traveled northeastward near the east side of Taiwan ( Fig. 1(b-1)). The northeastward movement of the moisture mass was established by the development of the subtropical ridge located at the bottom right side of the moisture mass. Hence, the moisture was concentrated in the midlatitude regions (30-60°N) of the Northern Hemisphere including South Korea and Japan because of the subtropical environment (approximately 20-30°N, 130-160°E) and an Okhotsk high (approximately 40-50°N, 140-160°E) located in the southeast and northeast of the Korea peninsula, respectively (called Changma front). As shown in Fig. 1(b-1) and (b-3), a large mass of moisture constantly traveled toward the region containing Taiwan, South Korea, and Japan and persisted through July 13, 2005 (approximately three weeks). This monsoon rainfall (a typical severe event of a non-TC rainfall) was an intense and prolonged event triggering landslides and floods with 345 casualties as well as severe social and economic losses.
As Fig. 1 indicates, the characteristics of TC and non-TC rainfall events present a clear difference in terms of their intensity and duration. Therefore, it is apparent that the future trends of TC and non-TC rainfall events are possibly different, and an extreme analysis associated with them must be performed individually.
Scientific RepoRts | 7: 5097 | DOI:10.1038/s41598-017-04629-1 Therefore, a non-parametric Mann-Kendall trend test 35 was used to assess the statistical significance of trends in different AMR series such as PM TOT , PM TC , and PM NTC . Figure 2(b) and (c) present the results of the Mann-Kendall test with a 10% significance level (p ≤ 0.1) for the PM TOT , PM TC , and PM NTC series affecting South Korea and Tokyo. In South Korea, the PM TOT series (the left side of Fig. 2(b)) tends to increase significantly over a part of South Korea including the Han River and Nakdong River basin (approximately 30% of the stations). For the PM TC (the middle map of Fig. 2(b)), a significant and increasing trend can be observed in a part of the Nakdong River basin (approximately 10% of the stations). The PM TC series exhibits stationary behavior for most areas of South Korea except for the Nakdong River basin.
In contrast, the PM NTC series shows a clearly increasing trend in the Han River and Nakdong River basin (38% of the stations). In particular, the PM NTC series shows a similar increasing trend as that observed in the PM TOT series in the northwestern area of South Korea such as the Han River and Geum River basins. This result indicates that the PM TC events have little effect on extreme rainfall; however, the PM NTC series contributes to most of the AMR series (PM TOT ) in the northwestern areas of South Korea (Han River: 87.3%; Geum River: 89.0%, the ratio of the PM NTC in the PM TOT series). In a part of the southern area of South Korea, the PM TOT shows different tendencies in terms of sign (increase or decrease) and strength compared with those of the PM TC and PM NTC because the trends may be offset by the behaviors of the two subcomponents (TC and non-TC). Therefore, the trend assessment with PM TOT may result in an unreliable future projection in the region where the AMR events (PM TOT ) comprise the two subcomponents.
In the Tokyo region, the PM TOT and PM TC series show no significant trends even though an increasing trend is seen in the PM TOT series. The PM TC series does not show a significant trend over the region. However, the PM NTC series increases significantly in the middle part of the Tokyo region (40% of the stations). Nonetheless, the PM NTC series shows a significant trend, and this result suggests that the PM TOT series for design rainfall can be treated as stationary. Therefore, the AMR series (PM TC and PM NTC ) must be treated differently when nonstationary processes are required and when extreme rainfall occurs due to two or more different weather systems.
The key features of the AMR series (PM TC and PM NTC ) over the study area that are identified in this study are as follows: First, the PM TC series exhibits nearly stationary behavior, while the PM NTC series appears to be nonstationary, showing noticeable trends. Second, although evidence for nonstationarity in the PM NTC series is clear, there is, however, evidence of stationarity in the PM TOT series. This may be interpreted as representing a stationary behavior buried in a nonstationary trend by compensating for the stationarity in the PM TC series.
For evaluating the suitability of five different probability distributions as normal, Gamma, logistic, Gumbel, and GEV distributions, the Akaike Information Criterion (AIC) 36 , Bayesian Information Criterion (BIC) 37 , and Kolmogorov-Smirnov (K-S) test 38,39 can be used. Calculating the AIC and BIC for five distributions, we found that the Gamma, gumbel, and GEV distributions are selected for a number of rainfall stations (Supplementary   S2). In the two study areas, the frequency analysis of extreme rainfall events often apply the Gumbel and GEV distributions in literature. In addition, Wi, et al. 29 suggested that the GEV distribution appropriately represents the AMP series for South Korea. In particular, the GEV distribution is more flexible relative to the Gumbel distribution because of using a larger number of parameters 40 . Overall, K-S test justified the use of the GEV distribution as an appropriate alternative for all stations in the two study areas. Therefore, the GEV distribution is used to fit the AMR series in the current study. Figure 3 presents the time series of the AMR series (PM TOT , PM TC , and PM NTC ) with linear slopes (left panels) and the rainfall quantile estimates with 95% confidence intervals (right panels) for the selected and stationary GEV models using the total rainfall (GEV TOT model) and integrated rainfall (IGEV model) for the selected sites of Busan (a), Hachioji (b), and Seogwipo (c). Details of the stationary or nonstationary GEV models are described in the 'Methods' section. To estimate confidence intervals for the parameters (Supplementary Fig. S3) and quantiles (Fig. 3) of selected stationary and nonstationary GEV models using the total rainfall and integrated rainfall, the bootstrap method was used. In the current study, 1000 bootstrap samples were generated from fitted selected GEV models. For example, in Fig. 3(a-2), the stationary model is used to fit the PM TOT series of Busan (station no. 22, see Fig. 2(a)) since the Mann-Kendall test shown in Fig. 2(b) has no significant trends. Meanwhile, the separate consideration of the PM TC and the PM NTC series by integrating the two series through equation (7) leads to the nonstationary model of IGEV 10 00 for reflecting the identified trends in the PM TC (i.e., stationarity) and the PM NTC (i.e., nonstationarity), as shown in Fig. 3(a-1).
The stationary model (GEV TOT 00 ) for the PM TOT series of the Busan and Hachioji stations is illustrated in Fig. 3(a-1) and (b-1). Note that one of the PM TC and PM NTC values is selected over time for each value of the PM TOT series such that the PM TOT value overlaps with one of the PM TC or the PM NTC values, as shown in Fig. 3(a1). The PM TOT series of the Busan (a-1) and Hachioji (b-1) stations show that the PM TC are generally selected for the PM TOT series during the early time periods (i.e., 1973-1995), while for the latter periods (i.e., 1996-2012), the PM NTC is more prevalent. As aforementioned, this is mainly due to different levels of nonstationarity over time in the two time series. For these reasons, the estimated design rainfalls may not be reliable as shown in the right panels of Fig. 3(a-2) and (b-2) without treating the AMR series separately. The IGEV model shows higher values than the stationary model for all the return periods, and more specifically, higher return periods shows higher differences. In this setting, the stationary model for critical hydraulic structures (e.g., dams and large river levees) associated with high return periods can lead to an underestimation of hydrologic risk.
As represented in Fig. 3(c), there are nonstationary behaviors of the PM TOT and PM NTC series, while the PM TC series presents stationary behavior for Seogwipo (i.e., station no. 29). Although the y-intercept of the linear fit from the PM TOT series is higher than that of the PM NTC , both have similar slopes to the PM TOT : 1.688/year and PM NTC : 1.657/year, which are statistically significant. Figure 3(c-2) shows the differences between quantiles of the GEV TOT 10 model and IGEV 10 00 model. For higher return periods, the quantiles of the IGEV 10 00 model tend to be larger   -1,b-1 and c-1). The straight lines indicate the linear fit to annual maximum daily rainfall. The right panels (a-2,b-2 and c-2) represent the quantiles with 95% confidence intervals estimated by the nonstationary GEV models between total rainfall (green dotted line with squares; (a-2 and b-2) are GEV TOT 00 models, and (c-2) is GEV TOT 10 ) and integrated rainfall (red dotted line with circles; IGEV 10 00 ); the gray dotteddashed lines with upward-pointing triangles in the c-2 indicate the quantiles estimated by the stationary GEV model (GEV TOT 00 ). Note that (1) the confidence intervals were estimated with the bootstrapping method with 1000 resampling; (2) there are some cases that TC rainfall does not occur in a certain year (for example, seogwipo station at 1996, 2005 at (c-1)) and the occurrence probability is adopted in these cases as shown in Eqs (6) and (7) with p o to take the non-occurrence of TC rainfall into account for quantile estimation. In other words, higher scale parameter ( Supplementary Fig. S3(b-2) and (c-2)) obtained from PM TC series may lead to a relatively higher quantile for higher return periods and vice versa while smaller location parameter of non-TC rainfall than total rainfall as shown in Fig. 3(c-2) may lead to relatively smaller quantiles during lower return periods. In terms of safety, it means that the proposed approach could be better than the existing nonstationary approach at estimating design rainfalls for dam structures with the high return periods.
To quantitatively identify the spatial features of design rainfalls with different return levels, the relative percent differences (RPD) were calculated as the differences of IGEV − GEV TOT divided by GEV TOT model and expressed as a percentage. Figure 4 presents the RPD of the quantiles estimated by the GEV models for the three different return levels (i.e., 10, 50, and 100 years). In South Korea, the quantiles corresponding to the 10-and 50-yr return periods estimated by the IGEV model are higher than those estimated by the GEV TOT model especially in the southeast region of South Korea ( Fig. 4(a-1) and (a-2)). In Tokyo, the quantiles estimated by the IGEV model shows higher values than those estimated by the GEV TOT model.
In general, the northern area of South Korea tends to be affected by extreme rainfall from only non-TC rainfall because numerous TCs move toward the southeast region of South Korea and the southern Sea of Japan. In the northwestern area including the Han River basin and Geum River basin, the quantiles estimated by the IGEV model are quantitatively similar to those estimated by the GEV TOT model. In contrast, in the southern area of South Korea and Tokyo, which are simultaneously influenced by the two subcomponents, the quantiles estimated by the IGEV and GEV TOT models show clear differences. These results indicate that the proposed approach can estimate relatively well by capturing nonstationary behaviors of non-TC rainfall better.

Discussion
In this study, we proposed a novel method of integrating the fitted AMR series of TC and non-TC rainfall by nonstationary models. In South Korea, the PM TOT and PM NTC series exhibited statistically significant increases over the Han River basin and the western Nakdong River basin. However, no significant trends were found in PM TC series, except for a statistically significant increasing trend at six stations in the Nakdong River basin. In Tokyo, the PM TOT and PM TC series showed no significant trends either. In contrast, the PM NTC observed a statistically significant increase over the middle areas of Tokyo. It indicated that the future long-term evolutions of the two subcomponents might be different from each other. Therefore, TC and non-TC rainfall events must be treated differently when nonstationary frequency analyses are required and when extreme rainfall occurs due to two or more different weather systems.
To treat the nonstationarity of the AMR series differently, we carried out a nonstationary frequency analysis of the different AMR series from TC and non-TC events individually and integrated the fitted PM TC and PM NTC series with the developed IGEV model in the current study. We found different features over South Korea and Tokyo according to extreme events such as the PM TC in stationary behavior and the PM NTC in nonstationary behavior. However, it was not possible to consider the existing frequency analysis to be different in the long-term trends of the AMR series from the two subcomponents. Therefore, the existing stationary approach might have led to imprecise inferences because it ignored the nonstationarity of the non-TC rainfall. It was concluded that the proposed approach can be a useful alternative in designing hydraulic structures and can ensure safety against flood damage in a nonstationary regime of extreme rainfall events.
Only the time variable was considered as a covariate in mapping extreme rainfall to highlight the need to handle the nonstationarity of the AMR series from the two subcomponents separately. The results of the proposed approach were expected to provide a more reliable estimate for designing flood-related hydraulic structures considering climate change and the management of water resources. One can extend the current study to include the variability change, but sufficient number of AMR data must be obtained in advance. Sufficient evidences such that different weather systems lead to different behaviors for overall trends of AMR series must be provided to apply the proposed method.

Methods
Integrating procedure for TC and non-TC rainfall. A radius-based TC rainfall extraction method (RTREM) [details in supplementary material] was developed to differentiate series of rainfall into TC and non-TC rainfall. Among TC rainfall values, X t was the annual maximum value at year t, and for non-TC rainfall, Y t was used. The time index t was dropped for simplicity, and X and Y were used for the AMR series of TC and non-TC rainfall, respectively. Eventually, an AMR value for a certain year (Z int ) was obtained by selecting the maximum value of X and Y.
int Note that Z int indicated the integration of the two AMR series (X and Y). From a general statistical derivation, it was easily shown that Z int int Since it was generally agreed upon that the AMR of the TC (X) and non-TC (Y) rainfall events were independent, equation (3) Vertically Integrated Moisture Flux (VIMF). The VIMF presented in Fig. 1  where W u and W v denoted the average wind vector in the x-and y-directions for the range 1000 and 850 hPa, respectively, and SH was the integrated specific humidity for the same range. The VIMF provided succinct information on the atmospheric moisture cycle.
Procedure for the quantile estimation. A detailed description of the modeling procedure is as follows ( Supplementary Fig. S1): (1) Extract TC rainfall events (P TC ) among the total rainfall events (P TOT ) every year with RTREM. Let the remainder be non-TC rainfall (P NTC ).