A nonparametric statistical framework using a kernel density estimator to approximate flood marginal distributions – a case study for the Kelantan River Basin in Malaysia

From the perspective of water resources operational planning in a river basin, the design and operation of flood management infrastructure often demands an accurate estimate of the flow exceedance probability for assessing hydrologic risk (i.e. Salvadori 2004; Xu et al. 2015). Flood frequency analysis (FFA) establishes an association between extreme event quantiles and their non-exceedance probabilities by fitting a univariate or multivariate probability distribution function (Cunnane 1988; Rao & Hameed 2000; Zhang 2005). A flood can be defined through three inter-correlated random vectors; that is, peak flow, flood volume and flood duration, which can lead univariate distribution analysis of return periods to underestimate or overestimate floods (Yue et al. 2001; Zhang & Singh 2006; Reddy & Ganguli 2012). This has led researchers to explore multivariate distribution analysis for the estimation of flood design quantiles under different notations of return periods based on joint distribution, conditional joint distribution or Kendall's distribution and the establishment of joint probability density functions (or PDFs) and joint cumulative distribution functions (or CDFs) for various possible combinations of the flood vectors (Yue 1999; Zhang & Singh 2006; Daneshkhan et al. 2016). Most of the earlier research applied traditional multivariate functions such as bivariate normal, gamma and lognormal distributions (i.e. Yue 1999, 2000, 2001). However, several modelling constraints and limitations have prompted the incorporation of multivariate copulas distribution analysis (i.e. De Michele & Salvadori 2003; Favre et al. 2004) which is frequently undertaken via a parametric framework, where both the marginal distribution of each targeted flood vector and their joint dependence structure are modeled under univariate parametric functions or copulas distribution (i.e. Favre et al. 2004; Zhang 2005; Veronika & Halmova 2014; Fan & Zheng 2016 and references therein). While some researchers have implemented copulas under semiparametric settings where marginal functions are approximated via nonparametric; that is, kernel density functions (KDE) or orthonormal series, their joint structures are still modelled under parametric settings (i.e. Karmakar & Simonovic 2008, 2009; Reddy & Ganguli 2012).

Selecting the most parsimonious univariate functions for defining flood marginal distributions is often a mandatory pre-requisite before establishing their joint dependence structure. Parametric functions have always imposed an assumption that the random samples are coming from known populations whose PDFs are pre-defined; that is, marginal distributions are assumed to follow a specific family of parametric probability functions. In reality, however, the best fitted marginal distributions may not be from the same probability distribution family (Adamowski 1985; Silverman 1986; Kim & Heo 2002; Botev et al. 2010). Dooge (1986) has already pointed out that no amount of statistical refinement can overcome the consequences of a lack of prior probability distribution information of the observed random samples. More especially, in the case of multimodal or skewed distributions, parametric distributions might be incompatible and lead to inconsistencies in the estimated quantiles. According to Bardsley (1988) and Bardsley & Manly (1987), the approximation of any distribution tail beyond the largest value under a parameter distribution can be a difficult task. In the last few decades, in the field of hydrologic or flood frequency analysis an attempt has been made to define bona fide density functions via kernel density estimators, or KDE, which are recognized as a flexible and very stable nonparametric data smoothing procedure for approximating or inferencing populations based on finite observational (Adamowski 1989, 2000; Guo 1991; Lall et al. 1993; Bowman & Azzalini 1997; Efromovich 1999; Kim et al. 2003; Ghosh & Mujumdar 2007). Alternate theoretical overviews for nonparametric setting have been described in the earlier literature such as Rosenblatt (1956), Parzen (1962) and Bartlett (1963). Nonparametric distributions do not require any prior distribution assumptions and have a higher level of flexibility in their univariate function in comparison with parametric density estimators (Adamowski 1989; Moon & Lall 1994).

