maxNR {maxLik} | R Documentation |
Unconstrained maximization based on Newton-Raphson method.
maxNR(fn, grad = NULL, hess = NULL, start, print.level = 0, tol = 1e-08, reltol=sqrt(.Machine$double.eps), gradtol = 1e-06, steptol = 1e-10, lambdatol = 1e-06, qrtol = 1e-10, iterlim = 150, constraints=NULL, activePar=rep(TRUE, length(start)), ...)
fn |
function to be maximized. It must return either a single
number or a numeric vector, which is summed (this is useful for BHHH
method). If the parameters are out of range, fn should
return NA . See details for constant parameters. |
grad |
gradient of the function. If NULL , numeric
gradient is used. It must return a gradient vector, or matrix where
columns correspond to individual parameters. The column sums
are treated as gradient components (this is useful for BHHH method).
|
hess |
Hessian matrix of the function. If missing, numeric Hessian, based on gradient, is used. |
start |
initial value for the parameter vector |
print.level |
this argument determines the level of printing which is done during the minimization process. The default value of 0 means that no printing occurs, a value of 1 means that initial and final details are printed and a value of 2 means that full tracing information for every iteration is printed. Higher values will result in even more details. |
tol |
stopping condition. Stop if the absolute difference
between successive iterations less than tol , return
code=2 . |
reltol |
Relative convergence tolerance. The algorithm stops if it is unable to increase the value by a factor of 'reltol * (abs(val) + reltol)' at a step. Defaults to 'sqrt(.Machine$double.eps)', typically about '1e-8'. |
gradtol |
stopping condition. Stop if the norm of the gradient
less than gradtol , return code=1 . |
steptol |
stopping/error condition. If the quadratic
approximation leads to lower function value instead of higher, or
NA , the step length
is halved and a new attempt is made. This procedure is repeated
until step < steptol , thereafter code=3 is returned. |
lambdatol |
control whether the Hessian is treated as negative
definite. If the
largest of the eigenvalues of the Hessian is larger than
-lambdatol , a suitable diagonal matrix is subtracted from the
Hessian (quadratic hill-climbing) in order to enforce nagetive
definiteness. |
qrtol |
QR-decomposition tolerance |
iterlim |
stopping condition. Stop if more than iterlim
iterations, return code=4 . |
constraints |
either NULL for unconstrained optimization
or a list with two components
eqA and eqB for equality-constrained optimization
A %*% theta + B = 0. The constrained
problem is forwarded
to sumt .
|
activePar |
index vector, which parameters are treated as free. Other parameters are treated as fixed constants. |
... |
further arguments to fn , grad and
hess . |
The algorithm can work with constant parameters and related changes of parameter values. Constant parameters are useful if a parameter value is converging toward the boundary of support, or for testing.
One way is to put
activePar
to non-NULL, specifying which parameters we want to
change. However, parameters can also be fixed in runtime by by
signaling with
fn
. This may be useful if an estimation converges toward the
edge of the parameter space possibly causing problems. The value of
fn
may have following attributes:
Hence, constVal
specifies which parameters are treated as
constants. newVal
allows one to overwrite the existing parameter
values, possibly the non-constant values as well. If the attribute
newVal
is present, the new function value need not to exceed
the previous one (maximization is not performed in that step).
list of class "maxim" with following components:
maximum |
fn value at maximum (the last calculated value
if not converged). |
estimate |
estimated parameter value. |
gradient |
last gradient value which was calculated. Should be close to 0 if normal convergence. |
hessian |
Hessian at the maximum (the last calculated value if not converged). |
code |
return code:
|
message |
a short message, describing code . |
last.step |
list describing the last unsuccessful step if
code=3 with following components:
|
activePar |
logical vector, which parameters are not constants. |
iterations |
number of iterations. |
type |
character string, type of maximization. |
constraints |
A list, describing the constrained optimization
(NULL if unconstrained). Includes the following components:
|
Ott Toomet otoomet@ut.ee
Greene, W., 2002, Advanced Econometrics Goldfeld, S. M. and Quandt, R. E., 1972, Nonlinear Methods in Econometrics. Amsterdam: North-Holland
nlm
for Newton-Raphson optimization,
optim
for different gradient-based optimization
methods.
## ML estimation of exponential duration model: t <- rexp(100, 2) loglik <- function(theta) sum(log(theta) - theta*t) ## Note the log-likelihood and gradient are summed over observations gradlik <- function(theta) sum(1/theta - t) hesslik <- function(theta) -100/theta^2 ## Estimate with numeric gradient and Hessian a <- maxNR(loglik, start=1, print.level=2) summary(a) ## You would probably prefer 1/mean(t) instead ;-) ## Estimate with analytic gradient and Hessian a <- maxNR(loglik, gradlik, hesslik, start=1) summary(a) ## ## Next, we give an example with vector argument: Estimate the mean and ## variance of a random normal sample by maximum likelihood ## Note: you might want to use maxLik instead ## loglik <- function(param) { mu <- param[1] sigma <- param[2] ll <- -0.5*N*log(2*pi) - N*log(sigma) - sum(0.5*(x - mu)^2/sigma^2) ll } x <- rnorm(1000, 1, 2) # use mean=1, stdd=2 N <- length(x) res <- maxNR(loglik, start=c(0,1)) # use 'wrong' start values summary(res) ### ### Now an example of constrained optimization ### ## We maximize exp(-x^2 - y^2) where x+y = 1 f <- function(theta) { x <- theta[1] y <- theta[2] exp(-(x^2 + y^2)) ## Note: you may want to use exp(- theta %*% theta) instead ;-) } ## use constraints: x + y = 1 A <- matrix(c(1, 1), 1, 2) B <- -1 res <- maxNR(f, start=c(0,0), constraints=list(eqA=A, eqB=B), print.level=1) print(summary(res))