#
Tps( ozone$x, ozone$y)-> obj # default is to find lambda by GCV
summary( obj)
gcv.Krig( obj)-> out
print( out$lambda.est) # results agree with Tps summary
# a simulation example
x<- seq( 0,1,,200)
f<- x**2*( 1-x)
f<- f/sqrt( var( f))
Tps( x,f)-> obj # create Krig object
set.seed(123) # let's all use the same seed
sigma<- .1
hold<- matrix( NA, ncol=4, nrow=100)
for( k in 1 :100){
# look at GCV estimates of lambda
# new data simulated
y<- f + rnorm(200)*sigma
# save GCV estimates
hold[k,]<- gcv.Krig(obj, y=y, give.warnings=FALSE)$lambda.est[1,]
}
plot( hold[,2], hold[,4], xlab="estimated eff. df", ylab="sigma hat")
yline( sigma, col=2)
# note some occaisional flaky behaviour with GCV ( eff df > 20!)
Run the code above in your browser using DataLab