GeoFit2: Fitting Gaussian and Non-Gaussian Random Fields with...

In GeoModels: Procedures for Gaussian and Non Gaussian Geostatistical (Large) Data Analysis

Fitting Gaussian and Non-Gaussian Random Fields with Automatic Marginal Initialization

Description

A univariate starting-value wrapper around GeoFit. The function optionally performs a preliminary fit under spatial independence to estimate available marginal parameters, uses those estimates to update the starting values, and then calls GeoFit for the requested final fit. If the preliminary independence fit is unavailable or fails, the original starting values are retained. Bivariate correlation models are not supported.

Usage

GeoFit2(data, coordx = NULL, coordy = NULL, coordz = NULL, coordt = NULL,
 coordx_dyn = NULL, copula = NULL, corrmodel = NULL,
 distance = "Eucl", fixed = NULL, anisopars = NULL,
 est.aniso = c(FALSE, FALSE), grid = FALSE,
 likelihood = "Marginal", lower = NULL, maxdist = Inf,
 neighb = NULL, p_neighb = 1, maxtime = Inf, memdist = TRUE,
 method = "cholesky", model = "Gaussian", n = 1,
 onlyvar = FALSE, optimizer = "Nelder-Mead", radius = 1,
 score = FALSE, sensitivity = FALSE, sparse = FALSE,
 start = NULL, thin_method = "bernoulli", type = "Pairwise",
 upper = NULL, varest = FALSE, weighted = FALSE, X = NULL,
 spobj = NULL, spdata = NULL, independence_start = TRUE,
 independence_optimizer = "Nelder-Mead",
 warn_independence_failure = TRUE, check.duplicates = FALSE)

Arguments

