Backtest Orbit Model

In this section, we will cover:

  • How to create a TimeSeriesSplitter

  • How to create a BackTester and retrieve the backtesting results

  • How to leverage the backtesting to tune the hyper-paramters for orbit models

[1]:
%matplotlib inline

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt

import orbit
from orbit.models.lgt import LGTMAP, LGTAggregated
from orbit.models.dlt import DLTMAP
from orbit.diagnostics.backtest import BackTester, TimeSeriesSplitter
from orbit.utils.dataset import load_iclaims
from orbit.utils.plot import get_orbit_style
plt.style.use(get_orbit_style())

import warnings
warnings.filterwarnings('ignore')
[2]:
print(orbit.__version__)
1.0.17
[3]:
# load log-transformed data
data = load_iclaims()
[4]:
data.shape
[4]:
(443, 7)

The way to gauge the performance of a time-series model is through re-training models with different historic periods and check their forecast within certain steps. This is similar to a time-based style cross-validation. More often, we called it backtest in time-series modeling.

The purpose of this notebook is to illustrate how to “backtest” a single model using BackTester

BackTester will compose a TimeSeriesSplitter within it, but TimeSeriesSplitter is useful as a standalone, in case there are other tasks to perform that requires splitting but not backtesting. TimeSeriesSplitter implemented each ‘slices’ as genertor, i.e it can be used in a for loop. You can also retrieve the composed TimeSeriesSplitter object from BackTester to utilize the additional methods in TimeSeriesSplitter

Currently, there are two schemes supported for the back-testing engine: expanding window and rolling window.

  • expanding window: for each back-testing model training, the train start date is fixed, while the train end date is extended forward.

  • rolling window: for each back-testing model training, the training window length is fixed but the window is moving forward.

Create a TimeSeriesSplitter

There two main way to splitting a timeseries: expanding and rolling. Expanding window has a fixed starting point, and the window length grows as we move forward in timeseries. It is useful when we want to incoporate all historical information. On the other hand, rolling window has a fixed window length, and the starting point of the window moves forward as we move forward in timeseries. Now, we will illustrate how to use TimeSeriesSplitter to split the claims timeseries.

Expanding window

[5]:
# configs
min_train_len = 380 # minimal length of window length
forecast_len = 20 # length forecast window
incremental_len = 20 # step length for moving forward
[6]:
ex_splitter = TimeSeriesSplitter(df=data,
                                 min_train_len=min_train_len,
                                 incremental_len=incremental_len,
                                 forecast_len=forecast_len,
                                 window_type='expanding',
                                 date_col='week')
[7]:
print(ex_splitter)

------------ Fold: (1 / 3)------------
Train start index: 0 Train end index: 379
Test start index: 380 Test end index: 399
Train start date: 2010-01-03 00:00:00 Train end date: 2017-04-09 00:00:00
Test start date: 2017-04-16 00:00:00 Test end date: 2017-08-27 00:00:00

------------ Fold: (2 / 3)------------
Train start index: 0 Train end index: 399
Test start index: 400 Test end index: 419
Train start date: 2010-01-03 00:00:00 Train end date: 2017-08-27 00:00:00
Test start date: 2017-09-03 00:00:00 Test end date: 2018-01-14 00:00:00

------------ Fold: (3 / 3)------------
Train start index: 0 Train end index: 419
Test start index: 420 Test end index: 439
Train start date: 2010-01-03 00:00:00 Train end date: 2018-01-14 00:00:00
Test start date: 2018-01-21 00:00:00 Test end date: 2018-06-03 00:00:00

We can visualize the splits, green is training window and yellow it the forecasting windown. The starting point is always 0 for three splits but window length increases from 380 to 420.

[8]:
_ = ex_splitter.plot()
../_images/tutorials_backtest_13_0.png

Rolling window

[9]:
# configs
min_train_len = 380 # in case of rolling window, this specify the length of window length
forecast_len = 20 # length forecast window
incremental_len = 20 # step length for moving forward
[10]:
roll_splitter = TimeSeriesSplitter(data,
                                   min_train_len=min_train_len,
                                   incremental_len=incremental_len,
                                   forecast_len=forecast_len,
                                   window_type='rolling', date_col='week')

We can visualize the splits, green is training window and yellow it the forecasting windown. The window length is always 380, while the starting point moves forward 20 weeks each steps.

[11]:
_ = roll_splitter.plot()
../_images/tutorials_backtest_18_0.png

Specifying number of splits

User can also define number of splits using n_splits instead of specifying minimum training length. That way, minimum training length will be automatically calculated.

[12]:
ex_splitter2 = TimeSeriesSplitter(data,
                                  min_train_len=min_train_len,
                                  incremental_len=incremental_len,
                                  forecast_len=forecast_len,
                                  n_splits=5,
                                  window_type='expanding', date_col='week')
[13]:
_ = ex_splitter2.plot()
../_images/tutorials_backtest_22_0.png

TimeSeriesSplitter as generator

TimeSeriesSplitter is implemented as a genetor, therefore we can call split() to loop through it. It comes handy even for tasks other than backtest.

