# Autoregression#

Here we fit NeuralProphet to data with 5-minute resolution (daily temperatures at Yosemite).

[1]:

if "google.colab" in str(get_ipython()):
!pip install git+https://github.com/ourownstory/neural_prophet.git # may take a while
#!pip install neuralprophet # much faster, but may not have the latest upgrades/bugfixes

import pandas as pd
from neuralprophet import NeuralProphet, set_log_level

# set_log_level("ERROR")

[2]:

data_location = "https://raw.githubusercontent.com/ourownstory/neuralprophet-data/main/datasets/"

[2]:

ds y
0 2017-05-01 00:00:00 27.8
1 2017-05-01 00:05:00 27.0
2 2017-05-01 00:10:00 26.8

## Next-step forecast#

Based on our first contact with the data, we set: * First, we disable weekly_seasonality, as nature does not follow the human week’s calendar. * Second, we increase n_changepoints, and increase changepoints_range, as we are doing short-term predictions.

Further, we can make use of the fact that tomorrow’s weather is most likely similar to yesterdays weather. This means that we will regress our time series on it’s most recent past values, also known as autoregression.

We can achieve this by setting n_lags to the desired number of past observations to regress over. This value is also known as the ‘AR order’.

Here, we predict the temperature in the next 5 minutes based on the last hour:

[3]:

m = NeuralProphet(
n_lags=12,
changepoints_range=0.95,
n_changepoints=30,
weekly_seasonality=False,
)
metrics = m.fit(df)

INFO - (NP.df_utils._infer_frequency) - Major frequency 5T corresponds to 99.995% of the data.
INFO - (NP.df_utils._infer_frequency) - Dataframe freq automatically defined as 5T
INFO - (NP.config.init_data_params) - Setting normalization to global as only one dataframe provided for training.
INFO - (NP.utils.set_auto_seasonalities) - Disabling yearly seasonality. Run NeuralProphet with yearly_seasonality=True to override this.
INFO - (NP.config.set_auto_batch_epoch) - Auto-set batch_size to 64
INFO - (NP.config.set_auto_batch_epoch) - Auto-set epochs to 88

INFO - (NP.utils_torch.lr_range_test) - lr-range-test results: steep: 5.60E-02, min: 7.11E-01

INFO - (NP.utils_torch.lr_range_test) - lr-range-test results: steep: 1.09E-01, min: 6.22E-01
INFO - (NP.forecaster._init_train_loader) - lr-range-test selected learning rate: 8.49E-02
Epoch[88/88]: 100%|█████████████████████████████████████████████| 88/88 [00:47<00:00,  1.85it/s, SmoothL1Loss=6.92e-5, MAE=0.302, RMSE=0.606, RegLoss=0]


Please note that a model with an autoregressive component can be harder to fit. The automatic selection of hyperparameters may not lead to ideal results. For best results, consider changing these manually (in order of importance):

• learning_rate

• epochs

• batch_size

The automatically set hyperparameters (printed out as ‘INFO’ level logs) can serve as a good starting point.

[4]:

forecast = m.predict(df)
fig = m.plot(forecast)

INFO - (NP.df_utils._infer_frequency) - Major frequency 5T corresponds to 99.995% of the data.
INFO - (NP.df_utils._infer_frequency) - Defined frequency is equal to major frequency - 5T
INFO - (NP.df_utils._infer_frequency) - Major frequency 5T corresponds to 99.995% of the data.
INFO - (NP.df_utils._infer_frequency) - Defined frequency is equal to major frequency - 5T


The predictions are now very precise, but this is not a big surprise, as we are just forecasting the very next 5 minutes.

When plotting the model parameters, the panel ‘AR weight’ displays the weight given to the 12 last observed values, these can be interpreted as our ‘AR coefficients’:

[5]:

# fig_comp = m.plot_components(forecast)
m = m.highlight_nth_step_ahead_of_each_forecast(1)  # temporary workaround to plot actual AR weights
fig_param = m.plot_parameters()


The ‘AR weight’ plot shows that the most recent observations are given significantly more weight compared to more distant observations.

## Multi-step forecast#

To predict multiple steps into the future, we could ‘unroll’ our single-step model, by predicting a step ahead, adding the forecasted value to our data, and then forecasting the next step until we reach the horizon we are interested in. However, there is a better way to do this: We can directly forecast multiple steps ahead with NeuralProphet.

