Apply and dynamically expand an interval using backtesting

Depending on the aim of generating a forecast, evaluating accurate confidence intervals could be a crucial task. Most classic econometric models, built upon assumptions about distributions of predictions and residuals, have a technique to do that in-built. When moving to machine learning to do time series, equivalent to with XGBoost or recurrent neural networks, it might be more complicated. A preferred technique is conformal intervals — a technique of quantifying uncertainty that makes no assumptions about prediction distributions.

The best implementation of this method is to coach a model and hold out a test set. If this test set is not less than 20 observations (assuming we wish 95% certainty), we will construct an interval by placing a plus/minus value on any point prediction that represents the ninety fifth percentile of the test-set residual absolute values. We then refit the model on the whole series and apply this plus/minus to all point predictions over the unknown horizon. This could be regarded as a naive conformal interval.

Scalecast is a forecasting library in Python that works well if you need to transform a series before applying an optimized machine or deep learning model on the time series, then easily revert the outcomes. While other libraries offer flexible and dynamic intervals for ML models, I’m unsure they’re built to effectively handle data that should be transformed after which reverted, especially with differencing. Please correct me if I’m fallacious.

I made scalecast specifically for this purpose. Transforming and reverting series is ridiculously easy with the package, including options to make use of cross validation to seek out optimal transformation combos. Nonetheless, applying a confidence interval, any interval, at a differenced level becomes complicated to map onto the series’ original level. In case you simply undifference the interval the identical way you’d the purpose predictions, it can more-than-likely be unrealistically wide. A suggestion to avoid that is to make use of models that don’t require stationary data — ARIMA, exponential smoothing, and others. But that’s no fun when you actually need to match ML models and your data will not be stationary.

The answer around this that scalecast uses is the naive conformal interval described above. If a series first, second, or seasonal difference is taken after which reverted, re-calculating test-set residuals and applying percentile functions on them is straightforward. I evaluated the efficacy of this interval using a measure called MSIS in a past post.

pip install --upgrade scalecast

But, this may very well be higher. In time series, it’s intuitive to consider that an interval will expand further out for some extent prediction that’s further away temporally from a baseline truth, once errors accumulate. It is less complicated to predict what I’ll do tomorrow than it’s a month from now. That intuitive concept is built into econometric methods but absent from the naive interval.

We could attempt to reconcile this problem in several ways, certainly one of them being conformalized quantile regression, equivalent to utilized by Neural Prophet. Which will make its technique to scalecast someday. But the tactic I’ll overview here involves using backtesting and applying percentiles depending on residuals from each backtest iteration. Versus employing assumptions, the tactic grounds every thing in some observed, empirical truth — a real relationship between the implemented model and the time series it’s being implemented on.

To do that, we want to separate our data into several training and test sets. Each test set will must be the identical length as our desired forecast horizon and the variety of splits should equal not less than one divided by alpha, where alpha is one minus the specified confidence level. Again, this may come out to twenty iterations for intervals of 95% certainty. Considering we want to iterate through the whole length of our forecast horizon 20 or more times, shorter series may run out of observations before this process finishes. A technique to mitigate this is that if we allow test sets to overlap. So long as test sets are not less than one commentary different from each other and no data leaks from any of the training sets, it must be okay. This will likely bias the interval towards more moderen observations, but the choice could be left open so as to add more room between training sets if the series comprises enough observations to achieve this. The method I explained is known as backtesting, but it might even be regarded as modified time-series cross validation, which is a standard technique to facilitate more accurate conformal intervals. Scalecast makes the means of obtaining this interval easy through pipelines and three utility functions.

Construct Full Model Pipeline

First we construct a pipeline. Assuming we wish differenced data and to forecast with the XGBoost model, the pipeline could be:

transformer = Transformer(['DiffTransform'])
reverter = Reverter(['DiffRevert'],base_transformer=transformer)
def forecaster(f):
f.add_ar_terms(100)
f.add_seasonal_regressors('month')
f.set_estimator('xgboost')
f.manual_forecast()
pipeline = Pipeline(
steps = [
('Transform',transformer),
('Forecast',forecaster),
('Revert',reverter)
]
)

