fields (version 16.3)

KrigFindLambda: Finds profile likelihood and GCV estimates of smoothing parameters for splines and Kriging.

Description

This is a secondary function that will use the computed Krig object and find various estimates of the smoothing parameter lambda. These are several different flavors of cross-validation, a moment matching strategy and the profile likelihood. This function can also be used independently with different data sets (the y's) if the covariates ( the x's) are the same and thus reduce the computation.

Usage

KrigFindLambda(
out, lambda.grid = NA, cost = 1, nstep.cv = 200, rmse
                 = NA, verbose = FALSE, tol = 1e-05, offset = 0, y =
                 NULL, give.warnings = TRUE)
gcv.sreg (
out, lambda.grid = NA, cost = 1, nstep.cv = 80, rmse =
                 NA, offset = 0, trmin = NA, trmax = NA, verbose =
                 FALSE, tol = 1e-05, give.warnings = TRUE)

Value

A list giving a summary of estimates and diagonostic details with the following components:

gcv.grid

A matrix describing results of the coarse search rows are values of lambda and the columns are lambda= value of smoothing parameter, trA=effective degrees of freedom, GCV=Usual GCV criterion, GCV.one=GCV criterion leave-one-out, GCV.model= GCV based on average response in the case of replicates, tauHat= Implied estimate of tau , -Log Profile= negative log of profiel likelihood for the lambda.

lambda.est

Summary table of all estimates Rows index different types of estimates: GCV, GCV.model, GCV.one, RMSE, pure error, -Log Profile and the columns are the estimated values for lambda, trA, GCV, tauHat.

Arguments

out

A Krig or sreg object.

lambda.grid

Grid of lambdas for coarse search. The default is equally spaced on effective degree of freedom scale.

cost

Cost used in GCV denominator

nstep.cv

Number of grid points in coarse search.

rmse

Target root mean squared error to match with the estimate of tau**2

verbose

If true prints intermediate results.

tol

Tolerance in delcaring convergence of golden section search or bisection search.

offset

Additional degrees of freedom to be added into the GCV denominator.

y

A new data vector to be used in place of the one associated with the Krig object (obj)

give.warnings

If FALSE will suppress warnings about grid search being out of range for various estimates based on GCV and REML.

trmin

Minimum value of lambda for grid search specified in terms of effective degrees of freedom.

trmax

Maximum value for grid search.

Author

Doug Nychka

Details

This function finds several estimates of the smoothing parameter using first a coarse grid search followed by a refinement using a minimization ( in the case of GCV or maximum likelihood) or bisection in the case of mathcing the rmse. Details of the estimators can be found in the help file for the Krig function.

The Krig object passed to this function has some matrix decompostions that facilitate rapid computation of the GCV and ML functions and do not depend on the independent variable. This makes it possible to compute the Krig object once and to reuse the decompostions for multiple data sets. (But keep in mind if the x values change then the object must be recalculated.) The example below show show this can be used for a simulation study on the variability for estimating the smoothing parameter.

See Also

Krig, Tps, predict.Krig

Examples

# 
Tps( ChicagoO3$x, ChicagoO3$y)-> obj # default is to find lambda by GCV
summary( obj)

KrigFindLambda( obj)-> out
print( out$lambda.est) # results agree with Tps summary

sreg( rat.diet$t, rat.diet$trt)-> out
gcv.sreg( out, tol=1e-10) # higher tolerance search for minimum 
if (FALSE) {
# a simulation example
x<- seq( 0,1,,150)
f<-  x**2*( 1-x)
f<- f/sqrt( var( f))

set.seed(123) # let's all use the same seed
tau<- .1
y<- f + rnorm( 150)*tau

Tps( x,y)-> obj # create Krig object

hold<- hold2<- matrix( NA, ncol=6, nrow=200)

for( k in 1:200){
# look at GCV estimates of lambda
# new data simulated
   y<- f + rnorm(150)*tau 
# save GCV estimates
lambdaTable<- KrigFindLambda(obj,  y=y, give.warnings=FALSE)$lambda.est
hold[k,]<-  lambdaTable[1,]
hold2[k,]<-  lambdaTable[6,]
}
matplot( cbind(hold[,2], hold2[,2]),cbind( hold[,4],hold2[,4]),
 xlab="estimated eff. df", ylab="tau hat", pch=16, col=c("orange3", "green2"), type="p")
yline( tau, col="grey", lwd=2)

}