linear.model.MLE: Maximum Likelihood estimation for the geostatistical linear...

In barryrowlingson/PrevMap: Geostatistical Modelling of Spatially Referenced Prevalence Data

Description

This function performs maximum likelihood estimation for the geostatistical linear Gaussian Model.

Usage

linear.model.MLE(formula, coords, data, kappa, fixed.rel.nugget = NULL,
  start.cov.pars, method = "BFGS", low.rank = FALSE, knots = NULL,
  messages = TRUE)

Arguments

formula	an object of class "formula" (or one that can be coerced to that class): a symbolic description of the model to be fitted.
coords	an object of class formula indicating the geographic coordinates.
data	a data frame containing the variables in the model.
kappa	shape parameter of the Matern covariance function.
fixed.rel.nugget	fixed value for the relative variance of the nugget effect; default is fixed.rel.nugget=NULL if this should be included in the estimation.
start.cov.pars	a vector of length two with elements corresponding to the starting values of phi and the relative variance of the nugget effect nu2, respectively, that are used in the optimization algorithm. If nu2 is fixed through fixed.rel.nugget, then start.cov.pars represents the starting value for phi only.
method	method of optimization. If method="BFGS" then the maxBFGS function is used; otherwise method="nlminb" to use the nlminb function. Default is method="BFGS".
low.rank	logical; if low.rank=TRUE a low-rank approximation of the Gaussian spatial process is used when fitting the model. Default is low.rank=FALSE.
knots	if low.rank=TRUE, knots is a matrix of spatial knots that are used in the low-rank approximation. Default is knots=NULL.
messages	logical; if messages=TRUE then status messages are printed on the screen (or output device) while the function is running. Default is messages=TRUE.

Details

This function estimates the parameters of a geostatistical linear Gaussian model, specified as

Y = d'β + S(x) + Z,

where Y is the measured outcome, d is a vector of coavariates, β is a vector of regression coefficients, S(x) is a stationary Gaussian spatial process and Z are independent zero-mean Gaussian variables with variance tau2. More specifically, S(x) has an isotropic Matern covariance function with variance sigma2, scale parameter phi and shape parameter kappa. In the estimation, the shape parameter kappa is treated as fixed. The relative variance of the nugget effect, nu2=tau2/sigma2, can be fixed though the argument fixed.rel.nugget; if fixed.rel.nugget=NULL, then the variance of the nugget effect is also included in the estimation.

Low-rank approximation. In the case of very large spatial data-sets, a low-rank approximation of the Gaussian spatial process S(x) can be computationally beneficial. Let (x_{1},…,x_{m}) and (t_{1},…,t_{m}) denote the set of sampling locations and a grid of spatial knots covering the area of interest, respectively. Then S(x) is approximated as ∑_{i=1}^m K(\|x-t_{i}\|; φ, κ)U_{i}, where U_{i} are zero-mean mutually independent Gaussian variables with variance sigma2 and K(.;φ, κ) is the isotropic Matern kernel (see matern.kernel). Since the resulting approximation is no longer a stationary process, the parameter sigma2 is adjusted by a factorconstant.sigma2. See adjust.sigma2 for more details on the the computation of the adjustment factor constant.sigma2 in the low-rank approximation.

Value

An object of class "PrevMap". The function summary.PrevMap is used to print a summary of the fitted model. The object is a list with the following components:

estimate: estimates of the model parameters; use the function coef.PrevMap to obtain estimates of covariance parameters on the original scale.

covariance: covariance matrix of the ML estimates.

log.lik: maximum value of the log-likelihood.

y: response variable.

D: matrix of covariates.

coords: matrix of the observed sampling locations.

method: method of optimization used.

kappa: fixed value of the shape parameter of the Matern function.

knots: matrix of the spatial knots used in the low-rank approximation.

const.sigma2: adjustment factor for sigma2 in the low-rank approximation.

fixed.rel.nugget: fixed value for the relative variance of the nugget effect.

call: the matched call.

Author(s)

Emanuele Giorgi [e-mail]

Peter J. Diggle [e-mail]

References

Higdon, D. (1998). A process-convolution approach to modeling temperatures in the North Atlantic Ocean. Environmental and Ecological Statistics 5, 173-190.

See Also

shape.matern, summary.PrevMap, coef.PrevMap, matern, matern.kernel, maxBFGS, nlminb.

Examples

data(loaloa)
# Empirical logit transformation
loaloa$logit <- log((loaloa$NO_INF+0.5)/(loaloa$NO_EXAM-loaloa$NO_INF+0.5))
fit.MLE <- linear.model.MLE(logit ~ 1,coords=~LONGITUDE+LATITUDE,
                data=loaloa, start.cov.pars=c(0.2,0.15),
                 kappa=0.5)
summary(fit.MLE)

# Low-rank approximation
data(data_sim)
n.subset <- 200
data_subset <- data_sim[sample(1:nrow(data_sim),n.subset),]

# Logit transformation
data_subset$logit <- log(data_subset$y+0.5)/
                                 (data_subset$units.m-data_subset$y+0.5)
knots <- as.matrix(expand.grid(seq(-0.2,1.2,length=8),seq(-0.2,1.2,length=8)))

fit <- linear.model.MLE(formula=logit~1,coords=~x1+x2,data=data_subset,
                             kappa=2,start.cov.pars=c(0.15,0.1),low.rank=TRUE,
                             knots=knots)
summary(fit,log.cov.pars=FALSE)

barryrowlingson/PrevMap documentation built on May 11, 2019, 6:24 p.m.