GeoFit2 {GeoModels}

Max-Likelihood-Based Fitting of Gaussian and non Gaussian RFs.

Description

Maximum weighted composite-likelihood fitting for Gaussian and some Non-Gaussian univariate spatial, spatio-temporal and bivariate spatial RFs. A first preliminary estimation is performed using independence composite-likelihood for the marginal parameters of the model. The estimates are then used as starting values in the second final estimation step. The function allows to fix any of the parameters and setting upper/lower bound in the optimization.

Usage

GeoFit2(data, coordx, coordy=NULL, coordt=NULL, coordx_dyn=NULL,
  copula=NULL,corrmodel,distance="Eucl",fixed=NULL,
     anisopars=NULL,est.aniso=c(FALSE,FALSE),GPU=NULL, 
     grid=FALSE, likelihood='Marginal',
     local=c(1,1),  lower=NULL,maxdist=Inf,neighb=NULL,
      maxtime=Inf, memdist=TRUE,method="cholesky", 
      model='Gaussian',n=1, onlyvar=FALSE ,
      optimizer='Nelder-Mead', parallel=FALSE, 
      radius=6371,  sensitivity=FALSE,sparse=FALSE,
       start=NULL, taper=NULL, tapsep=NULL, 
      type='Pairwise', upper=NULL, varest=FALSE, 
      vartype='SubSamp', weighted=FALSE, winconst=NULL, winstp=NULL, 
      winconst_t=NULL, winstp_t=NULL,X=NULL,nosym=FALSE,spobj=NULL,spdata=NULL)

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 spatial-temporal realisation) or an (d \times d \times t \times n )-array (a single spatial-temporal realisation on regular grid). For the description see the Section Details.

coordx

A numeric (d \times 2)-matrix (where d is the number of spatial sites) assigning 2-dimensions of spatial coordinates or a numeric d-dimensional vector assigning 1-dimension of spatial coordinates. Coordinates on a sphere for a fixed radius radius are passed in lon/lat format expressed in decimal degrees.

coordy

A numeric vector assigning 1-dimension of spatial coordinates; coordy is interpreted only if coordx is a numeric vector or grid=TRUE otherwise it will be ignored. Optional argument, the default is NULL then coordx is expected to be numeric a (d \times 2)-matrix.

coordt

A numeric vector assigning 1-dimension of temporal coordinates. Optional argument, the default is NULL then a spatial RF is expected.

coordx_dyn

A list of m numeric (d_t \times 2)-matrices containing dynamical (in time) spatial coordinates. Optional argument, the default is NULL

copula

String; the type of copula. It can be "Clayton" or "Gaussian"

corrmodel

String; the name of a correlation model, for the description see the Section Details.

distance

String; the name of the spatial distance. The default is Eucl, the euclidean distance. See the Section Details.

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 be not estimated.

anisopars

A list of two elements: "angle" and "ratio" i.e. the anisotropy angle and the anisotropy ratio, respectively.

est.aniso

A bivariate logical vector providing which anisotropic parameters must be estimated.

GPU

Numeric; if NULL (the default) no OpenCL computation is performed. The user can choose the device to be used. Use DeviceInfo() function to see available devices, only double precision devices are allowed

grid

Logical; if FALSE (the default) the data are interpreted as spatial or spatial-temporal realisations on a set of non-equispaced spatial sites (irregular grid).

likelihood

String; the configuration of the composite likelihood. Marginal is the default, see the Section Details.

local

Numeric; number of local work-items of the OpenCL setup

lower

An optional named list giving the values for the lower bound of the space parameter when the optimizer is L-BFGS-B or nlminb or optimize. The names of the list must be the same of the names in the start list.

maxdist

Numeric; an optional positive value indicating the maximum spatial distance considered in the composite or tapered likelihood computation. See the Section Details for more information.

neighb

Numeric; an optional positive integer indicating the order of neighborhood in the composite likelihood computation. See the Section Details for more information.

maxtime

Numeric; an optional positive integer indicating the order of temporal neighborhood in the composite likelihood computation.

memdist

Logical; if TRUE then all the distances useful in the composite likelihood estimation are computed before the optimization. FALSE is deprecated.

method

String; the type of matrix decomposition used in the simulation. Default is cholesky. The other possible choices is svd.

model

String; the type of RF and therefore the densities associated to the likelihood objects. Gaussian is the default, see the Section Details.

n

Numeric; number of trials in a binomial RF; number of successes in a negative binomial RF

onlyvar

Logical; if TRUE (and varest is TRUE) only the variance covariance matrix is computed without optimizing. FALSE is the default.

optimizer

String; the optimization algorithm (see optim for details). Nelder-Mead is the default. Other possible choices are nlm, BFGS, SANN, L-BFGS-B and nlminb. In these last two cases upper and lower bounds can be passed by the user. In the case of one-dimensional optimization, the function optimize is used.

