# Load ozone data set
data(ozone2)
x<-ozone2$lon.lat
y<- ozone2$y[16,]
# Find location that are not 'NA'.
# (LKrig is not set up to handle missing observations.)
good <- !is.na( y)
x<- x[good,]
y<- y[good]
LKinfo<- LKrigSetup( x,NC=20,nlevel=1, alpha=1, lambda= .3 , a.wght=5)
# BTW lambda is close to MLE
# Simulating this LKrig process
# simulate 4 realizations of process and plot them
# (these have unit marginal variance)
xg<- make.surface.grid(list( x=seq( -87,-83,,40), y=seq(36.5, 44.5,,40)))
out<- LKrig.sim(xg, LKinfo,M=4)
if (FALSE) {
set.panel(2,2)
for( k in 1:4){
image.plot( as.surface( xg, out[,k]), axes=FALSE) }
}
obj<- LKrig(x,y,LKinfo=LKinfo)
O3.cond.sim<- LKrig.sim.conditional( obj, M=3,nx=40,ny=40)
if (FALSE) {
set.panel( 2,2)
zr<- range( c( O3.cond.sim$draw, O3.cond.sim$ghat), na.rm=TRUE)
coltab<- tim.colors()
image.plot( as.surface( O3.cond.sim$x.grid, O3.cond.sim$ghat), zlim=zr)
title("Conditional mean")
US( add=TRUE)
for( k in 1:3){
image( as.surface( O3.cond.sim$x.grid, O3.cond.sim$g.draw[,k]),
zlim=zr, col=coltab)
points( obj$x, cex=.5)
US( add=TRUE)
}
set.panel()
}
Run the code above in your browser using DataLab