powell (version 1.0-0)

powell: powell

Description

Optimizes a function using Powell's UObyQA algorithm.

Usage

powell(par, fn, control = powell.control(), check.hessian = TRUE, ...)

Arguments

par

Starting values for objective function

fn

A function to be optimized. The function takes the parameters (par) as its first argument.

control

A list of control parameters

check.hessian

logical; if TRUE an eigenvalue decomposition is used to check the hessian for positive definiteness.

...

Additional arguments to be passed to fn

Value

A list with components

par

The final values of the parameters.

value

The final value of the function being optimized.

counts

The number of times the function is called.

hessian

A symmetric matrix of the estimated Hessian.

eigen.hessian

If check.hessian is TRUE the eigenvalues and eigenvectors; otherwise NULL.

convergence

0 if converged, 1 otherwise.

control

The input control parameters.

message

Information about the model fit. This will be non-null only if check.hessian is TRUE and the resulting hessian is not positive definite.

call

The original call to the optimizer.

Details

This function seeks the least value of a function of many variables, by a trust region method that forms quadratic models by interpolation. The algorithm is described in "UOBYQA: unconstrained optimization by quadratic approximation" by M.J.D. Powell, Report DAMTP 2000/NA14, University of Cambridge.

References

http://plato.asu.edu/topics/problems/nlounres.html

See Also

optim

Examples

set.seed(1)
fn <- function(beta, y, x, w) {
  # binomial deviance using double log link
  mu <- exp(x %*% beta)
  logLik <-  - y * mu + (w - y) * log(1 - exp(-mu))
  -2 * sum(logLik)
}
n <- 1000
beta <- c(-1, 0.5)
w <- rpois(n, 100)
x <- rep(c("A", "B"), length = n)
X <- model.matrix(~ x, data.frame(x))
y <- rbinom(n, w, exp(-exp(X %*% beta)))

x1 <- powell(beta, fn, y = y, x = X, w = w)
x2 <- optim(beta, fn, y = y, x = X, w = w, hessian = TRUE)
x3 <- glm(1 - y/w ~ x, data = data.frame(x, y, w),
          family = binomial("cloglog"), weights = w)

# compare coefficients from 3 fits
rbind(powell = x1$par, optim = x2$par, glm = coef(x3))
# compare standard errors from 3 fits
rbind(powell = sqrt(diag(2 * solve(x1$hessian))), 
      optim = sqrt(diag(2 * solve(x2$hessian))), 
      glm = sqrt(diag(vcov(x3))))