Population and sample L-moments

The theory and application of the L-moments, as an alternative system to describe the forms of the FDPs, were presented by ^{Hosking (1990)}. L-moments are linear combinations of probability weighted moments (MPP), developed by ^{Greenwood et al. (1979)}, and are defined, for a random variable x with accumulated FDP F (x) and function of quantiles x (F) as:

The expression of the moment: β_r = M_1,r,0 is:

The ordinary definition of moments is:

When compared with the previous equation, it follows that the conventional definition of moments involves successive powers of the function of quantiles x(F) and that the MPPs imply successive powers of F, and therefore can be considered as integrals of x(F) weighted by polynomials F^r , hence its name. MPPs significantly improve the sampling properties because extreme or scattered values do not influence them (^{Asquith, 2011}). The first three population L-moments correspond to (^{Hosking, 1990}; ^{Hosking and Wallis, 1997}):

In a sample of size n, with its elements in ascending order (x₁ ≤ x₂ ≤ …≤ x_n) the unbiased estimators of β_r are (^{Stedinger et al., 1993}; ^{Hosking and Wallis, 1997}):

The sampling estimators of λ _r will be l _r and the equations 4 to 6 define them.

Fit with L-moments of the distribution GVE₀

The GVE distribution is now the basic probabilistic model of the AHFs of maximum and minimum data, such as flows, rainfalls, levels, temperatures and winds. ^{Hosking and Wallis (1997)} and ^{El Adlouni et al. (2008)} highlighted that acceptance of the GVE function is due to the fact that when its form parameter is negative (k < 0) has its right tail thicker or denser than the other FDPs, which are usually used in AHFs; in addition, it has an upper limit when k > 0 and defines the Gumbel distribution of two parameters, when k = 0. Therefore, its application interval in x is: u + a/k ≤ x < ∞ if k < 0; - ∞ < x < ∞ if k = 0 and - ∞ < x ≤ u + a/k if k > 0.

The inverse solution x(F) of the GVE₀ distribution is as follows, with it, predictions associated with a certain probability of non-exceedance (F) are estimated and the equations that allow us to estimate their three fit parameters (u, a, k) corresponding to the location, scale and shape, with the L- moments method are:

where:

To estimate the Gamma function the Stirling formula was used (^{Davis, 1972}):

Verification of the record non-stationarity

A chronological series of maximum hydrological data is stationary if it is free of trends, jumps or periodicity. This implies that the statistical parameters of that series, such as the mean and variance, are basically constant over time; otherwise, the series will not be stationary (^{Cunderlik and Burn, 2003}; Khaliq et al., 2003). ^{Mudersbach and Jensen (2010 a)} indicated that the study of higher-order moments, such as asymmetry, is unnecessary because it is unlikely that this property will vary over time, but in a non-stationary record, the location and scale parameters of an FDP will surely change over time.

^{Campos-Aranda (2015)} estimated the linear trend with several methods and verified if it is statistically different from zero. The change in variance was tested with a simplified version of the procedure suggested by ^{Al Saji et al. (2015)}. This change consists of subtracting from each data the estimated linear trend and square them, to observe if they are relatively constant or have periods with high and low variability. This allowed us to conclude that the record is not stationary in the variance.

Fitting with mobile L-moments of the GVE₁₁distribution

The non-stationary GVE distribution, with its parameters of location (u) and scale (a) variables in time (t), is indicated in large records of maximum annual data with non-constant trend and variability (^{Cunderlik and Burn, 2003}). The linear or exponential variation of u and a is recommended because of its simple fitting and the ease of extrapolation, that is (^{Mudersbach and Jensen, 2010 a}, ^b):

In the mobile L-moments method a sliding time window is used, with length in years n that slides year after year through the record of maximum data of size n_t. The width of the window (n) will be selected so that in that amplitude of the record it can be considered stationary and then the number of windows (nv) to be processed will be:

For each section of the record of the nv mobile windows its location and scale parameters are calculated based on equations 11 to 15. The two chronological series obtained will be nonparametric functions that cannot be extrapolated to the future and that will be characterized by their average values: u_m , a_m and k_m . By adjusting to the chronological series of u the models defined by equations 16 and 18 and the series of a, those of equations 17 and 19, parametric descriptions of the mobile parameters are defined, whose quality of fit is measured with the mean standard error, expressed as:

In the previous expressions, the estimated values [u^t,a^t ] are calculated with equations 16 to 19, previous evaluation of their fit coefficients obtained by minimizing the squares of the residuals. For the linear and exponential variation of the location parameter [u (t)] the values are calculated with the following expressions:

In the above expressions "ln" is the natural logarithm and "exp" the exponential function, that is, the number e raised to the value of the parentheses. Changing in the equations 24 to 27, u_i by a_i and u_m by a_m , the coefficients β and δ of equations 17 and 19 are estimated.

For the predictions of a certain design return period (Tr) in years, estimated in a certain future time, commonly the years 2025, 2050 and 2100, equation 10 is applied with k_m and the following values of the probability of non-exceedance (F) and of the final time (T_f), in the location and scale parameters (equations 16 to 19), so there will be predictions with linear and exponential extrapolation:

where AFP and AFR are the final years of the prediction (2025, 2050 or 2100) and of the record.

Standard error of fit

In the mid-1970s the standard error of fit (EEA) was established as a quantitative statistical indicator of the quality of the fit achieved between data and predictions with the FDP that is tested, since it evaluates the standard deviation of those differences. In this study the FDP that is tested is the non-stationary distribution GVE₁₁ and its expression is as follows (^{Kite, 1977}, ^{Pandey and Nguyen, 1999}):

