parTLgld: Estimate the Parameters of the Generalized Lambda Distribution using Trimmed L-moments (t=1)

Description

This function estimates the parameters of the Generalized Lambda distribution given the trimmed L-moments (TL-moments) for $t=1$ of the data in a TL-moment object with a trim level of unity (trim=1). The relations between distribution parameters and TL-moments are seen under lmomTLgld. There are no simple expressions for the parameters in terms of the L-moments. Consider that multiple parameter solutions are possible with the Generalized Lambda distribution so some expertise with this distribution and other aspects is advised.

Usage

parTLgld(lmom, verbose=FALSE, initkh=NULL, eps=1e-3,
         aux=c("tau5", "tau6"), checklmom=TRUE, ...)

Value

An R

list is returned if result='best'.

type: The type of distribution: gld.
para: The parameters of the distribution.
delTau5: Difference between $\tilde{\tau}^{(1)}_5$ of the fitted distribution and true $\hat{\tau}^{(1)}_5$.
error: Smallest sum of square error found.
source: The source of the parameters: “parTLgld”.
rest: An R data.frame of other solutions if found.

The rest of the solutions have the following:

xi: The location parameter of the distribution.
alpha: The scale parameter of the distribution.
kappa: The 1st shape parameter of the distribution.
h: The 2nd shape parameter of the distribution.
attempt: The attempt number that found valid TL-moments and parameters of GLD.
delTau5: The absolute difference between $\hat{\tau}^{(1)}_5$ of data to $\tilde{\tau}^{(1)}_5$ of the fitted distribution.
error: The sum of square error found.
initial_k: The starting point of the $\kappa$ parameter.
initial_h: The starting point of the $h$ parameter.
valid.gld: Logical on validity of the GLD---TRUE by this point.
valid.lmr: Logical on validity of the L-moments---TRUE by this point.
lowerror: Logical on whether error was less than eps---TRUE by this point.

Arguments

lmom: A TL-moment object created by TLmoms.
verbose: A logical switch on the verbosity of output. Default is verbose=FALSE.
initkh: A vector of the initial guess of the $\kappa$ and $h$ parameters. No other regions of parameter space are consulted.
eps: A small term or threshold for which the square root of the sum of square errors in $\tau_3$ and $\tau_4$ is compared to to judge “good enough” for the alogrithm to order solutions based on smallest error as explained in next argument.
aux: Control the algorithm to order solutions based on smallest error in trimmed $\Delta \tau_5$ or $\Delta \tau_6$.
checklmom: Should the lmom be checked for validity using the are.lmom.valid function. Normally this should be left as the default and it is very unlikely that the L-moments will not be viable (particularly in the $\tau_4$ and $\tau_3$ inequality). However, for some circumstances or large simulation exercises then one might want to bypass this check.
...: Other arguments to pass.

Author

W.H. Asquith

Details

Karian and Dudewicz (2000) summarize six regions of the $\kappa$ and $h$ space in which the Generalized Lambda distribution is valid for suitably choosen $\alpha$. Numerical experimentation suggestions that the L-moments are not valid in Regions 1 and 2. However, initial guesses of the parameters within each region are used with numerous separate optim (the R function) efforts to perform a least sum-of-square errors on the following objective function. $$(\hat{\tau}^{(1)}_3 - \tilde{\tau}^{(1)}_3)^2 + (\hat{\tau}^{(1)}_4 - \tilde{\tau}^{(1)}_4)^2 \mbox{, }$$ where $\tilde{\tau}^{(1)}_r$ is the L-moment ratio of the data, $\hat{\tau}^{(1)}_r$ is the estimated value of the TL-moment ratio for the current pairing of $\kappa$ and $h$ and $\tau^{(1)}_r$ is the actual value of the L-moment ratio.

For each optimization a check on the validity of the parameters so produced is made--are the parameters consistent with the Generalized Lambda distribution and a second check is made on the validity of $\tau^{(1)}_3$ and $\tau^{(1)}_4$. If both validity checks return TRUE then the optimization is retained if its sum-of-square error is less than the previous optimum value. It is possible for a given solution to be found outside the starting region of the initial guesses. The surface generated by the $\tau^{(1)}_3$ and $\tau^{(1)}_4$ equations seen in lmomTLgld is complex; different initial guesses within a given region can yield what appear to be radically different $\kappa$ and $h$. Users are encouraged to “play” with alternative solutions (see the verbose argument). A quick double check on the L-moments (not TL-moments) from the solved parameters using lmomTLgld is encouraged as well.

References

Asquith, W.H., 2007, L-moments and TL-moments of the generalized lambda distribution: Computational Statistics and Data Analysis, v. 51, no. 9, pp. 4484--4496.

Karian, Z.A., and Dudewicz, E.J., 2000, Fitting statistical distributions---The generalized lambda distribution and generalized bootstrap methods: CRC Press, Boca Raton, FL, 438 p.

Examples

Run this code

# As of version 1.6.2, it is felt that in spirit of CRAN CPU
# reduction that the intensive operations of parTLgld() should
# be kept a bay.

if (FALSE) {
X <- rgamma(202,2) # simulate a skewed distribution
lmr <- TLmoms(X, trim=1) # compute trimmed L-moments
PARgldTL <- parTLgld(lmr) # fit the GLD

F <- pp(X) # plotting positions for graphing
plot(F,sort(X), col=8, cex=0.25)
lines(F, qlmomco(F,PARgldTL)) # show the best estimate
if(! is.null(PARgldTL$rest)) { 
  n <- length(PARgldTL$rest$xi)
  other <- unlist(PARgldTL$rest[n,1:4]) # show alternative
  lines(F, qlmomco(F,vec2par(other, type="gld")), col=2)
}
# Note in the extraction of other solutions that no testing for whether
# additional solutions were found is made. Also, it is quite possible
# that the other solutions "[n,1:4]" is effectively another numerical
# convergence on the primary solution. Some users of this example thus
# might not see two separate lines. Users are encouraged to inspect the
# rest of the solutions: print(PARgld$rest)

# For one run of the above example, the GLD results follow
#print(PARgldTL)
#$type
#[1] "gld"
#$para
#         xi       alpha       kappa           h
# 1.02333964 -3.86037875 -0.06696388 -0.22100601
#$delTau5
#[1] -0.02299319
#$error
#[1] 7.048409e-08
#$source
#[1] "pargld"
#$rest
#         xi     alpha       kappa          h attempt     delTau5        error
#1  1.020725 -3.897500 -0.06606563 -0.2195527       6 -0.02302222 1.333402e-08
#2  1.021203 -3.895334 -0.06616654 -0.2196020       4 -0.02304333 8.663930e-11
#3  1.020684 -3.904782 -0.06596204 -0.2192197       5 -0.02306065 3.908918e-09
#4  1.019795 -3.917609 -0.06565792 -0.2187232       2 -0.02307092 2.968498e-08
#5  1.023654 -3.883944 -0.06668986 -0.2198679       7 -0.02315035 2.991811e-07
#6 -4.707935 -5.044057  5.89280906 -0.3261837      13  0.04168800 2.229672e-10
}

if (FALSE) {
F <- seq(.01,.99,.01)
plot(F,qlmomco(F, vec2par(c( 1.02333964, -3.86037875,
                            -0.06696388, -0.22100601), type="gld")),
     type="l")
lines(F,qlmomco(F, vec2par(c(-4.707935, -5.044057,
                              5.89280906, -0.3261837), type="gld")))
}

Run the code above in your browser using DataLab