Adamowski (1985) retrieved the flood frequencies curve using nonparametric kernel estimators while Adamowski & Labatiuk (1987) compared the difference between real and synthetic data (derived from a Monte Carlo procedure) and identified the modelling efficiency of nonparametric density simulations. Similarly, Adamowski (1989) performed comparative assessments between parametric functions (i.e. Weibull and log-Pearson type III) and nonparametric density functions based on the Monte Carlo simulation statistics and found nonparametric density functions demonstrations performed better than the parametric functions. Adamowski & Feluch (1990) were the first to identify the efficiency of kernel estimators through selecting an appropriate kernel function for the extrapolation of distribution samples beyond the available record. Adamowski (1996) applied the nonparametric procedure to drought modelling while Kim & Heo (2002) compared nonparametric distributions from parametric functions for annual maximum flood samples. Lall et al. (1993) postulated the different orientations required to select an appropriate shape representation for kernel functions and their bandwidth for different distribution scenarios such as symmetrical, asymmetrical or skewed and mixed data structures. Kim et al. (2003) estimated the return periods of low-flow or drought characteristics using a nonparametric distribution framework. Kim et al. (2006) established bivariate drought characteristics using the Palmer drought index in a nonparametric procedure. Karmakar & Simonovic (2008) approximated the flood marginal densities by introducing both parametric and nonparametric functions. In this experiment, nonparametric based kernel density functions as well as orthonormal series functions, which are demonstrated earlier in the literature; for example, Bowman & Azzalini (1997) and Efromovich (1999), were used to establish univariate density functions. Santhosh & Srinivas (2013) demonstrated a flood frequency relationship by employing a bivariate diffusion process based on the adaptive kernel functions; that is, D-kernel function, and compared this against parametric copulas functions. This demonstration revealed that the nonparametric methodology is a data-driven approach that poses the higher degree of flexibility when establishing probability density functions. Lall et al. (1996), Zhang & Karunamuni (1998), Kim et al. (2003), Karmakar & Simonovic (2008) and Santhosh & Srinivas (2013) all pointed out the flexibility of an interactive sets of kernel functions for modelling extreme hydrological episodes such as Gaussian or Normal, Epanechnikov, Triangular and Quadratic functions.

Floods are becoming the most challenging hydrologic issue in the Kelantan River basin in Malaysia, and particularly during the period of wet monsoons (DID 2004; MMD 2007). Recent extreme weather across the basin includes the intense and prolonged precipitation in the year 2002, which caused flooding and affected a total area of 1,640 km² with a total population of 714,287, while in the year 2014 the worst flood ever recorded in history affected more than 200,000 people in the several parts of the east coast of the Kelantan River basin. The river is about 248 km long and drains a catchment area of 13,100 km² The catchment is 150 km long and 140 km wide. During the wet season, between mid-October and mid-January, this basin receives a rainfall of about 2,500 mm. According to the studies of Hussain & Ismail (2013), the Gulliemard Bridge, Lebir and Galas gauge stations have higher flood frequency than the Nenggiri gauge station. In 2018, Alamgir et al. (2018) performed a multivariate analysis of floods under a parametric copulas distribution framework for the different gauge stations of the river basin, while Nashwan et al. (2018) performed flood susceptibility assessments at the different gauge stations, which revealed that the downstream area is under the highest risk of devastating floods. A number of studies have also identified the negative consequences of land use changes on the catchment's response (i.e. Wan 1996; Jamaliah 2007).

All three flood characteristics; that is, peak flow, flood volume and flood duration, are important when formulating actions and measures to manage flood risk. For example, when designing retention basins or the spillways for reservoirs or any other infrastructures designs where flood storage is involved, the estimation of flood volume is required as well as the peak discharge in order to calculate the impact of inflow on the storage (Gaál et al. 2015). Similarly, estimating the joint behaviour of peak flow-flood volume and flood volume-flood duration when assessing flood control or diversion options (i.e. Xu et al. 2015). Therefore, estimating the multivariate designs and their associated return periods is an essential element of making informed risk-based decisions in this river basin.