where n _t and np are the number of record data and of parameters of fit, in the stationary case three (u, a, k ) and in the non-stationary case five (α₀, α₁, β₀, β₁, k _m , or, γ₀, γ₁, δ₀, δ₁, k _m); X_i are the maximum data ordered from lower to higher and X^i are the estimated values with x(F) or quantile function (equation 10) that uses the variable location and scale parameters (linearly or exponentially), for an estimated probability of non-exceedance with the Weibull formula (^{Benson, 1962}):

where m is the order number of the data, with 1 for the lowest and n_t for the highest. The evaluation of the EEA allows us to select between the linear and exponential models to extrapolate the design predictions.

Fitting of non-stationary distributions LOG₁₁ and PAG₁₁

FDP Generalized Logistics (LOG) and Generalized Pareto (PAG) are models used in the AHFs (^{Stedinger et al., 1993}; ^{Kim et al., 2015}), which can be fitted in their non-stationary versions, with the mobile L-moments method, by changing equations 11 to 14, to estimate their location and scale parameters and then make their predictions with equation 10, which corresponds to them. The formulas for the stationary distributions LOG₀ and PAG₀ were documented by ^{Campos-Aranda (2016)}. The mobile L-moments method can be applied using the higher-order L-moments (^{Wang, 1997}, ^{Campos-Aranda, 2016}) and obtaining the best fit in the recording sections, which define the sliding windows.

Records of annual maximum values to be processed

Two records of specialized literature were selected to illustrate the application of the non-stationary distribution GVE₁₁ fitted with the mobile L-moments method (^{Mudersbach y Jensen, 2010 a}, ^b). the data are approximate to those shown in the original graphs (Table 1).

The first record, taken from ^{Vogel et al. (2011)}, included 69 maximum annual expenditure data (1940 to 2008), at a hydrometric station of the Aberjona River, near the city of Boston, Massachusetts, USA, and whose drainage basin was impacted by urban development. This series had a behavior with upward linear trend and increasing variability at the end of the record (Figura 1).

The second record, taken from ^{Katz (2013)}, is an example of the drastic reduction of the maximum daily precipitation (PMD) in winter (May to October) recorded by the Manjimup station in the southwestern extreme of Australia. The record included 75 annual values (1930 to 2004), with a downward trend and greater variability at the beginning of the record (Figura 2).

Q Predictions in a station of the Aberjona River, USA

The data (Table 1) had an upward linear and significant trend, with statistics of the Student distribution of 3.850 and a critical value of 1.996 (p ≤ 0.05). The diagram of squared differences, between each data and the calculated linear trend showed that the variability was greater towards the end of the record (Figura 3).

The linear model with the time t of the fit parameters u and a, leads to lower standard errors of fit (EEA) compared to those of the exponential model (Cuadro 2). With a mobile window width of 20 data, the lowest EEA was obtained, so these results were adopted to obtain predictions for the future (Figura 4).

With the equations 11 to 15, the fit parameters of the stationary model (GVE₀) were obtained: u = 8.407, a = 4.481 and k = -0.3034, with an EEA = 1.556 m³·s^-1, and with equation 10 the predictions (Table 3). The non-stationary distribution GVE₁₁ provided higher predictions in all return periods and even at the end of the historical record (Table 3). In the extreme periods of return (500 and 1000 years) and towards the year 2100, the predictions of the non-stationary model were twice those of the stationary model that did not consider the upward trend or the greater variability towards the end of the record.

PMD predictions in the Manjimup station of Australia

The annual daily maximum precipitation data (PMD) (Table 1) have a significant linear downward trend (Figure 2), with statistics of the Student distribution of -4.085 and critical value of 1.993 with p ≤ 0.05. In the diagram of squared differences, between each data and linear trend, the variability was noticeably higher at the beginning of the record (Figure 5).

The exponential model with the time t of the fit parameters u and a, leads to lower standard errors of fit (EEA) compared with those of the linear model (Cuadro 4). As the size of the mobile window grows, lower EEAs are obtained, whereby those of the maximum proven value of n = 30 are adopted (Figura 6).

Based on equations 11 to 15, the fit parameters of the stationary model GVE₀ were obtained: u = 35.488, a = 8.431 and k = -0.1846, with an EEA = 2.881 mm. The respective predictions, obtained with equation 10, showed that the non-stationary distribution GVE₁₁ contributed lower predictions, from the end of the historical record and in all the return periods analyzed (Cuadro 5). In the extreme return periods analyzed (500 and 1000 years) and towards the year 2100, predictions of PMD tended to stabilize around 21 mm.

Applicability of the mobile L-moments method

The linear representation of both parameters and significantly better of the location parameter (u) is acceptable for the Aberjona River, USA (Figure 4). At the Manjimup station, Australia (Figure 6) the scale parameter representation (a) may be acceptable, but the one relative to the location parameter u may be not acceptable.

The higher statistical quality of the fit with the mobile L-moments method is numerically highlighted because of the similarity shown by the standard errors of fit (EEA) of the stationary (GVE₀) and non-stationary (GVE₁₁) models. In the case of the Aberjona River, USA, the EEA of the stationary and non-stationary GVE were 1.556 and 2.136 m³·s^-1, that is, with an order of similar magnitude. In contrast, at the Manjimup station, Australia, the EEA of both models (2.881 and 11.199 mm) were different and much greater that of the GVE₁₁ distribution, therefore, is not acceptable to the PMD data group.