data	A d-dimensional vector (a single spatial realisation) or a (d \times d)-matrix (a single spatial realisation on regular grid) or a (t \times d)-matrix (a single spatio-temporal realisation) or an (d \times d \times t \times n)-array (a single spatio-temporal realisation on regular grid). See Details for the accepted data layouts.
coordx	A numeric (d \times 2)-matrix or (d \times 3)-matrix. Coordinates on a sphere for a fixed radius radius are passed in lon/lat format expressed in decimal degrees. It may be omitted when coordinates are supplied through coordx_dyn or spobj.
coordy	A numeric vector giving 1-dimension of spatial coordinates; optional argument, default is NULL.
coordz	A numeric vector giving 1-dimension of spatial coordinates; optional argument, default is NULL.
coordt	A numeric vector assigning one dimension of the observation-time coordinates. Optional argument, default is NULL; if NULL, a spatial random field is expected. Temporal coordinates may be irregularly spaced; temporal lags are computed from the supplied coordinate values.
coordx_dyn	A list of m numeric (d_t \times 2)-matrices containing dynamical (in time) spatial coordinates. Optional argument, default is NULL.
copula	String; the type of copula. It can be "Clayton" or "Gaussian".
corrmodel	String; the name of a correlation model; see GeoCovmatrix for the list of implemented correlation models.
distance	String; the name of the spatial distance. Default is "Eucl" (Euclidean distance). See Details for the accepted options.
fixed	An optional named list giving the values of the parameters that will be considered as known values. The listed parameters for a given correlation function will not be estimated.
anisopars	A list of two elements: "angle" and "ratio", i.e. the anisotropy angle and the anisotropy ratio, respectively. Geometric anisotropy requires distance="Eucl".
est.aniso	A bivariate logical vector providing which anisotropic parameters must be estimated.
grid	Logical; if FALSE (default) the data are interpreted as spatial or spatio-temporal realisations on a set of non-equispaced spatial sites. If TRUE, the grid axes and data dimensions are handled by GeoFit and converted internally to the canonical explicit-coordinate layout before fitting.
likelihood	String; the configuration of the composite likelihood. "Marginal" is the default; see Details for the accepted options.
lower	An optional named list giving lower bounds for parameters when the optimizer is L-BFGS-B, nlminb, bobyqa or optimize. Names must match parameters that are being estimated.
maxdist	Numeric; an optional positive value indicating the maximum spatial distance considered in the composite likelihood computation. See Details for more information.
neighb	Numeric; an optional positive integer indicating the order of neighborhood in the composite likelihood computation. See Details for more information.
p_neighb	Numeric scalar in (0,1]. If 1 (default), no thinning is applied. For thin_method="bernoulli", it is the candidate-pair inclusion probability in the constant-weight design and controls the retained size in expectation. For thin_method="FixedBudget", it defines the exact global retained-pair budget K=\mathrm{round}(p_{neighb}d), where d is the candidate-pair count.
maxtime	Numeric; an optional non-negative maximum temporal-distance threshold, expressed in the same units as coordt, used to select pairwise composite-likelihood contributions.
memdist	Deprecated logical argument retained for backward compatibility. The selected pair structure is always precomputed and reused during composite-likelihood optimization. Supplying FALSE produces a warning and is treated as TRUE.
method	String; the type of matrix decomposition used in the likelihood computation. Default is "cholesky". Another possible choice is "svd" (when available).
model	String; the type of RF and therefore the densities associated to the likelihood objects. "Gaussian" is the default; see Details for the accepted options.
n	Positive integer scalar, or one positive integer per observation, for models using a Binomial/Negative-Binomial trial or success count.
onlyvar	Logical; if TRUE (and varest=TRUE) only the variance-covariance matrix is computed without optimizing. Default is FALSE.
optimizer	String; the optimization algorithm (see optim for details). "Nelder-Mead" is the default. Other possible choices are "nlm", "BFGS", "SANN", "L-BFGS-B", "nlminb", "bobyqa". In these last three cases upper and lower bounds can be passed by the user. In the one-dimensional case, optimize is used.
radius	Numeric; the radius of the sphere in the case of lon-lat coordinates. Default value is 1.
score	Logical; should score function be computed? Default is FALSE.
sensitivity	Logical; if TRUE then the sensitivity matrix is computed.
sparse	Logical; if TRUE then maximum likelihood is computed using sparse matrix algorithms (e.g., spam). It should be used with compactly supported covariance models. Default is FALSE.
start	An optional named list with initial values for parameters used by the numerical routines in the maximization procedure. Default is NULL; omitted parameters receive internal starting values (see Details).
thin_method	String; thinning scheme used when p_neighb < 1. Default is "bernoulli" (independent Bernoulli thinning). An alternative thinning scheme is "FixedBudget". The legacy name "TargetBalanced" is accepted as an alias for "FixedBudget".
type	String; the type of the likelihood objects. If "Pairwise" (default) then the marginal composite likelihood is formed by pairwise marginal likelihoods (see Details).
upper	An optional named list giving upper bounds for parameters when the optimizer is L-BFGS-B, nlminb, bobyqa or optimize. Names must match parameters that are being estimated.
varest	Logical; if TRUE the estimates' variances and standard errors are returned. For composite likelihood estimation it is deprecated. Use sensitivity=TRUE and update the object using GeoVarest. Default is FALSE.
weighted	Logical; if TRUE the likelihood objects are weighted; see Details for the accepted options. Default is FALSE.
X	Numeric; matrix of spatio(temporal) covariates in the linear mean specification.
spobj	An object of class sp or spacetime.
spdata	Character; the name of data in the sp or spacetime object.
independence_start	Logical; if TRUE (default), GeoFit2 attempts a preliminary marginal fit under spatial independence before the final fit. The preliminary fit is skipped when the requested fit already has type="Independence". Bivariate correlation models are not supported by GeoFit2; use GeoFit directly for those models.
independence_optimizer	String; optimizer used only for the preliminary independence fit. Default is "Nelder-Mead". The optimizer used for the final fit is still specified by optimizer.
warn_independence_failure	Logical; if TRUE (default), a warning is issued when the preliminary independence fit is unavailable or fails. In that case the supplied start values are used unchanged.
check.duplicates	Logical. If TRUE, perform a fast exact scan for duplicated observation locations at this user entry point. The default is FALSE, so no duplicate-location scan is imposed. Internal bootstrap, cross-validation, and refitting calls do not repeat the scan.

Details

GeoFit2 does not implement a second fitting engine. It is a wrapper around the canonical GeoFit function.

When independence_start=TRUE, the function first constructs a preliminary call to GeoFit with likelihood="Marginal" and type="Independence". The preliminary fit estimates the marginal parameters available for the selected model. When start is supplied, matching marginal entries and omitted mean coefficients are updated as before. When start=NULL, all eligible marginal estimates from the independence fit are used, while dependence parameters use the conservative automatic initial values prepared by GeoFit. Parameters supplied in fixed always take precedence.

