fit_mlgcp {stelfi}

Marked spatial log-Gaussian Cox process (mLGCP)

Description

Fit a marked LGCP using Template Model Builder (TMB) and the R_inla namespace for the SPDE-based construction of the latent field.

Usage

fit_mlgcp(
  locs,
  sf,
  marks,
  smesh,
  parameters = list(),
  methods,
  strfixed = matrix(1, nrow = nrow(locs), ncol = ncol(marks)),
  fields = rep(1, ncol(marks)),
  covariates,
  pp_covariates,
  marks_covariates,
  tmb_silent = TRUE,
  nlminb_silent = TRUE,
  ...
)

Arguments

locs

A data.frame of x and y locations, 2xn.

sf

An sf of type POLYGON specifying the spatial region of the domain.

marks

A matrix of marks for each observation of the point pattern.

smesh

A Delaunay triangulation of the spatial domain returned by fmesher::fm_mesh_2d().

parameters

a list of named parameters: log_tau, log_kappa, betamarks, betapp, marks_coefs_pp.

methods

An integer value:

• 0 (default), Gaussian distribution, parameter estimated is mean;

• 1, Poisson distribution, parameter estimated is intensity;

• 2, binomial distribution, parameter estimated is logit/probability;

• 3, gamma distribution, the implementation in TMB is shape-scale.

strfixed

A matrix of fixed structural parameters, defined for each event and mark. Defaults to 1. If mark distribution

• Normal, then this is the log of standard deviation;

• Poisson, then not used;

• Binomial, then this is the number of trials;

• Gamma, then this is the log of the scale.

fields

A binary vector indicating whether there is a new random field for each mark. By default, each mark has its own random field.

covariates

Covariate(s) corresponding to each area in the spatial mesh

pp_covariates

Which columns of the covariates apply to the point process

marks_covariates

Which columns of the covariates apply to the marks. By default, all covariates apply to the marks only.

tmb_silent

Logical, if TRUE (default) then TMB inner optimisation tracing information will be printed.

nlminb_silent

Logical, if TRUE (default) then for each iteration nlminb() output will be printed.

...

optional extra arguments to pass into stats::nlminb().

Details

The random intensity surface of the point process is (as fit_lgcp) \Lambda(\boldsymbol{x}) = \textrm{exp}(\boldsymbol{X}\beta + G(\boldsymbol{x}) + \epsilon), for design matrix \boldsymbol{X}, coefficients \boldsymbol{\beta}, and random error \epsilon.

Each mark, m_j, is jointly modelled and has their own random field M_j(s) = f^{-1}((\boldsymbol{X}\beta)_{m_j} + G_{m_j}(\boldsymbol{x}) + \alpha_{m_j}\; G(\boldsymbol{x}) + \epsilon_{m_j}) where \alpha_{.} are coefficient(s) linking the point process and the mark(s).

M_j(s) depends on the distribution of the marks. If the marks are from a Poisson distribution, it is the intensity (as with the point process). If the marks are from a Binomial distribution, it is the success probability, and the user must supply the number of trials for each event (via strfixed). If the marks are normally distributed then this models the mean, and the user must supply the standard deviation (via strfixed). The user can choose for the point processes and the marks to share a common GMRF, i.e. G_m(s) = G_{pp}(s); this is controlled via the argument fields.

Value

A list containing components of the fitted model, see TMB::MakeADFun. Includes

• par, a numeric vector of estimated parameter values;

• objective, the objective function; and

• gr, the TMB calculated gradient function.

References

Lindgren, F., Rue, H., and Lindström, J. (2011) An explicit link between Gaussian fields and Gaussian Markov random fields: the stochastic partial differential equation approach. Journal of the Royal Statistical Society: Series B (Statistical Methodology), 73: 423–498.

See Also

fit_lgcp

Examples

### ********************** ###
## A joint likelihood marked LGCP model
### ********************** ###
if(requireNamespace("fmesher")){
data(marked, package = "stelfi")
loc.d <- 3 * cbind(c(0, 1, 1, 0, 0), c(0, 0, 1, 1, 0))
domain <- sf::st_sf(geometry = sf::st_sfc(sf::st_polygon(list(loc.d))))
smesh <- fmesher::fm_mesh_2d(loc.domain = loc.d, offset = c(0.3, 1),
max.edge = c(0.3, 0.7), cutoff = 0.05)
locs <- cbind(x = marked$x, y = marked$y)
marks <- cbind(m1 = marked$m1) ## Gaussian mark
parameters <- list(betamarks = matrix(0, nrow = 1, ncol = ncol(marks)),
log_tau = rep(log(1), 2), log_kappa = rep(log(1), 2),
marks_coefs_pp = rep(0, ncol(marks)), betapp = 0)
fit <- fit_mlgcp(locs = locs, marks = marks,
sf = domain, smesh = smesh,
parameters = parameters, methods = 0,fields = 1)
}