The case study is divided into two parts, with each part covered in a full paper. The objective of this part of the case study is to identify the best marginal distribution of the flood characteristics. In this paper, the efficacy of a kernel density estimator is tested by assessing the adequacy of an interactive set of kernel functions for capturing the flood marginal density. Demonstration of the efficacy of copulas functions for establishing joint distribution functions whose marginal distributions are approximated with nonparametric distribution or kernel density functions (which are also called semiparametric copulas distribution settings) is introduced separately in the companion paper. An at-site event-based or block (annual) maxima methodology based on 50 years (from 1961 to 2016) of daily stream flow data collected at Gulliemard Bridge gauge station in the Kelantan River basin is described. A brief mathematical overview of the nonparametric procedure in estimating univariate margins via kernel density function is provided below. The sampling procedure for estimating multivariate flood characteristics, estimating marginal distribution functions by employing a variety of univariate kernel functions and also via parametric functions are discussed in the third section of this paper.

Rosenblatt (1956) introduced the concept of kernel estimators through smoothing kernel weights on each of the random observations. The univariate kernel density estimation (or KDE) is a nonparametric approach to approximate the PDFs, say of a given random observations of or flood characteristics such as peak flow, flood volume or flood duration, where inference about the flood populations is made based on finite samples. Therefore, univariate kernel functions are used to estimate the probability density of the random observations having the following statistical property.

The efficiency of the estimated kernel density depends upon two factors: (1) an appropriate choice of the kernel bandwidth and (2) selection of the kernel function. Bandwidth estimation procedures are discussed below. The kernel functions listed in Table 1 are ideally suited to unbounded observations but in actuality most hydrologic data including rainfall, streamflow or humidity and so on, are either upper or lower bounded. This can lead to boundary leakage problems when applying standard kernel functions in the hydrological data domain. This issue has been tackled in a number of ways including using a Method of reflection (Silverman 1986), Methods of transformations (Marron & Ruppert 1994), Beta kernels (Chen 2000) and Adaptive kernels (Botev et al. 2010). Such boundary leakage issues during kernel estimation are beyond the scope of this paper.

Monsoonal floods seem to have increased in the Kelantan River basin in Malaysia in the last few decades in terms of both frequency and magnitude (DID 2004; MMD 2007). It is the longest river of Kelantan state, which rises in the Tahan mountain range and flows to the South China Sea in the north-eastern part of Peninsular Malaysia between the geographical location of Lat 4° 30′ N to 6° 15′ N and Long 101°E to 102° 45′ E. The Galas River and the Lebir River are the two major tributaries of the Kelantan River. The land use in the upper catchment is forest while agriculture including paddy farms, rubber and oil-palm plantations are the major land-use activities in the middle and lower areas of the catchment. The precipitation for this region typically varies between 0 mm (in the dry period) to 1,750 mm (in the wet or north-eastern monsoonal period) (DID 2004).

Tests for stationarity or the existence of serial correlation (or autocorrelation) within the flood series is often a pre-requisite before introducing random samples into a univariate or multivariate framework (Daneshkhan et al. 2016). Ljung & Box (1978) based hypothesis testing, which is also known as Q-statistics, were performed on each time series of observations. As indicated in Table 4(a), the tests found negligible or zero first-order autocorrelations for each of the flood vector series for different lag sizes (i.e. lag 20, lag 10, lag 5). A nonparametric rank-based Mann-Kendall (or M-K) test (Mann 1945; Kendall 1975) was also performed to test for the existence of any monotonic trend within the historical flood series. As indicated in Table 4(b), the test found zero monotonic trend at the 5% or 0.05 level of significance within the flood vector series. Testing for the existence of a homogenous environment between any two given time points was also investigated for each flood vector through application of a Pettit test (Pettitt 1979), a Buishand test (Buishand 1982), von Neumann's test (Jaiswal et al. 2015) and by undertaking Alexanderson's SNHT based hypothesis testing (Alexandersson 1986). As demonstrated in Table 4(c), it was found that homogeneity existed within the flood vector series. It was concluded therefore that the time series of the flood vector do not exhibit any significant trend.

