2
$\begingroup$

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$)?

Improve this question
$\endgroup$
9
  • $\begingroup$ Upvoted the answer, but why would you fit a GAM to 20 data points anyway? In most cases it's hard to avoid parametric structure with such a small dataset. $\endgroup$
    – generic_user
    Oct 5, 2017 at 14:22
  • $\begingroup$ Are you interested in other approaches? You could simulate from the posterior distribution of the GAM and the conditional distribution of the response. If the model is a Gaussian GAM, simulate a bunch of draws (for each point you want to predict at) from the posterior of the model (if mgcv use the coefs(mod) for the mean vector and the vcov(mod) for $\Sigma$ of a multivariate normal, simulate from that. Plug each of those into rnorm() as the mu argument and for sigma use the residual standard error of residuals from model. With enough simulations per point use quantile() for CI $\endgroup$
    – Gavin Simpson
    Oct 5, 2017 at 15:42
  • 1
    $\begingroup$ You need to couple that matrix of draws from rcoefs <- mvrnom() with the $Xp$ matrix returned by Xp <- predict(model, newdata, type = 'lpmatrix'), where model is your fitted GAM, newdata is a data frame of new values you want predictions for - if you don't pass newdata then you'll get the $Xp$ for the observed data. Then compute Xp %*% 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 all i cols of the random draws... $\endgroup$
    – Gavin Simpson
    Oct 10, 2017 at 17:25
  • 1
    $\begingroup$ ...you have n.draws from the posterior distribution of the fitted values from the model. (Each i yields a new set of fitted values for ith draw from the model posterior.) You can then summarise the n.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 want vcov(model, unconditional =TRUE). $\endgroup$
    – Gavin Simpson
    Oct 10, 2017 at 17:29
  • 1
    $\begingroup$ I'll try to find some time later today to cook up an example and post it as an Answer. $\endgroup$
    – Gavin Simpson
    Oct 10, 2017 at 17:30

1 Answer 1

Reset to default
2
$\begingroup$

Let's assume for the sake of the argument that your sample points $\epsilon_i$ are i.i.d. You cannot assume they are normally distributed unless you have a good reason to. $n=20$ would be too few to assume that their sample mean is normally distributed. This combination leaves you already in a tough spot to compute your sample mean and standard error which you need for a confidence interval.

Furthermore, your sample points are actually not i.i.d., the cross validation introduces pseudo replication. Since you did LOOCV and thus without repetition, the pseudo replication comes only from the fact that each sample point is predicted through a model that was trained on almost the same training set. (Your "test-sets" which are individual records do not overlap.) There are corrected resampled tests which correct for cross validation induced pseudo replication. They make your standard error larger and thus your confidence interval wider.

The usual standard error would be:

$$\sigma/\sqrt{n}$$

The corrected resampled one is:

$$\sigma\times\sqrt{\frac{1}{n}+\frac{1}{1-n}}$$

Your confidence intervals most likely will be so wide as to be useless in this situation.

If you had slightly bigger data-sets $n>30$, you could do a confidence interval since you wouldn't need to assume normality anymore. If you have one like this, just try and see how large your confidence interval is with the resampled correction applied. Otherwise stay with the type of measures you mentioned.

Improve this answer
$\endgroup$
3
  • $\begingroup$ So, MAE is probably a more useful metric of model error? $\endgroup$
    – Dan
    Oct 5, 2017 at 13:58
  • $\begingroup$ If you had slightly bigger data-sets $n>30$, you could do a confidence interval since you wouldn't need to assume normality anymore. If you have one like this, just try and see how large your confidence interval is with the resampled correction applied. Otherwise stay with the type of measures you mentioned. $\endgroup$
    – David Ernst
    Oct 5, 2017 at 14:18
  • $\begingroup$ In the corrected standard error, the value under the square root is always negative, right? If so, isn't this a problem? $\endgroup$
    – MD004
    Jan 4, 2019 at 20:35

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.