parallel

Logical; if TRUE optmization is performed using optimParallel using the maximum number of cores, when optimizer is L-BFGS-B.FALSE is the default.

radius

Numeric; the radius of the sphere in the case of lon-lat coordinates. The default is 6371, the radius of the earth.

sensitivity

Logical; if TRUE then the sensitivy matrix is computed

sparse

Logical; if TRUE then maximum likelihood is computed using sparse matrices algorithms (spam packake).It should be used with compactly supported covariance models.FALSE is the default.

start

An optional named list with the initial values of the parameters that are used by the numerical routines in maximization procedure. NULL is the default (see Details).

taper

String; the name of the type of covariance matrix. It can be Standard (the default value) or Tapering for taperd covariance matrix.

tapsep

Numeric; an optional value indicating the separabe parameter in the space time adaptive taper (see Details).

type

String; the type of the likelihood objects. If Pairwise (the default) then the marginal composite likelihood is formed by pairwise marginal likelihoods (see Details).

upper

An optional named list giving the values for the upper bound of the space parameter when the optimizer is or L-BFGS-B or nlminb or optimize. The names of the list must be the same of the names in the start list.

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 the function GeoVarestbootstrap FALSE is the default.

vartype

String; (SubSamp the default) the type of method used for computing the estimates' variances, see the Section Details.

weighted

Logical; if TRUE the likelihood objects are weighted, see the Section Details. If FALSE (the default) the composite likelihood is not weighted.

winconst

Numeric; a bivariate positive vector for computing the spatial sub-window in the sub-sampling procedure. See Details for more information.

winstp

Numeric; a value in (0,1] for defining the the proportion of overlapping in the spatial sub-sampling procedure. The case 1 correspond to no overlapping. See Details for more information.

winconst_t

Numeric; a positive value for computing the temporal sub-window in the sub-sampling procedure. See Details for more information.

winstp_t

Numeric; a value in (0,1] for defining the the proportion of overlapping in the temporal sub-sampling procedure. The case 1 correspond to no overlapping. See Details for more information.

X

Numeric; Matrix of spatio(temporal)covariates in the linear mean specification.

nosym

Logical; if TRUE simmetric weights are not considered. This allows a faster but less efficient CL estimation.

spobj

An object of class sp or spacetime

spdata

Character:The name of data in the sp or spacetime object

Details

The function GeoFit2 is similar to the function GeoFit. However GeoFit2 performs a preliminary estimation using maximum indenpendence composite likelihood of the marginal parameters of the model and then use the obtained estimates as starting value in the global weighted composite likelihood estimation (that includes marginal and dependence parameters). This allows to obtain "good" starting values in the optimization algorithm for the marginal parameters.

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, if the full likelihood is considered then it coincides with the Akaike information criterion;

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 of data;

distance

The type of spatial distance;

fixed

A list of the fixed parameters;

iterations

The number of iteration 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 weigthed composite likelihood. If no spatial distance is specified then it is NULL;

maxtime

The maximum temporal distance used in the weigthed composite likelihood. If no spatial distance is specified then it is NULL;

message

Extra message passed from the numerical routines;

model

The density associated to the likelihood objects;

missp

True if a misspecified Gaussian model is ued in the composite likelihhod;

n

The number of trials in a binominal RF;the number of successes in a negative Binomial RFs;

neighb

The order of spatial neighborhood in the composite likelihood computation.

ns

The number of (different) location sites in the bivariate case;

nozero

In the case of tapered likelihood the percentage of non zero values in the covariance matrix. Otherwise is NULL.

numcoord

The number of spatial coordinates;

numtime

The number of the temporal realisations of the RF;

param

A list of the parameters' estimates;

radius

The radius of the sphere in the case of great circle distance;

stderr

The vector of standard errors for standard maximum likelihood estimation;

sensmat

The sensitivity matrix;

varcov

The matrix of the variance-covariance of the estimates;

varimat

The variability matrix;

vartype

The method used to compute the variance of the estimates;

type

The type of the likelihood objects.

winconst

The constant used to compute the window size in the spatial sub-sampling;

winstp

The step used for moving the window in the spatial sub-sampling;

winconst_t

The constant used to compute the window size in the spatio-temporal sub-sampling;

winstp_

The step used for moving the window in the spatio-temporal sub-sampling;

X

The matrix of covariates;

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

Examples

library(GeoModels)


###############################################################
############ Examples of spatial Gaussian RFs ################
###############################################################

################################################################
###
### 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 RFs with exponential correlation.
### 
###############################################################
fit1 <- GeoFit2(data=data,coordx=coords,corrmodel=corrmodel, 
                    optimizer="BFGS",neighb=3,likelihood="Conditional",
                    type="Pairwise", start=start,fixed=fixed,X=X)
print(fit1)




###############################################################
############ Examples of spatial non-Gaussian RFs #############
###############################################################


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