# Maximum Likelihood estimation for the geostatistical linear Gaussian model

### Description

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

### Usage

1 2 3 4 |

### Arguments

`formula` |
an object of class " |

`coords` |
an object of class |

`data` |
a data frame containing the variables in the model. |

`ID.coords` |
vector of ID values for the unique set of spatial coordinates obtained from |

`kappa` |
shape parameter of the Matern covariance function. |

`fixed.rel.nugget` |
fixed value for the relative variance of the nugget effect; default is |

`start.cov.pars` |
if |

`method` |
method of optimization. If |

`low.rank` |
logical; if |

`knots` |
if |

`messages` |
logical; if |

`profile.llik` |
logical; if |

`SPDE` |
logical; if |

`mesh` |
an object obtained as result of a call to the function |

`SPDE.analytic.hessian` |
logical; if |

### 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.

**Locations with multiple observations.**
If multiple observations are available at any of the sampled locations the above model is modified as follows. Let *Y_{ij}* denote the random variable associated to the measured outcome for the j-th individual at location *x_{i}*. The linear geostatistical model assumes the form

*Y_{ij} = d_{ij}'β + S(x_{i}) + Z{i} + U_{ij},*

where *S(x_{i})* and *Z_{i}* are specified as mentioned above, and *U_{ij}* are i.i.d. zer0-mean Gaussian variable with variance *ω^2*. his model can be fitted by specifing a vector of ID for the unique set locations thourgh the argument `ID.coords`

(see also `create.ID.coords`

).

**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 factor`constant.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.

`ID.coords`

: set of ID values defined through the argument `ID.coords`

.

`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.

`mesh`

: the mesh used in the SPDE approximation.

`call`

: the matched call.

### Author(s)

Emanuele Giorgi e.giorgi@lancaster.ac.uk

Peter J. Diggle p.diggle@lancaster.ac.uk

### 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

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 | ```
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)
``` |