The final model is then fitted by calling GeoFit with the original requested likelihood, likelihood-object type, optimizer, bounds, pair-selection settings, and the updated starting values. Consequently, the objective function, parameter constraints, and returned standard GeoFit components are the same as in a direct call to GeoFit; only the starting values may differ.

For Gaussian models, sill is the total marginal variance and nugget attenuates off-diagonal correlation. Therefore an independence estimate may initialize sill, while a user-supplied nugget is preserved unchanged and is never added to or subtracted from sill.

Bivariate correlation models are not supported by GeoFit2; use GeoFit directly. The preliminary initialization is skipped when type="Independence". If an independence likelihood is not implemented for the selected model, or if the preliminary optimization fails, the final fit still proceeds using the supplied starting values or, when start=NULL, the internal starting values prepared by GeoFit. This behavior can be controlled with warn_independence_failure.

Stochastic thinning of nearest-neighbor pairs in the final fit can be enabled through p_neighb < 1; thin_method specifies the thinning scheme.

Value

Returns an object of class GeoFit. An object of class GeoFit is a list containing at most the following components:

bivariate	Logical: TRUE if the Gaussian RF is bivariate, otherwise FALSE.
clic	The composite information criterion after a GeoVarest call; if the full likelihood is considered then it coincides with AIC.
coordx	A d-dimensional vector of spatial coordinates.
coordy	A d-dimensional vector of spatial coordinates.
coordt	A t-dimensional vector of temporal coordinates.
coordx_dyn	A list of dynamical (in time) spatial coordinates.
conf.int	Confidence intervals for standard maximum likelihood estimation.
convergence	A string that denotes if convergence is reached.
copula	The type of copula.
corrmodel	The correlation model.
data	The vector or matrix or array (or list) of data.
distance	The type of spatial distance.
fixed	A list of fixed parameters.
iterations	The number of iterations used by the numerical routine.
likelihood	The configuration of the composite likelihood.
logCompLik	The value of the log composite-likelihood at the maximum.
maxdist	The maximum spatial distance used in the weighted composite likelihood (or NULL).
maxtime	The maximum temporal-distance threshold used in the composite likelihood, expressed in the same units as coordt.
message	Extra message passed from the numerical routines.
model	The density associated to the likelihood objects.
missp	TRUE if a misspecified Gaussian model is used in the composite likelihood.
n	The number of trials in a binomial RF; the number of successes in a negative binomial random field.
neighb	The order of spatial neighborhood in the composite likelihood computation.
ns	The number of (different) location sites in the bivariate case.
numcoord	The number of spatial coordinates.
numtime	The number of temporal realisations of the random field.
param	A list of parameter estimates.
radius	The radius of the sphere in the case of great-circle distance.
stderr	Standard errors for standard maximum likelihood estimation.
sensmat	The sensitivity matrix.
varcov	The variance-covariance matrix of the estimates.
type	The type of the likelihood objects.
X	The matrix of covariates.
start_original	The named list of starting values supplied by the user.
start_used	The named list of starting values passed to the final GeoFit call after any successful marginal initialization.
independence_fit	The preliminary independence GeoFit object, or NULL when the preliminary fit was skipped or failed.
independence_start_message	A character string describing why the preliminary initialization was unavailable or failed, or NULL when no such message was generated.

Author(s)

Moreno Bevilacqua, [e-mail], https://sites.google.com/view/moreno-bevilacqua/home, Víctor Morales Oñate, [e-mail], https://sites.google.com/site/moralesonatevictor/, Christian Caamaño-Carrillo, [e-mail], https://www.researchgate.net/profile/Christian-Caamano

References

The methodological references for maximum weighted composite-likelihood fitting of Gaussian and non-Gaussian random fields are reported in GeoFit.

See Also

GeoCovmatrix for covariance matrix construction, GeoSim for simulation, GeoKrig for prediction, GeoVariogram and GeoWLS for variogram-based tools, GeoVarest for bootstrap variance estimation.

Examples

library(GeoModels)


###############################################################
############ Examples of spatial Gaussian random fields ################
###############################################################