We can set n_forecasts to the desired number of steps we would like to forecast (also known as ‘forecast horizon’). NeuralProphet will forecast n_forecasts steps into the future, at every single step. Thus, we have n_forecasts overlapping predictions of vaying age at every historic point.

When increasing the forecast horizon n_forecasts, we should also increase the number of past observations n_lags to at least the same value.

Here, we forecast the next 3 hours based on the last observed 6 hours, in 5-minute steps:

[6]:

m = NeuralProphet(
n_lags=3 * 12,
n_forecasts=2 * 12,
changepoints_range=0.95,
n_changepoints=30,
weekly_seasonality=False,
)
metrics = m.fit(df)

INFO - (NP.df_utils._infer_frequency) - Major frequency 5T corresponds to 99.995% of the data.
INFO - (NP.df_utils._infer_frequency) - Dataframe freq automatically defined as 5T
INFO - (NP.config.init_data_params) - Setting normalization to global as only one dataframe provided for training.
INFO - (NP.utils.set_auto_seasonalities) - Disabling yearly seasonality. Run NeuralProphet with yearly_seasonality=True to override this.
INFO - (NP.config.set_auto_batch_epoch) - Auto-set batch_size to 64
INFO - (NP.config.set_auto_batch_epoch) - Auto-set epochs to 88

INFO - (NP.utils_torch.lr_range_test) - lr-range-test results: steep: 9.56E-02, min: 1.39E+00

INFO - (NP.utils_torch.lr_range_test) - lr-range-test results: steep: 8.37E-02, min: 1.21E+00
INFO - (NP.forecaster._init_train_loader) - lr-range-test selected learning rate: 9.28E-02
Epoch[88/88]: 100%|██████████████████████████████████████████████| 88/88 [01:07<00:00,  1.31it/s, SmoothL1Loss=0.000723, MAE=1.31, RMSE=2.05, RegLoss=0]

[7]:

forecast = m.predict(df)
fig = m.plot(forecast)

INFO - (NP.df_utils._infer_frequency) - Major frequency 5T corresponds to 99.995% of the data.
INFO - (NP.df_utils._infer_frequency) - Defined frequency is equal to major frequency - 5T
INFO - (NP.df_utils._infer_frequency) - Major frequency 5T corresponds to 99.995% of the data.
INFO - (NP.df_utils._infer_frequency) - Defined frequency is equal to major frequency - 5T
WARNING - (NP.plotting.plot) - Legend is available only for the ten first handles


We see that our forecast sometimes miss-predicts for a day and then returns again to more accurate forecasts. Let’s have a closer look at that over-prediction on day 6 of our data:

[8]:

fig = m.plot(forecast[144 : 6 * 288])

WARNING - (NP.plotting.plot) - Legend is available only for the ten first handles


On day 6, the redorded temperature drops at nighttime. Based on the observed lower nighttime temperatures, our model predicts a lower daytime peak. However, the actual daytime temperatures are anomalously low, barely higher than during the night. Thus, leading to an overprediction.

We can again visualize the relative importance of the lags:

[9]:

# fig_comp = m.plot_components(forecast)
fig_param = m.plot_parameters()


Note that the ‘AR relevance’ is the relative importance of a lag, averaged over all n_forecasts.

## Reviewing a specific forecast step#

We can have a closer look at a particular forecast horizon by highlighting the n-th step ahead of each forecast. Here we focus on the temperature predicted 3 hours ahead (most distant prediction into the future). Let’s have a closer look at the weights for the different 3 hour ahead forecast:

[10]:

m = m.highlight_nth_step_ahead_of_each_forecast(2 * 12)
fig_param = m.plot_parameters()


Next, we replot the first 6 days and compare it to the forecast one-step-ahead. We observe the single-step ahead forecast to be much more accurate compared to the 3 hour ahead forecast. However, neither is able to forecast the anomaly on day 6.

[11]:

fig = m.plot(forecast[144 : 6 * 288])
fig = m.plot(forecast[144 : 6 * 288])


## Plotting the most recent forecast#

When we are more concerned with our actual forecast than the model fit, we can plot the last most recent forecast:

[12]:

future = m.make_future_dataframe(df, n_historic_predictions=12)
forecast = m.predict(future)
m = m.highlight_nth_step_ahead_of_each_forecast(None)  # reset highlight
fig = m.plot_last_forecast(forecast)

