Predicition interval for LOESS smoothed data Is the concept of a prediction interval for non-parameteric fits such as LOESS meaningful? I can't see why it wouldn't be - and the fact that it is implemented in R also implies it is.
https://mclust-org.github.io/msir/reference/loess.sd.html
But I can't find any reference on how to derive it (including in the book referenced in the link above) however, which makes me think I am missing something - unlike for parameteric modelling where the prediction interval is heavily documented.
I would like to derive prediction intervals for a various fitting techniques, including LOESS, which is why I'd like to understand it better- guidance welcome
 A: For reference, the specific function you've linked basically calculates the non-parametric variance function of the Loess fit.
The intervals constructed are referred to as variance smooths by Cook & Weisberg (1999, pp. 52 - 53) for the Loess, and are calculated as:
$$\hat{E}(y | x) \pm \sqrt{\widehat{Var}(y | x)}$$
Where $\hat{E}(y | x)$ is the estimated mean function and $\sqrt{\widehat{Var}(y | x)}$ is the estimated standard deviation function, calculated by a procedure as follows (outlined in Section 2.5.2 of the same text, pp. 36 - 37):

*

*Calculate $E(y | x)$ using Loess with a given smoothing parameter, $\alpha$


*Create a new variable representing the squared deviation of each observation from the mean function
$$\hat{e}^2 = (y - \hat{E}(y | x))^2$$


*Calculate $\hat{E}(\hat{e}^2 | x)$ using Loess with the same $\alpha$ - this is our estimate of $\widehat{Var}(y | x)$ and conceptually is like smoothing a scatterplot of the squared deviations on $x$


*Take the square root of the variance function to calculate the standard deviation function
If you dig into the loess.sd function in msir, you'll find this is basically what it does (also allowing you to scale the standard deviation function using the nsigma argument).
But I'll work through an example manually using R to demonstrate. I'm using the cholesterol data from the bootstrap package which evaluated the effects of cholostyramine compliance on blood cholesterol levels.
# Load the Data
cholesterol <- bootstrap::cholost
colnames(cholesterol) <- c("compliance", "improvement")
with(cholesterol, plot(compliance, improvement))

# 1. Fit loess to find mean function E(y | x)
alpha <- 0.6
m1 <- loess(improvement ~ compliance, 
            data = cholesterol,
            span = alpha)

# Generate and plot predictions from mean function
new_dat <- data.frame(compliance = seq(0, 100))
preds <- predict(m1, new_dat)
lines(seq(0, 100), preds)

# 2. Find squared deviations
cholesterol$e <- residuals(m1)^2

# 3. Fit loess to find variance function E(e^2 | x)
m2 <- loess(e ~ compliance, 
            data = cholesterol,
            span = alpha)

# 4. Estimate Std. Deviation function across values of x
sd <- sqrt(predict(m2, new_dat))

# Find upper and lower limits and plot
ll <- preds - sd
ul <- preds + sd
lines(new_dat$compliance, ll, lty = 2)
lines(new_dat$compliance, ul, lty = 2)

Which graphically looks like this:


References
Cook, R. D., & Weisberg, S. (1999). Applied regression including computing and graphics. John Wiley & Sons.