It will be important to notice that this framework can be applied to deep learning models, classic econometric models, RNNs, and even naive models. For any model you need to apply on time series, this may work.

Next, we fit_predict() the pipeline, generating 24 future observations:

f = Forecaster(
y=starts, # an array of observed values
current_dates=starts.index, # an array of dates
future_dates=24, # 24-length forecast horizon
test_length=24, # 24-length test-set for confidence intervals
cis=True, # generate naive intervals for comparison with the final result
)
f = pipeline.fit_predict(f)

Backtest Pipeline and Construct Residual Matrix

Now, we do the backtest. For 95% intervals, this implies not less than 20 train/test splits iteratively moving backward through probably the most recent observations. That is probably the most computationally expensive a part of the method, depending on what number of models we wish to send through our pipeline (we will place more by expanding the forecaster() function), if we wish to optimize each model’s hyperparameters, and if we wish to make use of multivariate processes, it might take some time. On my Macbook, this easy pipeline takes somewhat over a minute to backtest with 20 iterations.

backtest_results = backtest_for_resid_matrix(
f, # a number of positional Forecaster objects can go here
pipeline=pipeline, # each univariate and multivariate pipelines supported
alpha = 0.05, # 0.05 for 95% cis (default) 
bt_n_iter = None, # by default uses the minimum required: 20 for 95% cis, 10 for 90%, etc.
jump_back = 1, # space between training sets, default 1
)

The backtest results from this function could be multipurpose. We will use them to report average errors from our model(s), we will use them to glean insight about many out-of-sample predictions, or we will use them to generate intervals. To generate intervals, we do:

backtest_resid_matrix = get_backtest_resid_matrix(backtest_results)

This creates one matrix for every evaluated model with a row that represents each backtest iteration and a column for every forecast step. The worth in each cell is a prediction error (residual).

Applying a column-wise percentile function, we will generate plus/minus values to seek out the ninety fifth percentile of absolutely the residuals along each forecast step. On average, this must be a bigger value the further out the forecast goes. In our example, the plus/minus is 15 for step 1, 16 for step 4, and 46 for step 24 (the last step). Not all values are progressively larger than the last, but they often are.

Construct Backtested Intervals

We then overwrite the stale naive intervals with the brand new dynamic ones.

overwrite_forecast_intervals(
f, # a number of positional Forecaster objects can go here
backtest_resid_matrix=backtest_resid_matrix,
models=None, # if a couple of models are within the matrix, subset down here
alpha = .05, # 0.05 for 95% cis (default)
)

Voila! We’ve got an assumption-free and dynamic conformal interval built for our time series model.

How a lot better is that this interval than the default? Using MSIS, a measure not many learn about or use, we will rating each obtained interval before and after this process. We can even use the coverage rate from each interval (the percent of actual observations fell throughout the interval). We’ve put aside a separate slice of the info, which overlaps with none of our previously evaluated test sets, for just this purpose. The naive interval looks like:

This ended up being an accurate forecast with a decent interval. It contained 100% of the particular observations and scored an MSIS of 4.03. From my limited usage of MSIS, I feel anything under 5 will likely be pretty good. We apply the dynamic interval and get the next:

This is sweet. We’ve got an expanding interval that’s tighter on average than the default interval. The MSIS rating improved barely to three.92. The bad news: 3 out of the 24 test-set observations fall out of this latest interval’s range for a coverage rating of 87.5%. For a 95% interval, that is probably not ideal.

After all, this is only one example so we should always watch out about drawing too broad conclusions. I do feel confident that the backtested interval will almost all the time expand out, which makes it more intuitive than the default interval. It probably on average is more accurate as well. It just costs more computational power to acquire.

Along with gaining latest intervals, we also obtained backtest information. Over 20 iterations, we observed the next error metrics:

We will feel higher about reporting these than the errors from only one test set.

Thanks for following along! In case you like this content, it could mean loads to me when you gave scalecast a star on GitHub and take a look at the full notebook that accompanies this text. The information used is publicly available through FRED.