# \donttest{
if(require(fda)){
## load the data
data("CanadianWeather", package = "fda")
## use data on a daily basis
canada <- with(CanadianWeather,
list(temp = t(dailyAv[ , , "Temperature.C"]),
l10precip = t(dailyAv[ , , "log10precip"]),
l10precip_mean = log(colMeans(dailyAv[ , , "Precipitation.mm"]), base = 10),
lat = coordinates[ , "N.latitude"],
lon = coordinates[ , "W.longitude"],
region = factor(region),
place = factor(place),
day = 1:365, ## corresponds to t: evaluation points of the fun. response
day_s = 1:365)) ## corresponds to s: evaluation points of the fun. covariate
## center temperature curves per day
canada$tempRaw <- canada$temp
canada$temp <- scale(canada$temp, scale = FALSE)
rownames(canada$temp) <- NULL ## delete row-names
## fit the model
mod <- FDboost(l10precip ~ 1 + bolsc(region, df = 4) +
bsignal(temp, s = day_s, cyclic = TRUE, boundary.knots = c(0.5, 365.5)),
timeformula = ~ bbs(day, cyclic = TRUE, boundary.knots = c(0.5, 365.5)),
data = canada)
mod <- mod[75]
#### create folds for 3-fold bootstrap: one weight for each curve
set.seed(124)
folds_bs <- cv(weights = rep(1, mod$ydim[1]), type = "bootstrap", B = 3)
## compute out-of-bag risk on the 3 folds for 1 to 75 boosting iterations
cvr <- applyFolds(mod, folds = folds_bs, grid = 1:75)
## compute out-of-bag risk and coefficient estimates on folds
cvr2 <- validateFDboost(mod, folds = folds_bs, grid = 1:75)
## weights per observation point
folds_bs_long <- folds_bs[rep(1:nrow(folds_bs), times = mod$ydim[2]), ]
attr(folds_bs_long, "type") <- "3-fold bootstrap"
## compute out-of-bag risk on the 3 folds for 1 to 75 boosting iterations
cvr3 <- cvrisk(mod, folds = folds_bs_long, grid = 1:75)
## plot the out-of-bag risk
oldpar <- par(mfrow = c(1,3))
plot(cvr); legend("topright", lty=2, paste(mstop(cvr)))
plot(cvr2)
plot(cvr3); legend("topright", lty=2, paste(mstop(cvr3)))
## plot the estimated coefficients per fold
## more meaningful for higher number of folds, e.g., B = 100
par(mfrow = c(2,2))
plotPredCoef(cvr2, terms = FALSE, which = 1)
plotPredCoef(cvr2, terms = FALSE, which = 3)
## compute out-of-bag risk and predictions for leaving-one-curve-out cross-validation
cvr_jackknife <- validateFDboost(mod, folds = cvLong(unique(mod$id),
type = "curves"), grid = 1:75)
plot(cvr_jackknife)
## plot oob predictions per fold for 3rd effect
plotPredCoef(cvr_jackknife, which = 3)
## plot coefficients per fold for 2nd effect
plotPredCoef(cvr_jackknife, which = 2, terms = FALSE)
par(oldpar)
}
# }
Run the code above in your browser using DataLab