Universal Kriging
In KRIG: Spatial Statistic with Kriging

An example presenting an application of Universal Kriging.

Loading the library

library( KRIG )

Defining objective function to approach.

m<-100
XLim<-c( -3, 6 )
x<-seq( XLim[1], XLim[2], length.out = m )
f<-function(x){
  return( 1 + exp(-x^2) + 0.4 * exp(-(x-2)^2 ) - 0.5 * exp( -(x-4)^2 ) ) 
}
z<-sapply( x, f )

Sampling function values.

n<-10
X<-matrix( runif( n, XLim[1], XLim[2] ), n, 1 )

Z<-matrix( sapply( X, f ), n, 1 )

m<-100
YLim<-c(-3.5,6.5)
Y<-matrix( seq( YLim[1], YLim[2], length.out = m ), m, 1 )

Computing constraint matrices for universal Kriging

G<-rbind( t( X ), t( X * X ), t( X * X * X ) )

g<-rbind( t( Y ), t( Y * Y ), t( Y * Y * Y ) )

Kernel definition.

s<-0.1
t<-1

Kern<-function( x, y ) {
  h<-sqrt( sum( ( x - y )^2 ) )
  return( gaussian_kernel( h, s, t ) )
}

Computing Kriging.

K = Kov( X, X, Kern, TRUE )
k = Kov( Y, X, Kern );

KRIG<-Krig( Z = Z,
            K = K, 
            k = k,
            G = G,
            g = g,
            type = "universal", 
            cinv = 'syminv' )

Verification that the constraint is satisfied.

round( max( G %*% KRIG$L - g ), 10 )

Plotting results.

ymin<-min( z, KRIG$Z[,1], Z[,1] )
ymax<-max( z, KRIG$Z[,1], Z[,1] )
plot( x, z, type = 'l', lwd = 2, col='gold', ylim = c( ymin, ymax ), xlim = YLim )
points( Y, KRIG$Z[,1], col='darkgreen', type = 'l', lwd = 2 )
points( X, Z[,1], col = 'dodgerblue3', pch = 16, cex = 1.2 )