mleOptimWrapper: Optimization wrapper for likelihood-based procedures

In egarpor/sdetorus: Statistical Tools for Toroidal Diffusions

Optimization wrapper for likelihood-based procedures

Description

A convenient wrapper to perform local optimization of the likelihood function via nlm and optim including several practical utilities.

Usage

mleOptimWrapper(minusLogLik, region = function(pars) list(pars = pars,
  penalty = 0), penalty = 1e+10, optMethod = "Nelder-Mead", start,
  lower = rep(-Inf, ncol(start)), upper = rep(Inf, ncol(start)),
  selectSolution = "lowestLocMin", checkCircular = TRUE, maxit = 500,
  tol = 1e-05, verbose = 0, eigTol = 1e-04, condTol = 10000, ...)

Arguments

minusLogLik	function computing the minus log-likelihood function. Must have a single argument containing a vector of length p.
region	function to impose a feasibility region via a penalty. See details.
penalty	imposed penalty if value is not finite.
optMethod	one of the following strings: "nlm", "Nelder-Mead", "BFGS", "CG", "L-BFGS-B", "SANN", or "Brent".
start	starting values, a matrix with p columns, with each entry representing a different starting value.
lower, upper	bound for box constraints as in method "L-BFGS-B" of optim.
selectSolution	which criterion is used for selecting a solution among possible ones, either "lowest", "lowestConv" or "lowestLocMin". "lowest" returns the solution with lowest value in the minusLogLik function. "lowestConv" restricts the search of the best solution among the ones for which the optimizer has converged. "lowestLocMin" in addition imposes that the solution is guaranteed to be a local minimum by examining the Hessian matrix.
checkCircular	logical indicating whether to automatically treat the variables with lower and upper entries equal to -pi and pi as circular. See details.
maxit	maximum number of iterations.
tol	tolerance for convergence (passed to reltol, pgtol or gradtol).
verbose	an integer from 0 to 2 if optMethod = "Nelder-Mead" or from 0 to 4 otherwise giving the amount of information displayed.
eigTol, condTol	eigenvalue and condition number tolerance for the Hessian in order to guarantee a local minimum. Used only if selectSolution = "lowestLocMin".
...	further arguments passed to the optMethod selected. See options in nlm or optim.

Details

If checkCircular = TRUE, then the corresponding lower and upper entries of the circular parameters are set to -Inf and Inf, respectively, and minusLogLik is called with the principal value of the circular argument.

If no solution is found satisfying the criterion in selectSolution, NAs are returned in the elements of the main solution.

The Hessian is only computed if selectSolution = "lowestLocMin".

Region feasibility can be imposed by a function with the same arguments as minusLogLik that resets pars in to the boundary of the feasibility region and adds a penalty proportional to the violation of the feasibility region. Note that this is not the best procedure at all to solve the constrained optimization problem, but just a relatively flexible and quick approach (for a more advanced treatment of restrictions, see optimization-focused packages). The value must be a list with objects pars and penalty. By default no region is imposed, i.e., region = function(pars) list("pars" = pars, "penalty" = 0). Note that the Hessian is computed from the unconstrained problem, hence localMinimumGuaranteed might be FALSE even if a local minimum to the constrained problem was found.

Value

A list containing the following elements:

• par: estimated minimizing parameters

• value: value of minusLogLik at the minimum.

• convergence: if the optimizer has converged or not.

• message: a character string giving any additional information returned by the optimizer.

• eigHessian: eigenvalues of the Hessian at the minimum. Recall that for the same solution slightly different outputs may be obtained according to the different computations of the Hessian of nlm and optim.

• localMinimumGuaranteed: tests if the Hessian is positive definite (all eigenvalues larger than the tolerance eigTol and condition number smaller than condTol).

• solutionsOutput: a list containing the complete output of the selected method for the different starting values. It includes the extra objects convergence and localMinimumGuaranteed.

Examples

# No local minimum
head(mleOptimWrapper(minusLogLik = function(x) -sum((x - 1:4)^2),
                     start = rbind(10:13, 1:2), selectSolution = "lowest"))
head(mleOptimWrapper(minusLogLik = function(x) -sum((x - 1:4)^2),
                     start = rbind(10:13, 1:2),
                     selectSolution = "lowestConv"))
head(mleOptimWrapper(minusLogLik = function(x) -sum((x - 1:4)^2),
                     start = rbind(10:13, 1:2),
                     selectSolution = "lowestLocMin"))

# Local minimum
head(mleOptimWrapper(minusLogLik = function(x) sum((x - 1:4)^2),
                     start = rbind(10:13), optMethod = "BFGS"))
head(mleOptimWrapper(minusLogLik = function(x) sum((x - 1:4)^2),
                     start = rbind(10:13), optMethod = "Nelder-Mead"))

# Function with several local minimum and a 'spurious' one
f <- function(x)  0.75 * (x[1] - 1)^2 -
                  10 / (0.1 + 0.1 * ((x[1] - 15)^2 + (x[2] - 2)^2)) -
                  9.5 / (0.1 + 0.1 * ((x[1] - 15)^2 + (x[2] + 2)^2))
plotSurface2D(x = seq(0, 20, l = 100), y = seq(-3, 3, l = 100), f = f)
head(mleOptimWrapper(minusLogLik = f,
                     start = rbind(c(15, 2), c(15, -2), c(5, 0)),
                     selectSolution = "lowest"))
head(mleOptimWrapper(minusLogLik = f,
                     start = rbind(c(15, 2), c(15, -2), c(5, 0)),
                     selectSolution = "lowestConv"))
head(mleOptimWrapper(minusLogLik = f,
                     start = rbind(c(15, 2), c(15, -2), c(5, 0)),
                     selectSolution = "lowestLocMin", eigTol = 0.01))

# With constraint region
head(mleOptimWrapper(minusLogLik = function(x) sum((x - 1:2)^2),
                     start = rbind(10:11),
                     region = function(pars) {
                       x <- pars[1]
                       y <- pars[2]
                       if (y <= x^2) {
                         return(list("pars" = pars, "penalty" = 0))
                       } else {
                        return(list("pars" = c(sqrt(y), y),
                                    "penalty" = y - x^2))
                       }
                     }, lower = c(0.5, 1), upper = c(Inf, Inf),
                optMethod = "Nelder-Mead", selectSolution = "lowest"))
head(mleOptimWrapper(minusLogLik = function(x) sum((x - 1:2)^2),
                     start = rbind(10:11), lower = c(0.5, 1),
                     upper = c(Inf, Inf),optMethod = "Nelder-Mead"))

egarpor/sdetorus documentation built on March 4, 2024, 1:23 a.m.