Table 1 identified some standard univariate kernel functions, which are adopted when defining a best-fit flood marginal distribution. The bandwidth of the candidate kernel functions was estimated using the optimal bandwidth algorithm given in Equation (6) (Azzalini 1981; Silverman 1986). According to Kim et al. (2006), the nonparametric density approximations do not facilitate a closed form of the PDF and CDF, thus CDFs were estimated through an empirical procedure that is based on numerical integration (Kim & Heo 2002).

Some frequently applied parametric family functions such as the Log Pearson type III distribution (Bobee 1974), the Lognormal-2P function (Yue 2000), the Weibull-3P distribution (Heo et al. 2001), the Johnson SB-4P distribution (Keshtkaran & Torabihaghighi 2011), the Gamma-3P (Xu et al. 2015) and the Inverse Gaussian-3P functions (Daneshkhan et al. 2016) were also tested and compared to the nonparametric estimates. The vectors of unknown statistical parameters or model parameters were estimated using the Maximum likelihood estimators (MLE) and Methods of Moments estimators (MOM) and their estimated values are listed in Table 5.

Table 6 presents the estimated values of MSE, RMSE, AIC, BIC and the HQC statistics of all the nonparametric kernel distribution functions fitted to the flood vectors. It was found that the Triweight kernel outperformed other functions as it gave the lowest values of the fitness test statistics for all three flood distribution series; that is, MSE (0.00022), RMSE (0.01483), AIC (−419.0760), BIC (−417.164) and HQC (−418.348) for peak flow, MSE (0.00016), RMSE (0.01287), AIC (−433.27), BIC (−431.36) and HQC (−432.55) for flood volume series and MSE (0.00048), RMSE (0.02208), AIC (−379.27), BIC (−377.26) and HQC (−378.54) for flood duration. The performance of the Biweight and Triangular functions was also more effective than the targeted parametric functions. While the Epanechnikov kernel was less effective than the other candidate kernel functions, it was still better than the parametric functions as revealed in Table 6. Based on analytically based fitness measures, it was concluded that it is likely that Triweight kernel function is the best-fitted distribution for defining the marginal distribution of peak flows, flood volumes and flood durations in the Kelantan River basin.

A qualitative approach based on a graphically based visual inspection was also conducted for each flood vector for the probability density plot, the cumulative density plot, the probability- probability (or p-p) plot, the K-S test comparison cumulative fraction plot and the K-S test comparison percentile plot as illustrated in Figures 5 and 6(a)–6(c). It is noted that the Kolmogorov-Smirnov test (K-S) is a nonparametric distribution-free test that seeks to investigate the largest vertical gap between cumulative empirical and theoretical probabilities and also has the advantage of not assuming the distribution of data (Xu et al. 2015).

It was concluded that these plots clearly indicate the effectiveness of a nonparametric kernel structure and support the adoption of a Triweight kernel function for defining univariate flood marginal distributions.

Floods are becoming the most challenging hydrologic issue in the Kelantan River basin in Malaysia, and particularly during the period of wet monsoons. All three flood characteristics; that is, peak flow, flood volume and flood duration, are important when formulating actions and measures to manage flood risk. Therefore, estimating the multivariate designs and their associated return periods is an essential element of making informed risk-based decisions in this river basin.

In this paper, the efficacy of a kernel density estimator is tested by assessing the adequacy of an interactive set of kernel functions for capturing the flood marginal density of 50 years (from 1961 to 2016) of daily stream flow data collected at Gulliemard Bridge gauge station in the Kelantan River basin.

Tests for stationarity or existence of serial correlation (or autocorrelation) within the flood series is often a pre-requisite before introducing the random samples into a univariate or a multivariate framework. It was found that homogeneity existed within the flood vector series. It was concluded therefore that the time series of the flood vectors do not exhibit any significant trend.

Based on analytically based fitness measures, it was concluded that it is likely that Triweight kernel function is the best-fitted distribution for defining the marginal distribution of peak flows, flood volumes and flood durations in the Kelantan River basin.

Special thanks are extended to the Drainage and Irrigation Department (DID), Malaysia, for supplying streamflow data of the Kelantan river basin.