# mle2par: Use Maximum Likelihood to Estimate Parameters of a... In lmomco: L-Moments, Censored L-Moments, Trimmed L-Moments, L-Comoments, and Many Distributions

## Description

This function uses the method of maximum likelihood (MLE) to estimate the parameters of a distribution.

MLE is a straightforward optimization problem that is formed by maximizing the sum of the logarithms of probability densities. Let Θ represent a vector of parameters for a candidate fit to the specified probability density function g(x|Θ) and x_i represent the observed data for a sample of size n. The objective function is

\mathcal{L}(Θ) = -∑_{i=1}^{n} \log\, g(x_i|Θ)\mbox{,}

where the Θ for a maximized {-}\mathcal{L} (note the 2nd negation for the adjective “maximized”, optim() defaults as a minimum optimizer) represents the parameters fit by MLE. The initial parameter estimate by default will be seeded by the method of L-moments.

## Usage

 1 2 3 mle2par(x, type, para.int=NULL, silent=TRUE, null.on.not.converge=TRUE, ptransf= function(t) return(t), pretransf=function(t) return(t), ...) 

## Arguments

 x A vector of data values. type Three character (minimum) distribution type (for example, type="gev"). para.int Initial parameters as a vector Θ or as an lmomco parameter “object” from say vec2par. If a vector is given, then internally vec2par is called with distribution equal to type. silent A logical to silence the try() function wrapping the optim() function. null.on.not.converge A logical to trigging simple return of NULL if the optim() function returns a nonzero convergence status. ptransf An optional parameter transformation function (see Examples) that is useful to guide the optimization run. For example, suppose the first parameter of a three parameter distribution resides in the positive domain, then ptransf(t) =  function(t) c(log(t[1]), t[2], t[3]). pretransf An optional parameter retransformation function (see Examples) that is useful to guide the optimization run. For example, suppose the first parameter of a three parameter distribution resides in the positive domain, then pretransf(t) =  function(t) c(exp(t[1]), t[2], t[3]). ... Additional arguments for the optim() function and other uses.

## Value

An R list is returned. This list should contain at least the following items, but some distributions such as the revgum have extra.

 type The type of distribution in three character (minimum) format. para The parameters of the distribution. source Attribute specifying source of the parameters. AIC The Akaike information criterion (AIC). optim The returned list of the optim() function.

## Note

During the optimization process, the function requires evaluation at the initial parameters. The following error rarely will be seen:

 1 2  Error in optim(para.int$para, afunc) : function cannot be evaluated at initial parameters  if Inf is returned on first call to the objective function. The silent by default though will silence this error. Alternative starting parameters might help. This function is not built around subordinate control functions to say keep the parameters within distribution-specific bounds. However, in practice, the L-moment estimates should already be fairly close and the optimizer can take it from there. More sophisticated MLE for many distributions is widely available in other R packages. The lmomco package uses its own probability density functions. ## Author(s) W.H. Asquith ## See Also lmom2par, mps2par ## Examples   1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 ## Not run: # This example might fail on mle2par() or mps2par() depending on the values # that stem from the simulation. Trapping for a NULL return is not made here. father <- vec2par(c(37,25,114), type="st3"); FF <- nonexceeds(); qFF <- qnorm(FF) X <- rlmomco(78, father) # rerun if MLE and MPS fail to get a solution plot(qFF, qlmomco(FF, father), type="l", xlim=c(-3,3), xlab="STANDARD NORMAL VARIATE", ylab="QUANTILE") # parent (black) lines(qFF, qlmomco(FF, lmom2par(lmoms(X), type="gev")), col=2) # L-moments (red) lines(qFF, qlmomco(FF, mps2par(X, type="gev")), col=3) # MPS (green) lines(qFF, qlmomco(FF, mle2par(X, type="gev")), col=4) # MLE (blue) points(qnorm(pp(X)), sort(X)) # the simulated data ## End(Not run) ## Not run: # REFLECTION SYMMETRY set.seed(451) X <- rlmomco(78, vec2par(c(2.12, 0.5, 0.6), type="pe3")) # MLE and MPS are almost reflection symmetric, but L-moments always are. mle2par( X, type="pe3")$para # 2.1796827 0.4858027 0.7062808 mle2par(-X, type="pe3")$para # -2.1796656 0.4857890 -0.7063917 mps2par( X, type="pe3")$para # 2.1867551 0.5135882 0.6975195 mps2par(-X, type="pe3")$para # -2.1868252 0.5137325 -0.6978034 parpe3(lmoms( X))$para # 2.1796630 0.4845216 0.7928016 parpe3(lmoms(-X))$para # -2.1796630 0.4845216 -0.7928016 ## End(Not run) ## Not run: Ks <- seq(-1,+1,by=0.02); n <- 100; MLE <- MPS <- rep(NA, length(Ks)) for(i in 1:length(Ks)) { sdat <- rlmomco(n, vec2par(c(1,0.2,Ks[i]), type="pe3")) mle <- mle2par(sdat, type="pe3")$para[3] mps <- mps2par(sdat, type="pe3")$para[3] MLE[i] <- ifelse(is.null(mle), NA, mle) # A couple of failures expected as NA's. MPS[i] <- ifelse(is.null(mps), NA, mps) # Some amount fewer failures than MLE. } plot( MLE, MPS, xlab="SKEWNESS BY MLE", ylab="SKEWNESS BY MPS")# ## End(Not run) ## Not run: # Demonstration of parameter transformation and retransformation set.seed(9209) # same seed used under mps2par() in parallel example x <- rlmomco(500, vec2par(c(1,1,3), type="gam")) # 3-p Generalized Gamma guess <- lmr2par(x, type="gam", p=3) # By providing a 3-p guess the 3-p # Generalized Gamma will be triggered internally. There are problems passing # "p" argument to optim() if that function is to pick up the ... argument. mle2par(x, type="gam", para.int=guess, silent=FALSE, ptransf= function(t) { c(log(t[1]), log(t[2]), t[3])}, pretransf=function(t) { c(exp(t[1]), exp(t[2]), t[3])})$para # Reports: mu sigma nu for some simulated data. # 1.0341269 0.9731455 3.2727218 ## End(Not run) ## Not run: # Demonstration of parameter estimation with tails of density zero, which # are intercepted internally to maintain finiteness. We explore the height # distribution for male cats of the cats dataset from the MASS package and # fit the generalized lambda. The log-likelihood is shown by silent=FALSE # to see that the algorithm converges slowly. and shows how we can control # the relative tolerance of the optim() function as shown below and # investigate the convergence by reviewing the five fits to the data. FF <- nonexceeds(sig6=TRUE) library(MASS); data(cats); x <- cats$Hwt[cats$Sex == "M"] p2 <- mle2par(x, type="gld", silent=FALSE, control=list(reltol=1E-2)) p3 <- mle2par(x, type="gld", silent=FALSE, control=list(reltol=1E-3)) p4 <- mle2par(x, type="gld", silent=FALSE, control=list(reltol=1E-4)) p5 <- mle2par(x, type="gld", silent=FALSE, control=list(reltol=1E-5)) p6 <- mle2par(x, type="gld", silent=FALSE, control=list(reltol=1E-6)) plot( FF, quagld(FF, p2), type="l", col=1); points(pp(x), sort(x)) lines(FF, quagld(FF, p3), col=2); lines(FF, par2qua(FF, p4), col=3) lines(FF, quagld(FF, p5), col=4); lines(FF, par2qua(FF, p6), col=6) # ## End(Not run) 

lmomco documentation built on Sept. 21, 2018, 6:40 p.m.