################################################################
###
### Example 1 : Maximum pairwise conditional likelihood fitting 
### of a Gaussian RF with Matern correlation
###
###############################################################
model="Gaussian"
# Define the spatial-coordinates of the points:
set.seed(3)
N=400 # number of location sites
x <- runif(N, 0, 1)
set.seed(6)
y <- runif(N, 0, 1)
coords <- cbind(x,y)

# Define spatial matrix covariates
X=cbind(rep(1,N),runif(N))

# Set the covariance model's parameters:
corrmodel <- "Matern"
mean <- 0.2
mean1 <- -0.5
sill <- 1
nugget <- 0
scale <- 0.2/3
smooth=0.5
param<-list(mean=mean,mean1=mean1,sill=sill,nugget=nugget,scale=scale,smooth=smooth)

# Simulation of the spatial Gaussian RF:
data <- GeoSim(coordx=coords,model=model,corrmodel=corrmodel, param=param,X=X)$data

fixed<-list(nugget=nugget,smooth=smooth)
start<-list(mean=mean,mean1=mean1,scale=scale,sill=sill)

################################################################
###
### Maximum pairwise likelihood fitting of
### Gaussian random fields with exponential correlation.
### 
###############################################################
fit1 <- GeoFit2(data=data,coordx=coords,corrmodel=corrmodel, 
 neighb=3,likelihood="Conditional",
 type="Pairwise", start=start,fixed=fixed,X=X)
print(fit1)




###############################################################
############ Examples of spatial non-Gaussian random fields #############
###############################################################


################################################################
###
### Example 2. Maximum pairwise likelihood fitting of 
### a LogGaussian RF with Generalized Wendland correlation
### 
###############################################################
set.seed(524)
# Define the spatial-coordinates of the points:
N=500
x <- runif(N, 0, 1)
y <- runif(N, 0, 1)
coords <- cbind(x,y)
X=cbind(rep(1,N),runif(N))
mean=1; mean1=2 # regression parameters
nugget=0
sill=0.5
scale=0.2
smooth=0

model="LogGaussian"
corrmodel="GenWend"
param=list(mean=mean,mean1=mean1,sill=sill,scale=scale,
 nugget=nugget,power2=4,smooth=smooth)
# Simulation of a non stationary LogGaussian RF:
data <- GeoSim(coordx=coords, corrmodel=corrmodel,model=model,X=X,
 param=param)$data

fixed<-list(nugget=nugget,power2=4,smooth=smooth)
start<-list(mean=mean,mean1=mean1,scale=scale,sill=sill)
I=Inf
lower<-list(mean=-I,mean1=-I,scale=0,sill=0)
upper<-list(mean= I,mean1= I,scale=I,sill=I)

# Maximum pairwise composite-likelihood fitting of the RF:
fit <- GeoFit2(data=data,coordx=coords,corrmodel=corrmodel, model=model,
 neighb=3,likelihood="Conditional",type="Pairwise",X=X,
 optimizer="nlminb",lower=lower,upper=upper,
 start=start,fixed=fixed)
print(unlist(fit$param))


################################################################
###
### Example 3. Maximum pairwise likelihood fitting of
### SinhAsinh random fields with Wendland0 correlation
###
###############################################################
set.seed(261)
model="SinhAsinh"
# Define the spatial-coordinates of the points:
x <- runif(500, 0, 1)
y <- runif(500, 0, 1)
coords <- cbind(x,y)

corrmodel="Wend0"
mean=0;nugget=0
sill=1
skew=-0.5
tail=1.5
power2=4
c_supp=0.2

# model parameters
param=list(power2=power2,skew=skew,tail=tail,
 mean=mean,sill=sill,scale=c_supp,nugget=nugget)
data <- GeoSim(coordx=coords, corrmodel=corrmodel,model=model, param=param)$data

plot(density(data))
fixed=list(power2=power2,nugget=nugget)
start=list(scale=c_supp,skew=skew,tail=tail,mean=mean,sill=sill)
# Maximum pairwise likelihood:
fit1 <- GeoFit2(data=data,coordx=coords,corrmodel=corrmodel, model=model,
 neighb=3,likelihood="Marginal",type="Pairwise",
 start=start,fixed=fixed)
print(unlist(fit1$param))

GeoModels documentation built on Sept. 23, 2026, 5:07 p.m.