Say that I have some observations, $y_1, y_2, ...y_n$ that are described using a generalized additive model (GAM). A Leave-One-Out Cross-Validation (LOOCV) is then performed where each observation $y_i$ is removed and then the GAM is refit to the remaining observations and the omitted observation, $\hat{y}_i$, is predicted. This gives errors $\epsilon_i=(\hat{y}_i - y_i)$ that can be used to calculate useful metrics of my model's predictive strength. For example, I could calculate mean absolute error (MAE):
$\text{MAE} = \frac{1}{n}\sum_{i=1}^n|\epsilon_i|$
So, that's cool, but rather than a point estimate of the error like this I'm curious if I'm able to calculate a confidence interval for the error based on these calculated errors? If so, is this also possible for small data sets (e.g., $n<20$)?
coefs(mod)
for the mean vector and thevcov(mod)
for $\Sigma$ of a multivariate normal, simulate from that. Plug each of those intornorm()
as themu
argument and forsigma
use the residual standard error of residuals from model. With enough simulations per point usequantile()
for CI $\endgroup$rcoefs <- mvrnom()
with the $Xp$ matrix returned byXp <- predict(model, newdata, type = 'lpmatrix')
, wheremodel
is your fitted GAM,newdata
is a data frame of new values you want predictions for - if you don't passnewdata
then you'll get the $Xp$ for the observed data. Then computeXp %*% t(rcoefs[, i])
(i.e. you do a matrix multiplication of the $Xp$ with a row vector of coefficients (hence the transpose,t()
). If you do that for alli
cols of the random draws... $\endgroup$n.draws
from the posterior distribution of the fitted values from the model. (Eachi
yields a new set of fitted values for ith draw from the model posterior.) You can then summarise then.draws
values for each observation to say get the 0.025 and 0.975 probability quantiles which gives you a 95% interval on the fitted/predicted values. If your response isn't Gaussian, then you'll need to simulate from the required distribution using the fitted values as means/parameters for the relevant distribution. Note you probably wantvcov(model, unconditional =TRUE)
. $\endgroup$