# Regarding Hyndman's approach to estimating prediction intervals for forecasts generated by neural networks

I'm currently looking for ways to estimate prediction intervals from an LSTM generated forecast. Several advanced methods are suggested in the literature (e.g. SQF-RNN), but as a first pass, I'm looking for a "quick and dirty" method.

I came across the approach that Hyndman et al. use for the Neural Network forecasting models in the forecast package (NNETAR method), which are described in this post. The approach can be outlined as:

• Assume that your trained neural network is equivalent to $$y_{t+1}=f(y_t,...,y_{t-n})+\epsilon_t$$.
• Generate random samples of $$\epsilon_t$$, and plug those samples into the recursive forecasting equation, so that instead of $$\hat{y}_{t+1}=f(y_t,...,y_{t-n})$$, we have $$y^{*}_{t+1}=f(y_t,...,y_{t-n})+\epsilon^{*}_{t+1}$$.
• Apply this formula recursively all the way up to the forecast horizon $$h$$, each time using a different sample of $$\epsilon^*_{t+k}$$, so that we end up with $$[y^{*}_{t+1},...,y^{*}_{t+h}]$$. This represent a single sample path.
• Repeat this process enough times, until we have a large enough population of sample paths to estimate the forecast distribution.
• Calculate the prediction intervals from the percentiles of the forecast distribution.

In order to generate samples of $$\epsilon_t$$, Hyndman suggests either assuming that $$\epsilon_t$$ follows a normal distribution, or randomly drawing from the historical errors $$[\epsilon_t, \epsilon_{t-1}, \epsilon_{t-2},...]$$ (i.e. the model residuals).

This is where I find the approaches problematic:

• In the first case, the assumption of normality is a pretty strong assumption to make in a lot (most?) of time series cases. In fact, by assuming $$\epsilon_t \sim \mathcal{N}(0,\sigma^2)$$, aren't we effectively creating a self-fulfilling prophecy w/r regards to the subsequent distribution of the sample paths?
• In the second case, the errors $$\epsilon^*_{t+k}$$ are being sampled directly from the residuals, which at first pass makes slightly better sense than the first approach, since it doesn't assume normality (even if it does still erase any time dependencies in the variance). But then I remembered that with neural networks, the out of sample (test) error is almost always significantly larger than the in-sample (train) error - so this approach leads to a serious risk of underestimating the forecast intervals.

Are there any facts that I am missing about the two approaches, such that one or both of them, are actually reasonable approximations, and not just ad-hoc stop-gaps?

On top of these two concerns, I noticed when generating prediction intervals using NNETAR, that the forecast intervals don't get wider as the forecast horizon gets bigger (usually a cardinal sin in the demand forecasting world that I spend most of my time in).

Since the two approaches for sampling $$\epsilon^*_{t+k}$$ are ad-hoc as it is, shouldn't be possible to tweak the above approach in order to have more realistic forecast intervals that widen as the forecast horizon increases? Why is it not possible, following the formula for multi-step prediction intervals for naive forecasts, to use something like:

$$\epsilon^{*}_{t+k} = \epsilon^{*}_k \sqrt{k}$$, with $$\epsilon^{*}_k$$ the value obtained form sampling the distribution of $$\epsilon_t$$ (using either of the sampling approaches) and $$\epsilon^{*}_{t+k}$$ a weighted estimate of $$\epsilon_{t+k}$$ that takes into account how far into the future $$k$$ is?

• Your approach in the last paragraph would be valid if the residuals have no autocorrelation, but in general they will - hence Hyndman's approach is a form of bootstrap which attempts to maintain the autocorrelation of the residuals. There are approximations such as yours which can be used if you "know" the ACF of the residuals - but of course normally you don't. See github.com/bashtage/arch/issues/323 for a discussion I had recently on implementing a bootstrap which includes the variance of the estimator and the residuals. – David Waterworth Feb 11 at 5:39
• @DavidWaterworth I'm not following: since Bootstraping is sampling with replacement, if anything it would muddy up any auto-correlations in the variance as opposed to maintaining them? – Skander H. Feb 11 at 6:08
• There are bootstrap methods specifically for timeseries i.e. arch.readthedocs.io/en/latest/bootstrap/… – David Waterworth Feb 11 at 21:19
• Also I'd add that I believe Hyndmanns method (and bootstrap in general) is intended to correct for heteroskedastic iid errors. Using a timeseries bootstrap extends to autocorrelated heteroskedastic errors. If neither are present the naive formula is valid. – David Waterworth Feb 12 at 22:36