Spectral Log-Likelihood Function and Derivatives

Share:

Description

These functions evaluate the negative of the spectral log-likelihood function of a linear Gaussian state space model and its first and second order derivatives.

Usage

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
mloglik.fd(x, model, 
  barrier = list(type = c("1", "2"), mu = 0), inf = 99999, xreg)

mcloglik.fd(x, model, xreg = NULL,
  barrier = list(type = c("1", "2"), mu = 0), inf = 99999)

mloglik.fd.deriv(model, xreg = NULL,
  gradient = TRUE, hessian = TRUE, infomat = TRUE, modcovgrad = TRUE,
  barrier = list(type = c("1", "2"), mu = 0),
  version = c("2", "1"))

mcloglik.fd.deriv(model, xreg = NULL, 
  gradient = TRUE, hessian = TRUE, infomat = TRUE)

mloglik.fd.grad(x, model, xreg = NULL,
  barrier = list(type = c("1", "2"), mu = 0),
  inf)

mcloglik.fd.grad(x, model, xreg = NULL, inf, barrier)

Arguments

x

a numeric vector containing the parameters of the model. This is an auxiliary argument so that this function can be used as input to optim.

model

an object of class stsm.

xreg

optional list containing constant terms. See details.

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.

gradient

logical. If TRUE, first order derivatives of the negative of the spectral log-likelihood function are returned.

hessian

logical. If TRUE, second order derivatives of the negative of the spectral log-likelihood function are returned.

infomat

logical. If TRUE, the information matrix of the spectral log-likelihood are returned.

modcovgrad

logical. If TRUE, a mixture of the analytical expressions for the Hessian and the outer product of the gradient is used. This option is experimental and may be removed in future versions of the package.

version

a character indicating whether implementation "2" or "2" (the default) should be used.They yield the same result but are kept for debugging and comparison of timings. This argument may be removed in future versions.

Details

The spectral log-likelihood of a linear Gaussian state space model is given by (Harvery, 1989 Section 4.3):

logLik = -0.5 log(2π) - 0.5 ∑_{j=0}^{n-1} log\,g(λ[j]) - π ∑_{j=0}^{n-1} I(λ[j])/g(λ[j])

where λ[j] is a frequency defined as λ[j] = 2π j/n; I(λ[j]) is the periodogram at frequency λ[j] and g(λ[j]) is the spectral generating function of the model at frequency λ[j].

The derivation of the spectral likelihood function defined above relies on the assumption that the process is circular (its covariance matrix is circulant). If the process is not circular the value of the likelihood is an approximation.

First and second order derivatives are computed by means of their analytical expressions. The first order derivatives of the spectral log-likelihood with respect to parameter θ are given by:

(d\,logLik)/(d\,θ) = 0.5 ∑_{j=0}^{n-1} ( (2π I(λ[j])) / (g(λ[j])) - 1 ) (1/g(λ[j])) d g(λ[j])/d θ

Second order derivatives are given by:

(d^2\, logLik)(d\,θ θ') = ∑_{j=0}^{n-1} ( (2π I(λ[j]))/(g(λ[j])) - 1 ) 1/(2g(λ[j])) (d^2 g(λ[j]))/(d θ d θ') -

2 ∑_{j=0}^{n-1} ( (4π I(λ[j]))(g(λ[j])) - 1 ) ( \frac{1}{2g(λ[j])} )^2 (d g(λ[j]))/(d θ) (d g(λ[j]))/(d θ')

The argument x is an auxiliary vector that is necessary in some contexts. For example, the input to function optim must contain as first argument the vector of parameters where optimization is performed. If it is not required or is redundant information contained in model@pars it can be set to NULL.

The functions mcloglik.fd, mcloglik.fd.deriv and mcloglik.fd.grad use the expressions for the spectral log-likelihood function where the parameter specified in model@cpar is concentrated out of the likelihood function.

For further information about the barrier term see Bounds on parameters and barrier term in the details section in maxlik.fd.scoring.

Arguments inf and barrier are not used by mloglik.fd.grad and mcloglik.fd.grad but they are needed in maxlik.fd.optim, where this function is passed as the gradient to be used by optim including the arguments inf and barrier.

Argument xreg. It is an optional list of constant terms. It is used by maxlik.fd.optim when analytical derivatives are employed and by maxlik.fd.scoring. It avoids computing some constant terms each time the function mloglik.fd.grad is called.

The list xreg should contain an element called dxreg, the external regressors differeneced by means of the differencing filter that renders stationarity in the model and the element fft.xreg, the Fourier transform of each regressor in dxreg.

The list xreg is not used by mloglik.fd. It is necessary to define this argument in the prototype of the function because when this function is passed to optim along with mloglik.fd.grad, the argument xreg is passed to mloglik.fd when it is defined in optim as an argument to be passed to mloglik.fd.grad.

Argument xreg is not currently implemented in functions with concentration of a parameter, mloglik.fd, mcloglik.fd.deriv and mcloglik.fd.grad.

Note: modcovgrad is not available when external regressors are defined in the input model, model.

Value

mloglik.fd returns a numeric value of the negative of the spectral log-likelihood evaluated at the parameter values defined in the input model model or at x if this argument is not NULL. If the value happens to be NA or non-finite the value of argument inf is returned. This function is suited to be passed as the objective function to optim.

mloglik.fd.deriv returns a list containing a vector of the first order derivatives of the negative of the spectral likelihood and a matrix for the second order derivatives. Those derivative terms that are not requested by setting the corresponding argument to FALSE are set to NULL in the output list.

mloglik.fd.grad returns a numeric vector containing the gradient. This function is suited to be passed as the gradient function to optim.

mcloglik.fd, mcloglik.fd.deriv and mcloglik.fd.grad return the value of the same information as the other functions but for the concentrated likelihood function.

Note

mcloglik.fd.deriv is not currently implemented for model with non-null model@transPars or with a barrier term.

References

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

See Also

barrier.eval, logLik, maxlik.fd, stsm.

Examples

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
# define the local level model for Nile time series
pars <- c("var1" = 11000, "var2" = 1700)
m <- stsm.model(model = "local-level", y = Nile, pars = pars)
# 'mloglik.fd' returns the negative of the log-likelihood function
mloglik.fd(model = m)
# 'logLik' returns the value of the log-likelihood function
logLik(object = m, domain = "frequency")

# compare analytical and numerical derivatives
# more tests in file 'test-derivatives-mloglik-fd.R' in the 
# folder 'inst' of the source package
system.time(da <- mloglik.fd.deriv(m, gradient = TRUE, hessian = TRUE))
dgn <- numDeriv::grad(func = mloglik.fd, x = m@pars, model = m)
dhn <- numDeriv::hessian(func = mloglik.fd, x = m@pars, model = m)
all.equal(as.vector(da$gradient), dgn)
all.equal(da$hessian, dhn)

# the same as above for the local level plus seasonal model and 
# a sample simulated series
data("llmseas")
m <- stsm.model(model = "llm+seas", y = llmseas)
system.time(a <- mloglik.fd.deriv(model = m, gradient = TRUE, hessian = TRUE))
system.time(g <- numDeriv::grad(func = mloglik.fd, x = m@pars, model = m))
system.time(h <- numDeriv::hessian(func = mloglik.fd, x = m@pars, model = m))
all.equal(a$gradient, g, check.attributes = FALSE)
all.equal(a$hessian, h, check.attributes = FALSE)

Want to suggest features or report bugs for rdrr.io? Use the GitHub issue tracker.