R: Optimization wrapper for likelihood-based procedures

mleOptimWrapper {sdetorus}

R Documentation

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"))

[Package sdetorus version 0.1.10 Index]