mle {spam}

Maximum likelihood estimates

Description

Maximum likelihood estimates of a simple spatial model

Usage

neg2loglikelihood.spam(y, X, distmat, Covariance,
                 beta, theta, Rstruct = NULL, cov.args = NULL, ...)

neg2loglikelihood(y, X, distmat, Covariance,
                 beta, theta, cov.args = NULL, ...)

neg2loglikelihood.nomean(y, distmat, Covariance,
                 theta, cov.args = NULL, ...)

mle.spam(y, X, distmat, Covariance,
     beta0, theta0, thetalower, thetaupper, optim.control=NULL,
     Rstruct = NULL, hessian = FALSE, cov.args = NULL, ...)

mle(y, X, distmat, Covariance,
     beta0, theta0, thetalower, thetaupper, optim.control=NULL,
     hessian = FALSE, cov.args = NULL, ...)

mle.nomean.spam(y, distmat, Covariance,
     theta0, thetalower, thetaupper, optim.control=NULL,
     Rstruct = NULL, hessian = FALSE, cov.args = NULL, ...) 

mle.nomean(y, distmat, Covariance,
     theta0, thetalower, thetaupper, optim.control=NULL,
     hessian = FALSE, cov.args = NULL, ...) 

Arguments

y

data vector of length n.

X

the design matrix of dimension n x p.

distmat

a distance matrix. Usually the result of a call to nearest.dist.

Covariance

function defining the covariance. See example.

beta

parameters of the trend (fixed effects).

theta

parameters of the covariance structure.

Rstruct

the Cholesky structure of the covariance matrix.

beta0, theta0

inital values.

thetalower, thetaupper

lower and upper bounds of the parameter theta.

optim.control

arguments passed to optim.

hessian

Logical. Should a numerically differentiated Hessian matrix be returned?

cov.args

additional arguments passed to Covariance.

...

additional arguments passed to chol.

Details

We provide functions to calculate the negative-2-log-likelihood and maximum likelihood estimates for the model

y ~ N_n( X beta, Sigma(h;theta) )

in the case of a sparse or ordinary covariance matrices.

In the case of the *.spam versions, the covariance function has to return a spam object. In the other case, the methods are correctly overloaded and work either way, slightly slower than the *.spam counterparts though.

When working on the sphere, the distance matrix has to be transformed by

h -> R / 2 sin(h/2)

where R is the radius of the sphere.

The covariance function requires that the first argument is the distance matrix and the second the parameters. One can image cases in which the covariance function does not take the entire distance matrix but only some partial information thereof. (An example is the use of a kronecker type covariance structure.) In case of a sparse covariance construction where the argument Rstruct is not given, the first parameter element needs to be the range parameter. (This results from the fact, that a sparse structure is constructed that is independent of the parameter values to exploit the fast Choleski decomposition.)

In the zero-mean case, the neg2loglikelihood is calculated by setting the parameters X or beta to zero.

Value

The negative-2-loglikelihood or the output from the function optim.

Author(s)

Reinhard Furrer

See Also

covmat, rmvnorm.spam

Examples

# True parameter values:
truebeta <- c(1,2,.2)    # beta = (intercept, linear in x, linear in y)
truetheta <- c(.5,2,.02) # theta = (range, sill, nugget)



# We now define a grid, distance matrix, and a sample:
x <- seq(0,1,l=5)
locs <- expand.grid( x, x)
X <- as.matrix( cbind(1,locs))  # design matrix

distmat <- nearest.dist( locs, upper=NULL) # distance matrix
Sigma <- cov.sph( distmat, truetheta)    # true covariance matrix


set.seed(15)
y <- c(rmvnorm.spam(1,X %*% truebeta,Sigma)) # construct sample

# Here is the negative 2 log likelihood:
neg2loglikelihood.spam( y, X, distmat, cov.sph,
                       truebeta, truetheta)

# We pass now to the mle:
res <- mle.spam(y, X, distmat, cov.sph,
         truebeta, truetheta,thetalower=c(0,0,0),thetaupper=c(1,Inf,Inf))

# Similar parameter estimates here, of course:
mle.nomean.spam(y-X%*%res$par[1:3], distmat, cov.sph,
         truetheta, thetalower=c(0,0,0), thetaupper=c(1,Inf,Inf))