GeoKrig {GeoModels} | R Documentation |
Spatial (bivariate) and spatio temporal optimal linear prediction for Gaussian and non Gaussian random fields.
Description
For a given set of spatial location sites (and temporal instants), the function computes optimal linear prediction and associated mean square error for the Gaussian and non Gaussian case.
Usage
GeoKrig(estobj=NULL,data, coordx, coordy=NULL, coordt=NULL,
coordx_dyn=NULL, corrmodel,distance="Eucl",
grid=FALSE, loc, maxdist=NULL, maxtime=NULL,
method="cholesky", model="Gaussian", n=1,nloc=NULL,mse=FALSE,
lin_opt=TRUE, param, anisopars=NULL,radius=6371, sparse=FALSE,
taper=NULL,tapsep=NULL, time=NULL, type="Standard",type_mse=NULL,
type_krig="Simple",weigthed=TRUE,which=1,
copula=NULL, X=NULL,Xloc=NULL,Mloc=NULL,spobj=NULL,spdata=NULL)
Arguments
estobj |
An object of class Geofit that includes information about data, model and estimates. |
data |
A |
coordx |
A numeric ( |
coordy |
A numeric vector giving 1-dimension of
spatial coordinates used for prediction; |
coordt |
A numeric vector giving 1-dimension of
temporal coordinates used for prediction; the default is |
coordx_dyn |
A list of |
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 |
grid |
Logical; if |
lin_opt |
Logical;If TRUE (default) then optimal (pairwise) linear kriging is computed. Otherwise optimal (pairwise) kriging is computed in the mean square sense. |
loc |
A numeric ( |
maxdist |
Numeric; an optional positive value indicating the maximum spatial compact support in the case of covariance tapering kriging. |
maxtime |
Numeric; an optional positive value indicating the maximum temporal compact support in the case of covasriance tapering kriging. |
method |
String; the type of matrix decomposition used in the simulation. Default is |
n |
Numeric; the number of trials in a binomial random fields.
Default is |
nloc |
Numeric; the number of trials of the locations sites to be predicted in a binomial random fields type II.
Default is |
mse |
Logical; if |
model |
String; the type of RF and therefore the densities associated to the likelihood
objects. |
param |
A list of parameter values required for the correlation model.See the Section Details. |
anisopars |
A list of two elements: "angle" and "ratio" i.e. the anisotropy angle and the anisotropy ratio, respectively. |
radius |
Numeric: the radius of the sphere if coordinates are passed in lon/lat format; |
sparse |
Logical; if |
taper |
String; the name of the taper correlation function, see the Section Details. |
tapsep |
Numeric; an optional value indicating the separabe parameter in the space time quasi taper (see Details). |
time |
A numeric ( |
type |
String; if |
type_mse |
String; if |
type_krig |
String; the type of kriging. If |
weigthed |
Logical; if |
which |
Numeric; In the case of bivariate (tapered) cokriging it indicates which variable to predict. It can be 1 or 2 |
copula |
String; the type of copula. It can be "Clayton" or "Gaussian" |
X |
Numeric; Matrix of spatio(temporal)covariates in the linear mean specification. |
Xloc |
Numeric; Matrix of spatio(temporal)covariates in the linear mean specification associated to predicted locations. |
Mloc |
Numeric; Vector of spatio(temporal) estimated means associated to predicted locations. |
spobj |
An object of class sp or spacetime |
spdata |
Character:The name of data in the sp or spacetime object |
Details
Best linear unbiased predictor and associated mean square error is computed
for Gaussian and some non Gaussian cases.
Specifically, for a spatial or spatio-temporal or spatial bivariate dataset, given a set of spatial locations and
temporal istants and a correlation model
corrmodel
with some fixed parameters and given the type of RF (model
) the function computes
simple kriging, for the specified spatial locations
loc
and temporal instants time
, providing also the respective mean square error.
For the choice of the spatial or spatio temporal correlation model see details in GeoCovmatrix
function.
The list param
specifies mean and covariance parameters, see CorrParam
and GeoCovmatrix
for details. The type_krig
parameter indicates the type of kriging. In the
case of simple kriging, the known mean can be specified by the parameter
mean
in the list param
(See examples).
Value
Returns an object of class Kg
.
An object of class Kg
is a list containing
at most the following components:
bivariate |
|
coordx |
A |
coordy |
A |
coordt |
A |
corrmodel |
String: the correlation model; |
covmatrix |
The covariance matrix if |
data |
The vector or matrix or array of data used for prediction |
distance |
String: the type of spatial distance; |
grid |
|
loc |
A ( |
n |
The number of trial for Binomial RFs |
nozero |
In the case of tapered simple kriging the percentage of non zero values in the covariance matrix. Otherwise is NULL. |
numcoord |
Numeric:he number |
numloc |
Numeric: the number |
numtime |
Numeric: the number |
numt |
Numeric: the number |
model |
The type of RF, see |
param |
Numeric: The covariance parameters; |
pred |
A ( |
radius |
Numeric: the radius of the sphere if coordinates are pssed in lon/lat format; |
spacetime |
|
tapmod |
String: the taper model if |
time |
A |
type |
String: the type of kriging (Standard or Tapering). |
type_krig |
String: the type of kriging. |
mse |
A ( |
Author(s)
Moreno Bevilacqua, moreno.bevilacqua89@gmail.com,https://sites.google.com/view/moreno-bevilacqua/home, Víctor Morales Oñate, victor.morales@uv.cl, https://sites.google.com/site/moralesonatevictor/, Christian", Caamaño-Carrillo, chcaaman@ubiobio.cl,https://www.researchgate.net/profile/Christian-Caamano
References
Gaetan, C. and Guyon, X. (2010) Spatial Statistics and Modelling. Spring Verlang, New York.
See Also
Examples
library(GeoModels)
################################################################
########### Examples of spatial kriging ############
################################################################
################################################################
###
### Example 1. Spatial kriging of a
### Gaussian random fields with Gen wendland correlation.
###
################################################################
model="Gaussian"
set.seed(79)
x = runif(300, 0, 1)
y = runif(300, 0, 1)
coords=cbind(x,y)
# Set the exponential cov parameters:
corrmodel = "GenWend"
mean=0; sill=5; nugget=0
scale=0.2;smooth=0;power2=4
param=list(mean=mean,sill=sill,nugget=nugget,scale=scale,smooth=smooth,power2=power2)
# Simulation of the spatial Gaussian random field:
data = GeoSim(coordx=coords, corrmodel=corrmodel,
param=param)$data
## estimation with pairwise likelihood
fixed=list(nugget=nugget,smooth=0,power2=power2)
start=list(mean=0,scale=scale,sill=1)
I=Inf
lower=list(mean=-I,scale=0,sill=0)
upper=list(mean= I,scale=I,sill=I)
# Maximum pairwise likelihood fitting :
fit = GeoFit(data, coordx=coords, corrmodel=corrmodel,model=model,
likelihood='Marginal', type='Pairwise',neighb=3,
optimizer="nlminb", lower=lower,upper=upper,
start=start,fixed=fixed)
# locations to predict
xx=seq(0,1,0.03)
loc_to_pred=as.matrix(expand.grid(xx,xx))
## first option
#param=append(fit$param,fit$fixed)
#pr=GeoKrig(loc=loc_to_pred,coordx=coords,corrmodel=corrmodel,
# model=model,param=param,data=data,mse=TRUE)
## second option using object GeoFit
pr=GeoKrig(fit,loc=loc_to_pred,mse=TRUE)
colour = rainbow(100)
opar=par(no.readonly = TRUE)
par(mfrow=c(1,3))
quilt.plot(coords,data,col=colour)
# simple kriging map prediction
image.plot(xx, xx, matrix(pr$pred,ncol=length(xx)),col=colour,
xlab="",ylab="",
main=" Kriging ")
# simple kriging MSE map prediction variance
image.plot(xx, xx, matrix(pr$mse,ncol=length(xx)),col=colour,
xlab="",ylab="",main="Std error")
par(opar)
################################################################
###
### Example 2. Spatial kriging of a Skew
### Gaussian random fields with Matern correlation.
###
################################################################
model="SkewGaussian"
set.seed(79)
x = runif(300, 0, 1)
y = runif(300, 0, 1)
coords=cbind(x,y)
# Set the exponential cov parameters:
corrmodel = "Matern"
mean=0
sill=2
nugget=0
scale=0.1
smooth=0.5
skew=3
param=list(mean=mean,sill=sill,nugget=nugget,scale=scale,smooth=smooth,skew=skew)
# Simulation of the spatial skew Gaussian random field:
data = GeoSim(coordx=coords, corrmodel=corrmodel,model=model,
param=param)$data
fixed=list(nugget=nugget,smooth=smooth)
start=list(mean=0,scale=scale,sill=1,skew=skew)
I=Inf
lower=list(mean=-I,scale=0,sill=0,skew=-I)
upper=list(mean= I,scale=I,sill=I,skew=I)
# Maximum pairwise likelihood fitting :
fit = GeoFit2(data, coordx=coords, corrmodel=corrmodel,model=model,
likelihood='Marginal', type='Pairwise',neighb=3,
optimizer="nlminb", lower=lower,upper=upper,
start=start,fixed=fixed)
# locations to predict
xx=seq(0,1,0.03)
loc_to_pred=as.matrix(expand.grid(xx,xx))
## optimal linear kriging
pr=GeoKrig(fit,loc=loc_to_pred,mse=TRUE)
colour = rainbow(100)
opar=par(no.readonly = TRUE)
par(mfrow=c(1,3))
quilt.plot(coords,data,col=colour)
# simple kriging map prediction
image.plot(xx, xx, matrix(pr$pred,ncol=length(xx)),col=colour,
xlab="",ylab="",
main=" Kriging ")
# simple kriging MSE map prediction variance
image.plot(xx, xx, matrix(pr$mse,ncol=length(xx)),col=colour,
xlab="",ylab="",main="Std error")
par(opar)
################################################################
###
### Example 3. Spatial kriging of a
### Gamma random field with mean spatial regression
###
###############################################################
set.seed(312)
model="Gamma"
corrmodel = "GenWend"
# Define the spatial-coordinates of the points:
NN=300
coords=cbind(runif(NN),runif(NN))
## matrix covariates
a0=rep(1,NN)
a1=runif(NN,0,1)
X=cbind(a0,a1)
##Set model parameters
shape=2
## regression parameters
mean = 1;mean1= -0.2
# correlation parameters
nugget = 0;power2=4
scale = 0.3;smooth=0
## simulation
param=list(shape=shape,nugget=nugget,mean=mean,mean1=mean1,
scale=scale,power2=power2,smooth=smooth)
data = GeoSim(coordx=coords,corrmodel=corrmodel, param=param,
model=model,X=X)$data
#####starting and fixed parameters
fixed=list(nugget=nugget,power2=power2,smooth=smooth)
start=list(mean=mean,mean1=mean1, scale=scale,shape=shape)
## estimation with pairwise likelihood
fit2 = GeoFit(data=data,coordx=coords,corrmodel=corrmodel,X=X,
neighb=3,likelihood="Conditional",type="Pairwise",
start=start,fixed=fixed, model = model)
# locations to predict with associated covariates
xx=seq(0,1,0.03)
loc_to_pred=as.matrix(expand.grid(xx,xx))
NP=nrow(loc_to_pred)
a0=rep(1,NP)
a1=runif(NP,0,1)
Xloc=cbind(a0,a1)
#optimal linear kriging
pr=GeoKrig(fit2,loc=loc_to_pred,Xloc=Xloc,sparse=TRUE,mse=TRUE)
## map
opar=par(no.readonly = TRUE)
par(mfrow=c(1,3))
quilt.plot(coords,data,main="Data")
map=matrix(pr$pred,ncol=length(xx))
mapmse=matrix(pr$mse,ncol=length(xx))
image.plot(xx, xx, map,
xlab="",ylab="",main="Kriging ")
image.plot(xx, xx, mapmse,
xlab="",ylab="",main="MSE")
par(opar)
################################################################
########### Examples of spatio temporal kriging ############
################################################################
################################################################
###
### Example 4. Spatio temporal simple kriging of n locations
### sites and m temporal instants for a Gaussian random fields
### with estimated double Wendland correlation.
###
###############################################################
model="Gaussian"
# Define the spatial-coordinates of the points:
x = runif(300, 0, 1)
y = runif(300, 0, 1)
coords=cbind(x,y)
times=1:4
# Define model correlation modek and associated parameters
corrmodel="Wend0_Wend0"
param=list(nugget=0,mean=0,power2_s=4,power2_t=4,
scale_s=0.2,scale_t=2,sill=1)
# Simulation of the space time Gaussian random field:
set.seed(31)
data=GeoSim(coordx=coords,coordt=times,corrmodel=corrmodel,sparse=TRUE,
param=param)$data
# Maximum pairwise likelihood fitting of the space time random field:
start = list(scale_s=0.15,scale_t=2,sill=1,mean=0)
fixed = list(nugget=0,power2_s=4,power2_t=4)
fit = GeoFit(data, coordx=coords, coordt=times, model=model, corrmodel=corrmodel,
likelihood='Conditional', type='Pairwise',start=start,fixed=fixed,
neighb=3,maxtime=1)
# locations to predict
xx=seq(0,1,0.04)
loc_to_pred=as.matrix(expand.grid(xx,xx))
# Define the times to predict
times_to_pred=2
pr=GeoKrig(fit,loc=loc_to_pred,time=times_to_pred,sparse=TRUE,mse=TRUE)
opar=par(no.readonly = TRUE)
par(mfrow=c(1,3))
zlim=c(-2.5,2.5)
colour = rainbow(100)
quilt.plot(coords,data[2,] ,col=colour,main = paste(" data at Time 2"))
image.plot(xx, xx, matrix(pr$pred,ncol=length(xx)),col=colour,
main = paste(" Kriging at Time 2"),ylab="")
image.plot(xx, xx, matrix(pr$mse,ncol=length(xx)),col=colour,
main = paste("Std err Time at time 2"),ylab="")
par(opar)
################################################################
###
### Example r. Spatial bivariate simple cokriging of n locations
### sites for a bivariate Gaussian random fields
### with estimated Matern correlation.
###
###############################################################
#set.seed(6)
#NN=1500 # number of spatial locations
#x = runif(NN, 0, 1);
#y = runif(NN, 0, 1)
#coords=cbind(x,y)
## setting parameters
#mean_1 = 2; mean_2= -1
#nugget_1 =0;nugget_2=0
#sill_1 =0.5; sill_2 =1;
### correlation parameters
#CorrParam("Bi_Matern")
#scale_1=0.2/3; scale_2=0.15/3; scale_12=0.5*(scale_2+scale_1)
#smooth_1=smooth_2=smooth_12=0.5
#pcol = -0.4
#param= list(nugget_1=nugget_1,nugget_2=nugget_2,
# sill_1=sill_1,sill_2=sill_2,
# mean_1=mean_1,mean_2=mean_2,
# smooth_1=smooth_1, smooth_2=smooth_2,smooth_12=smooth_12,
# scale_1=scale_1, scale_2=scale_2,scale_12=scale_12,
# pcol=pcol)
## simulation
#data = GeoSim(coordx=coords, corrmodel="Bi_Matern",model=model,param=param)$data
#fixed=list(mean_1=mean_1,mean_2=mean_2, nugget_1=nugget_1,nugget_2=nugget_2,
# smooth_1=smooth_1, smooth_2=smooth_2,smooth_12=smooth_12)
#start=list( sill_1=sill_1,sill_2=sill_2,
# scale_1=scale_1,scale_2=scale_2,scale_12=scale_12, pcol=pcol)
## estimation with maximum likelihood
#fit = GeoFit(data=data,coordx=coords, corrmodel="Bi_Matern",
#likelihood="Full",type="Standard",optimizer="BFGS",
#likelihood="Marginal",type="Pairwise",optimizer="BFGS",neighb=5,
#start=start,fixed=fixed)
###### co-kriging for the fist component ##############
#xx=seq(0,1,0.022)
#loc_to_pred=as.matrix(expand.grid(xx,xx))
#pr1 = GeoKrig(fit,which=1,mse=TRUE,loc=loc_to_pred)
#opar=par(no.readonly = TRUE)
#par(mfrow=c(1,2))
#zlim=c(-2.5,2.5)
#colour = rainbow(100)
#quilt.plot(coords,data[1,] ,col=colour,main = paste(" Fist component"))
#quilt.plot(loc_to_pred,pr1$pred,col=colour,
# main = paste(" Kriging first component"),ylab="")
#par(opar)