INFO - (NP.df_utils._infer_frequency) - Major frequency 5T corresponds to 99.995% of the data.
INFO - (NP.df_utils._infer_frequency) - Defined frequency is equal to major frequency - 5T
INFO - (NP.df_utils._infer_frequency) - Major frequency 5T corresponds to 98.611% of the data.
INFO - (NP.df_utils._infer_frequency) - Defined frequency is equal to major frequency - 5T
INFO - (NP.df_utils._infer_frequency) - Major frequency 5T corresponds to 98.611% of the data.
INFO - (NP.df_utils._infer_frequency) - Defined frequency is equal to major frequency - 5T


The last few historic forecasts can be included to assess how the forecast has changed over time. Here, we focus on the 3-hour ahead forecast, given over the last 2 hours.

[13]:

m = m.highlight_nth_step_ahead_of_each_forecast(2 * 12)
fig = m.plot_last_forecast(forecast, include_previous_forecasts=4)


## Larger forecast horizon#

For predictions further into the future, you could reduce the resulution of the data. Using a 5-minute resolution may be useful for a high-resolution short-term forecast, but counter-productive for a long-term forecast. As we only have a limited amount of data (approx 2 months), we want to avoid over-specifying the model.

As an example: If we set the model to forecast 24 hours into the future (nforecasts=24*12) based on the last day’s temperatures (n_lags=24*12), the number of parameters of our AR component grows to 24*12*24*12 = 82,944. However, we only have about 2*30*24*12 = 17,280 samples in our dataset. The model would be overspecified.

If we first downsample our data to hourly data, we reduce our dataset to 2*30*24=1440 and our model parameters to 24*24=576. Thus, we are able to fit the model. However, it would be better to collect more data.

[14]:

df.loc[:, "ds"] = pd.to_datetime(df.loc[:, "ds"])
df_hourly = df.set_index("ds", drop=False).resample("H").mean().reset_index()
len(df_hourly)

[14]:

1561

[15]:

df_hourly.head(3)

[15]:

ds y
0 2017-05-01 00:00:00 24.891667
1 2017-05-01 01:00:00 16.741667
2 2017-05-01 02:00:00 11.733333
[16]:

m = NeuralProphet(
n_lags=24,
n_forecasts=24,
changepoints_range=0.95,
n_changepoints=30,
weekly_seasonality=False,
)
metrics = m.fit(df_hourly)

INFO - (NP.df_utils._infer_frequency) - Major frequency H corresponds to 99.936% of the data.
INFO - (NP.df_utils._infer_frequency) - Dataframe freq automatically defined as H
INFO - (NP.config.init_data_params) - Setting normalization to global as only one dataframe provided for training.
INFO - (NP.utils.set_auto_seasonalities) - Disabling yearly seasonality. Run NeuralProphet with yearly_seasonality=True to override this.
INFO - (NP.config.set_auto_batch_epoch) - Auto-set batch_size to 32
INFO - (NP.config.set_auto_batch_epoch) - Auto-set epochs to 171

INFO - (NP.utils_torch.lr_range_test) - lr-range-test results: steep: 8.52E-02, min: 4.99E-01

INFO - (NP.utils_torch.lr_range_test) - lr-range-test results: steep: 7.25E-02, min: 6.87E-01
INFO - (NP.forecaster._init_train_loader) - lr-range-test selected learning rate: 8.52E-02
Epoch[171/171]: 100%|████████████████████████████████████████████| 171/171 [00:20<00:00,  8.52it/s, SmoothL1Loss=0.00202, MAE=2.13, RMSE=3.4, RegLoss=0]

[17]:

forecast = m.predict(df_hourly)
fig = m.plot(forecast)
# fig_param = m.plot_parameters()

INFO - (NP.df_utils._infer_frequency) - Major frequency H corresponds to 99.936% of the data.
INFO - (NP.df_utils._infer_frequency) - Defined frequency is equal to major frequency - H
INFO - (NP.df_utils._infer_frequency) - Major frequency H corresponds to 99.937% of the data.
INFO - (NP.df_utils._infer_frequency) - Defined frequency is equal to major frequency - H
WARNING - (NP.plotting.plot) - Legend is available only for the ten first handles


Finally, we plot the most recent and last ten historic 24-hour forecasts, marking the 24-th hour ahead with an ‘x’.

[18]:

future = m.make_future_dataframe(df_hourly)
forecast = m.predict(future)

INFO - (NP.df_utils._infer_frequency) - Major frequency H corresponds to 99.936% of the data.