maxlik-fd: Maximization of the Spectral Likelihood Function
In stsm: Structural Time Series Models

Description Usage Arguments Details Value References See Also Examples

Maximize the spectral log-likelihood function of a structural time series model by means of a scoring algorithm or a general purpose optimization algorithm available in optim.

maxlik.fd.optim(m, 
  barrier = list(type = c("1", "2"), mu = 0), inf = 99999, 
  method = c("BFGS", "L-BFGS-B", "Nelder-Mead", "CG", "SANN"), 
  gr = c("analytical", "numerical"), optim.control = list())

maxlik.fd.scoring(m, step = NULL, 
  information = c("expected", "observed", "mix"),
  ls = list(type = "optimize", tol = .Machine$double.eps^0.25, cap = 1),
  barrier = list(type = c("1", "2"), mu = 0), 
  control = list(maxit = 100, tol = 0.001, trace = FALSE, silent = FALSE), 
  debug = FALSE)
  
maxclik.fd.scoring(m, step = NULL, information = c("expected", "observed"),
  ls = list(type = "optimize", tol = .Machine$double.eps^0.25, cap = 1),
  barrier = list(type = c("1", "2"), mu = 0), 
  control = list(maxit = 100, tol = 0.001, trace = FALSE, silent = FALSE))

`m`	an object of class `stsm`.

barrier

a list defining a barrier term to penalize parameter values close to the bounds m@lower and m@upper.

inf

a numeric indicating the value to be returned if the value of the log-likelihood function happens to be NA or non-finite at some iteration of the optimization procedure.

method

character indicating the method to be used by optim.

`gr`	character indicating whether numerical or analytical derivatives should be used.

optim.control

a list of control parameters passed to optim.

step

if it is a numeric it stands for a fixed step size, otherwise an automatic procedure is used to choose the step size.

information

the type of information about second order derivatives used to project the gradient in the scoring algorithm.

`ls`	control parameters for the line search procedure used to chose the step size in the scoring algorithm.

control

a list of control parameters for the scoring algorithm.

debug

logical. If TRUE, tracing information is printed for debugging purposes of the scoring algorithm.

The matrix used to project the gradient may be based on expected or observed information. The former, information = "expected", uses the analytical expression of the information matrix. The latter, information = "observed", uses the analytical expression of the Hessian as in a Newton-Raphson procedure. The option information = "mix" uses a mixture of both expressions; in simulations it performed similar to the information matrix, this option may be removed in future versions of the package.

maxclik.fd.scoring maximizes the concentrated likelihood function. The parameter to be concentrated must be defined in the slot cpar of the input model m, see stsm. maxlik.fd.optim detects whether cpar is defined in the input model. In the scoring algorithm, if m@cpar is not NULL maxclik.fd.scoring should be used.

Bounds on parameters and barrier term. The lower and upper bounds within which the L-BFGS-B algorithm conducts the search are taken from the slots lower and upper defined in the input object m. As an alternative to the L-BFGS-B procedure, a barrier term can be passed as argument. The barrier term is added to the likelihood function and acts as a penalization for parameter values close to the bounds. For details about the barrier term see barrier.eval. In the scoring algorithm, if step = NULL the procedure automatically searches the optimum step size that is compatible with the bounds on the parameters. See step.maxsize and further details below.

Control parameters for the scoring algorithm. maxit, maximum number of iterations (the default is 100); tol, tolerance to assess convergence of the algorithm (the default is 0.001); trace, logical, if TRUE, the values of the parameters at each iteration of the procedure are stored and returned in a matrix; silent, logical, if FALSE, a warning is printed if convergence is not achieved.

Choice of the step size in the scoring algorithm. If step is a numerical value, the step size is fixed to that value at all iterations of the algorithm. Otherwise, the choice of the step size in the scoring algorithm is made by means of a line search procedure specified in the argument ls. See step.maxsize for a description of the elements that can be passed through the argument ls.

External regressors If external regressors are included in the model m, starting values for their coefficients are obtained in a linear regression of the differenced series on the differenced regressors. The values in the slot pars are therefore overwritten and not used as initial values.

A list of class stsmFit with components:

`call`	an object of class `call` specifying the arguments passed to the function.
`init`	initial parameter values.
`pars`	parameter values at the local optimum.
`m`	the `stsm` model object updated with the optimal parameter values.
`loglik`	the value of the log-likelihood function at the local optimum.
`convergence`	convergence code returned by `optim`; for the scoring algorithm a logical indicating whether convergence was achieved.
`iter`	for `maxlik.fd.optim` it is a two-element vector with the number of calls made by `optim` to the objective function and to the gradient; for `maxlik.fd.scoring` it is the number of iterations employed by the scoring algorithm.
`message`	an empty character or a character message giving some additional information about the optimization process.
`Mpars`	a matrix or `NULL`. If `control$trace = TRUE` in the scoring algorithm, the path to the local optimum is traced storing by rows the parameter values at each iteration.
`steps`	a vector or `NULL`. If `control$trace = TRUE` in the scoring algorithm, the step size used at each iteration is stored in this vector.
`ls.iter`	a vector containing the number of iterations employed by the line search procedure at each step of the scoring algorithm. It is `NULL` if `ls$type = "optimize"`.
`ls.counts`	a two-element vector containing the total number of calls to the objective function and the gradient made by the line search procedure in all the iterations. (It is `NULL` if `ls$type = "optimize"`.)

Note: if m@transPars is not NULL, the elements init and pars are in terms of the auxiliary set of parameters. If the output is stored for example in an object called res, get.pars(res$model) will return the actual variance parameters.

The version based on optim, maxlik.fd.optim, returns also the element hessian containing the numerically differentiated Hessian matrix at the local optimum. Note that, if the model is parameterized in terms of an auxiliary set of parameters and gradient = "numerical" is used, the Hessian returned by optim is defined with respect to the auxiliary set of parameters, not the variances.

Harvey, A. C. (1989). Forecasting, Structural Time Series Models and the Kalman Filter. Cambridge University Press.

Nocedal, J. and Wright, J. W. (2006). Numerical Optimization. Springer-Verlag.

barrier.eval, mloglik.fd, stsm, optim.

# fit the local level plus seasonal model to a 
# sample simulated series
# further examples and options can be explored in the 
# script files 'sim-llmseas-ml-fd.R' and 'sim-llmseas-mcl-fd.R'
# available in the 'inst' folder of the source package
data("llmseas")

# initial parameters and 'stsm' model
initpars <- c(var1 = 1, var2 = 1, var3 = 1)
m <- stsm.model(model = "llm+seas", y = llmseas, pars = initpars)

# Newton-Raphson algorithm (analytical Hessian)
res1 <- maxlik.fd.scoring(m = m, step = NULL, 
  information = "observed", control = list(maxit = 100, tol = 0.001))
res1

# Scoring algorithm (information matrix)
res2 <- maxlik.fd.scoring(m = m, step = NULL, 
  information = "expected", control = list(maxit = 100, tol = 0.001))
res2

# wrapper function for 'optim()' in the 'stats' package
res3 <- maxlik.fd.optim(m, method = "L-BFGS-B", gr = "analytical")
res3

# concentrating one of the parameters
# the model must be first defined accordingly, here 'var1', i.e.,
# the variance of the disturbance in the observation equation
# is concentrated, its standard error is reported as 'NA'
mc <- stsm.model(model = "llm+seas", y = llmseas, 
  pars = initpars[-1], cpar = initpars[1])
res4 <- maxclik.fd.scoring(m = mc, step = NULL, 
    information = "observed", control = list(maxit = 100, tol = 0.001))
res4