[14]:
 for train_df, test_df, scheme, key in roll_splitter.split():
        print('Initial Claim slice {} rolling mean:{}'.format(key, train_df['claims'].mean()))
Initial Claim slice 0 rolling mean:12.712341458856958
Initial Claim slice 1 rolling mean:12.67144746151623
Initial Claim slice 2 rolling mean:12.64657747604057

Create a BackTester

Now, we are ready to do backtest, first let’s initialize a DLT model and a BackTester. You pass in TimeSeriesSplitter parameters to BackTester.

[15]:
# instantiate a model
dlt = DLTMAP(
    date_col='week',
    response_col='claims',
    regressor_col=['trend.unemploy', 'trend.filling', 'trend.job'],
    seasonality=52,
)
[16]:
# configs
min_train_len = 100
forecast_len = 20
incremental_len = 100
window_type = 'expanding'

bt = BackTester(
    model=dlt,
    df=data,
    min_train_len=min_train_len,
    incremental_len=incremental_len,
    forecast_len=forecast_len,
    window_type=window_type,
)

Backtest fit and predict

The most expensive portion of backtesting is fitting the model iteratively. Thus, we separate the api calls for fit_predict and score to avoid redundant computation for multiple metrics or scoring methods

[17]:
bt.fit_predict()

Once fit_predict() is called, the fitted models and predictions can be easily retrieved from BackTester. Here the data is grouped by the date, split_key, and whether or not that observation is part of the training or test data

[18]:
predicted_df = bt.get_predicted_df()
predicted_df.head()
[18]:
date split_key training_data actuals prediction
0 2010-01-03 0 True 13.386595 13.386595
1 2010-01-10 0 True 13.624218 13.646964
2 2010-01-17 0 True 13.398741 13.371691
3 2010-01-24 0 True 13.137549 13.151830
4 2010-01-31 0 True 13.196760 13.187538

Users might find this useful for any custom computations that may need to be performed on the set of predicted data. Note that the columns are renamed to generic and consistent names.

Sometimes, it might be useful to match the data back to the original dataset for ad-hoc diagnostics. This can easily be done by merging back to the orignal dataset

[19]:
predicted_df.merge(data, left_on='date', right_on='week')
[19]:
date split_key training_data actuals prediction week claims trend.unemploy trend.filling trend.job sp500 vix
0 2010-01-03 0 True 13.386595 13.386595 2010-01-03 13.386595 0.219882 -0.318452 0.117500 -0.417633 0.122654
1 2010-01-03 1 True 13.386595 13.386595 2010-01-03 13.386595 0.219882 -0.318452 0.117500 -0.417633 0.122654
2 2010-01-03 2 True 13.386595 13.386595 2010-01-03 13.386595 0.219882 -0.318452 0.117500 -0.417633 0.122654
3 2010-01-03 3 True 13.386595 13.386595 2010-01-03 13.386595 0.219882 -0.318452 0.117500 -0.417633 0.122654
4 2010-01-10 0 True 13.624218 13.646964 2010-01-10 13.624218 0.219882 -0.194838 0.168794 -0.425480 0.110445
... ... ... ... ... ... ... ... ... ... ... ... ...
1075 2017-12-17 3 False 12.568616 12.566123 2017-12-17 12.568616 0.298663 0.248654 -0.216869 0.434042 -0.482380
1076 2017-12-24 3 False 12.691451 12.675590 2017-12-24 12.691451 0.328516 0.233616 -0.358839 0.430410 -0.373389
1077 2017-12-31 3 False 12.769532 12.782966 2017-12-31 12.769532 0.503457 0.069313 -0.092571 0.456087 -0.553539
1078 2018-01-07 3 False 12.908227 12.799673 2018-01-07 12.908227 0.527849 0.051295 0.029532 0.471673 -0.456456
1079 2018-01-14 3 False 12.777193 12.536322 2018-01-14 12.777193 0.465717 0.032946 0.006275 0.480271 -0.352770

1080 rows × 12 columns

Backtest Scoring

The main purpose of BackTester are the evaluation metrics. Some of the most widely used metrics are implemented and built into the BackTester API.

The default metric list is smape, wmape, mape, mse, mae, rmsse.

[20]:
bt.score()
[20]:
metric_name metric_values is_training_metric
0 smape 0.006617 False
1 wmape 0.006613 False
2 mape 0.006599 False
3 mse 0.012875 False
4 mae 0.084002 False
5 rmsse 0.809899 False

It is possible to filter for only specific metrics of interest, or even implement your own callable and pass into the score() method. For example, see this function that uses last observed value as a predictor and computes the mse. Or naive_error which computes the error as the delta between predicted values and the training period mean.

Note these are not really useful error metrics, just showing some examples of callables you can use ;)

[21]:
def mse_naive(test_actual):
    actual = test_actual[1:]
    predicted = test_actual[:-1]
    return np.mean(np.square(actual - predicted))

def naive_error(train_actual, test_predicted):
    train_mean = np.mean(train_actual)
    return np.mean(np.abs(test_predicted - train_mean))
[22]:
bt.score(metrics=[mse_naive, naive_error])
[22]:
metric_name metric_values is_training_metric
0 mse_naive 0.019628 False
1 naive_error 0.229262 False

