estimateCV.STmodel: Cross-Validated Estimation and Prediction

Description

Functions that perform cross-validated parameter estimation and prediction for the spatio-temporal model.

Usage

## S3 method for class 'STmodel':
estimateCV(object, x, Ind.cv,
    control = list(trace = 3), ...)
  ## S3 method for class 'STmodel':
predictCV(object, x, Ind.cv, ...,
    silent = TRUE)

Arguments

object

STmodel object for which to perform cross-validation.

Either a vector or matrix of starting point(s) for the optimisation, see estimate.STmodel; or a matrix with parameters, the i:th row being used for prediction of the i:th cross-validation set

Ind.cv

Ind.cv defines the cross-validation scheme. Either a (number or observations) - by - (groups) logical matrix or an integer valued vector with length equal to (number or observations). See further

control

A list of control parameters for the optimisation. See optim for details; setting trace=0 eliminates all ouput.

...

All additional parameters for estimate.STmodel or predict.STmodel. For predict.STmodel<

silent

Show status after each iteration?

Value

Either a estCVSTmodel object with elements:
statusData.frame with convergence information and best function value for each cross-validation group.
Ind.cvThe cross-validation grouping.
x.fixedFixed parameters in the estimation, see estimate.STmodel.
par.all, par.covMatrices with estimated parameters for each cross-validation group.
par.all.sd, par.cov.sdStandard deviations computed from the Hessian/information matrix for set of estimated parameters.
res.allEstimation results for each cross-validation group, contains the res.best and summary fields from the estimate.STmodel results.
Or a predCVSTmodel object with elements:
optsCopy of the opts field in the output from predict.STmodel.
Ind.cvThe cross-validation grouping.
pred.obsA data.frame with a copy of observations from object$obs, predictions (for different model components), variances, and residuals. Variance field will be missing if pred.var=FALSE.
pred.allA list with time-by-location data.frames containing predictions and variances for all space-time locations as well as predictions and variances for the beta-fields. Unobserved points are NA for the option only.obs=TRUE.

Details

For predictCV.STmodel the parameters used to compute predictions for the left out observations can be either a single vector or a matrix. For a single vector the same parameter values will be used for all cross-validation predictions; for a matrix the parameters in x[,i] will be used for the predictions of the i:th cross-validation set (i.e. for Ind.cv[,i]). Suitable matrices are provided in the output from estimateCV.STmodel.

The cross-validation groups are given by Ind.cv. Ind.cv should be either a (number of observations) - by - (groups) logical matrix or an integer valued vector with length equal to (number of observations). If a matrix then each column defines a cross-validation set with the TRUE values marking the observations to be left out. If a vector then 1:s denote observations to be dropped in the first cross-validation set, 2:s observations to be dropped in the second set, etc. Observations marked by values <=0< code=""> are never dropped. See createCV for details.

Examples

Run this code

##load data
data(mesa.model)

################
## estimateCV ##
################
##create the CV structure defining 10 different CV-groups
Ind.cv <- createCV(mesa.model, groups=10, min.dist=.1)

##create vector of initial values
dim <- loglikeSTdim(mesa.model)
x.init <- cbind(rep(2,dim$nparam.cov),c(rep(c(1,-3),dim$m+1),-3))

##estimate different parameters for each CV-group
  est.cv.mesa <- estimateCV(mesa.model, x.init, Ind.cv)
##lets load precomputed results instead
data(CV.mesa.model)

##examine the estimation results
print( est.cv.mesa )
##estimated parameters for each CV-group
est.cv.mesa$par.cov

##boxplot of the different estimates from the CV
par(mfrow=c(1,1), mar=c(12,2.5,2,.5), las=2)
boxplot( est.cv.mesa, plot.type="cov")

###############
## predictCV ##
###############
##Do cross-validated predictions using the just estimated parameters
  pred.cv.mesa <- predictCV(mesa.model, est.cv.mesa, est.cv.mesa$Ind.cv)
##lets load precomputed results instead
data(CV.mesa.model)

##examine the prediction results
names( pred.cv.mesa )
print( pred.cv.mesa )
##compute CV-statistics
pred.cv.stats <- summary( pred.cv.mesa, LTA=TRUE)
print( pred.cv.stats )
 
##Plot observations with CV-predictions and prediction intervals
##TODO TODO TODO
#par(mfcol=c(4,1),mar=c(2.5,2.5,2,.5))
#plotCV(mesa.data.res$pred.cv,  1, mesa.data.model)
#plotCV(mesa.data.res$pred.cv, 10, mesa.data.model)
#plotCV(mesa.data.res$pred.cv, 17, mesa.data.model)
#plotCV(mesa.data.res$pred.cv, 22, mesa.data.model)

##Residual QQ-plots
#par(mfrow=c(1,2),mar=c(3,2,1,1),pty="s")
##Raw Residual QQ-plot
#CVresiduals.qqnorm(pred.cv.stats$res)
##Normalized Residual QQ-plot
#CVresiduals.qqnorm(pred.cv.stats$res.norm, norm=TRUE)

##A faster option is to only consider the observations and not to compute
  ##variances
  pred.cv.fast <- predictCV(mesa.model, est.cv.mesa, est.cv.mesa$Ind.cv,
                            only.obs=TRUE, pred.var=FALSE)
  print( pred.cv.fast )
  summary( pred.cv.fast )

Run the code above in your browser using DataLab