It doesn’t take additional time to refit and predict the model, since the results are stored when fit_predict() is called. Check docstrings for function criteria that is required for it to be supported with this api.

In some cases, we may want to evaluate our metrics on both train and test data. To do this you can call score again with the following indicator

[23]:
bt.score(include_training_metrics=True)
[23]:
metric_name metric_values is_training_metric
0 smape 0.006617 False
1 wmape 0.006613 False
2 mape 0.006599 False
3 mse 0.012875 False
4 mae 0.084002 False
5 rmsse 0.809899 False
6 smape 0.002739 True
7 wmape 0.002742 True
8 mape 0.002738 True
9 mse 0.003117 True
10 mae 0.035042 True

Backtest Get Models

In cases where BackTester doesn’t cut it or for more custom use-cases, there’s an interface to export the TimeSeriesSplitter and predicted data, as shown earlier. It’s also possible to get each of the fitted models for deeper diving

[24]:
fitted_models = bt.get_fitted_models()
[25]:
model_1 = fitted_models[0]
model_1.get_regression_coefs()
[25]:
regressor regressor_sign coefficient
0 trend.unemploy Regular -0.045665
1 trend.filling Regular -0.109453
2 trend.job Regular 0.027940

BackTester composes a TimeSeriesSplitter within it, but TimeSeriesSplitter can also be created on its own as a standalone object. See section below on TimeSeriesSplitter for more details on how to use the splitter.

All of the additional TimeSeriesSplitter args can also be passed into BackTester on instantiation

[26]:
ts_splitter = bt.get_splitter()
_ = ts_splitter.plot()
../_images/tutorials_backtest_49_0.png

Hyperparameter Tunning

After seeing the results fromt the backtest, users may wish to fine tune the hyperparmeters. Orbit also provide a grid_search_orbit utilities for parameter searching. It uses Backtester under the hood so users can compare backtest metrics for different paramters combination.

[27]:
from orbit.diagnostics.backtest import grid_search_orbit
# defining the search space for level smoothing paramter and seasonality smooth paramter
param_grid = {
    'level_sm_input': [0.3, 0.5, 0.8],
    'seasonality_sm_input': [0.3, 0.5, 0.8],
}
[28]:
# configs
min_train_len = 380 # in case of rolling window, this specify the length of window length
forecast_len = 20 # length forecast window
incremental_len = 20 # step length for moving forward
best_params, tuned_df = grid_search_orbit(param_grid,
                                        model=dlt,
                                        df=data,
                                        min_train_len=min_train_len,
                                        incremental_len=incremental_len,
                                        forecast_len=forecast_len,
                                        metrics=None, criteria=None, verbose=True)
  0%|          | 0/9 [00:00<?, ?it/s]
tuning hyper-params {'level_sm_input': 0.3, 'seasonality_sm_input': 0.3}
 11%|█         | 1/9 [00:04<00:35,  4.48s/it]
tuning metric:0.0049017
tuning hyper-params {'level_sm_input': 0.3, 'seasonality_sm_input': 0.5}
 22%|██▏       | 2/9 [00:08<00:31,  4.46s/it]
tuning metric:0.0040405
tuning hyper-params {'level_sm_input': 0.3, 'seasonality_sm_input': 0.8}
 33%|███▎      | 3/9 [00:13<00:27,  4.55s/it]
tuning metric:0.0036082
tuning hyper-params {'level_sm_input': 0.5, 'seasonality_sm_input': 0.3}
 44%|████▍     | 4/9 [00:17<00:22,  4.46s/it]
tuning metric:0.0079585
tuning hyper-params {'level_sm_input': 0.5, 'seasonality_sm_input': 0.5}
 56%|█████▌    | 5/9 [00:22<00:17,  4.41s/it]
tuning metric:0.0063208
tuning hyper-params {'level_sm_input': 0.5, 'seasonality_sm_input': 0.8}
 67%|██████▋   | 6/9 [00:26<00:13,  4.37s/it]
tuning metric:0.0050748
tuning hyper-params {'level_sm_input': 0.8, 'seasonality_sm_input': 0.3}
 78%|███████▊  | 7/9 [00:30<00:08,  4.40s/it]
tuning metric:0.010801
tuning hyper-params {'level_sm_input': 0.8, 'seasonality_sm_input': 0.5}
 89%|████████▉ | 8/9 [00:35<00:04,  4.50s/it]
tuning metric:0.011081
tuning hyper-params {'level_sm_input': 0.8, 'seasonality_sm_input': 0.8}
100%|██████████| 9/9 [00:39<00:00,  4.44s/it]
tuning metric:0.010595

[29]:
tuned_df.head() # backtest output for each parameter searched
[29]:
level_sm_input seasonality_sm_input metrics
0 0.3 0.3 0.004902
1 0.3 0.5 0.004041
2 0.3 0.8 0.003608
3 0.5 0.3 0.007958
4 0.5 0.5 0.006321
[30]:
best_params # output best parameters
[30]:
[{'level_sm_input': 0.3, 'seasonality_sm_input': 0.8}]