Probability: Forecasting with Structural Time Series

[cmpgervais]

Referencing: https://github.com/tensorflow/probability/issues/343#issuecomment-480952321

Trying to clarify my understanding of best practices for forecasting with STS. So I'm good up until:

q_samples_demand_ = sess.run({k: q.sample(50)
                                for k, q in variational_posteriors.items()})

From this I can get an estimate of the inferred parameters:

Inferred parameters:
observation_noise_scale: 5.289973735809326 +- 0.02318570390343666
hour_of_day_effect/_drift_scale: 5.65079402923584 +- 0.018391167744994164
day_of_week_effect/_drift_scale: 2.2549798488616943 +- 0.02803918346762657
month_of_year_effect/_drift_scale: 0.3897028863430023 +- 0.037868596613407135
x_effect/_weights: [ 5.2341375  -2.3535817   0.51498556] +- [0.17526183 0.05769397 0.78600216]
autoregressive/_coefficients: [0.88036] +- [0.0036582]
autoregressive/_level_scale: 9.465166091918945 +- 0.018752284348011017

Then, before actually building the sts.forecast() object, I should create a new model object with the future regressor values in the design_matrix:

def build_model_forecast(observed_time_series):

    hour_of_day_effect = sts.Seasonal(
        num_seasons = 24,
        num_steps_per_season = 1,
        observed_time_series = observed_time_series,
        name = 'hour_of_day_effect')

    day_of_week_effect = sts.Seasonal(
        num_seasons = 7,
        num_steps_per_season = 24,
        observed_time_series = observed_time_series,
        name = 'day_of_week_effect')

    month_of_year_effect = tfp.sts.Seasonal(
        num_seasons=12,
        num_steps_per_season= 30 * 24,
        drift_scale_prior=tfd.LogNormal(loc=-1., scale=0.1),
        initial_effect_prior=tfd.Normal(loc=0., scale=5.),
        name='month_of_year_effect'
    )

    x_effect = sts.LinearRegression(
        design_matrix = tf.stack([load, wind, outage], axis = -1),
        weights_prior = tfd.Normal(loc=0., scale=1.),
        name='x_effect')

    autoregressive = sts.Autoregressive(
        order=1,
        observed_time_series=observed_time_series,
        name='autoregressive')

    model = sts.Sum([hour_of_day_effect,
                   day_of_week_effect,
                   month_of_year_effect,
                   x_effect,
                   autoregressive],
                   observed_time_series=observed_time_series)

    return model

Then at forecast time I would pass the original model parameter samples, along with the model containing future regressor values:

demand_forecast_dist = tfp.sts.forecast(
    model=build_model_forecast(),
    observed_time_series=target_training_data,
    parameter_samples=q_samples_demand_,
    num_steps_forecast=num_forecast_steps)

Does that seem about right? It feels a bit awkward, but other than appending the original model design matrix directly, I haven't found any alternatives.

[skeydan]

Hi,

I think it's okay to define the model using the complete design matrix, but for training, make use of the training data only (for loss calculation). At least that's what I did here: https://blogs.rstudio.com/tensorflow/posts/2019-06-25-dynamic_linear_models_tfprobability/

[trainorpj]

This is interesting, and you've stated the problem nicely. One suggestion I have in order to prevent data leakage (e.g. test data in your training set) is to define a general build_model, which accepts a x_effect_design_matrix parameter. It would look like:

def build_model(observed_time_series, x_effect_design_matrix):

    hour_of_day_effect = sts.Seasonal(...)

    day_of_week_effect = sts.Seasonal(...)

    month_of_year_effect = tfp.sts.Seasonal(...)

    autoregressive = sts.Autoregressive(...)

    x_effect = sts.LinearRegression(
        design_matrix = x_effect_design_matrix,
        weights_prior = tfd.Normal(loc=0., scale=1.),
        name='x_effect')

    model = sts.Sum([hour_of_day_effect,
                   day_of_week_effect,
                   month_of_year_effect,
                   x_effect,
                   autoregressive],
                   observed_time_series=observed_time_series)

    return model

Then you could use it like:

# get data
training_time_series, training_x_effect = get_training_data()
forecast_time_series, forecast_x_effect = get_forecast_data()

# train model
training_model = build_model(training_time_series, training_x_effect)
q_samples_ = train_model(training_model)

# forecast
forecast_model = build_model(forecast_time_series, forecast_x_effect)
forecast_dist = tfp.sts.forecast(
    model=forecast_model,
    observed_time_series=forecast_time_series,
    parameter_samples=q_samples_,
    num_steps_forecast=num_forecast_steps)

Do you think that would work for you?

[cmpgervais]

Thanks @trainorpj! I haven't gotten around to finishing the inference part yet, but that's a great suggestion - will definitely try it out! Would be interesting to see this functionality built into TFP natively, I suspect others may run into the same scenario.

[CodeMySky]

Hi, if I try to build a second model during inference and pass in the previous model params. The second model params will actually get a '_1' suffix so they will not be matched to the naming of params for the first model. Any suggestions?

[cmpgervais]

@CodeMySky can you post your code so we can try to reproduce? I didn't have any issues when re-building the model during inference, but not sure what version I was working with at the time.

[CodeMySky]

In the Structural_Time_Series_Modeling_Case_Studies_Atmospheric_CO2_and_Electricity_Demand.ipynb

If I run the original code, which is:

demand_forecast_dist = tfp.sts.forecast(
    model=demand_model,
    observed_time_series=demand_training_data,
    parameter_samples=q_samples_demand_,
    num_steps_forecast=num_forecast_steps)

It works fine, if I create a new model, it will raise an error:

demand_model_new = build_model(demand)
demand_forecast_dist = tfp.sts.forecast(
    model=demand_model_new,
    observed_time_series=demand_training_data,
    parameter_samples=q_samples_demand_,
    num_steps_forecast=num_forecast_steps)

The error is:

<ipython-input-12-4660e0488b68> in <module>
     10     observed_time_series=demand_training_data,
     11     parameter_samples=q_samples_demand_,
---> 12     num_steps_forecast=24 * 7 * 2)

/anaconda3/lib/python3.7/site-packages/tensorflow_probability/python/sts/forecast.py in forecast(model, observed_time_series, parameter_samples, num_steps_forecast)
    304         tf.shape(input=observed_time_series))[-2]
    305     observed_data_ssm = model.make_state_space_model(
--> 306         num_timesteps=num_observed_steps, param_vals=parameter_samples)
    307     (_, _, _, predictive_means, predictive_covs, _, _
    308     ) = observed_data_ssm.forward_filter(observed_time_series, mask=mask)

/anaconda3/lib/python3.7/site-packages/tensorflow_probability/python/sts/structural_time_series.py in make_state_space_model(self, num_timesteps, param_vals, initial_state_prior, initial_step)
    157         param_map=self._canonicalize_param_vals_as_map(param_vals),
    158         initial_state_prior=initial_state_prior,
--> 159         initial_step=initial_step)
    160 
    161   def prior_sample(self,

/anaconda3/lib/python3.7/site-packages/tensorflow_probability/python/sts/sum.py in _make_state_space_model(self, num_timesteps, param_map, initial_step, initial_state_prior, constant_offset)
    524 
    525     # List the model parameters in canonical order
--> 526     param_vals_list = [param_map[p.name] for p in self.parameters]
    527     observation_noise_scale = param_vals_list[0]
    528 

/anaconda3/lib/python3.7/site-packages/tensorflow_probability/python/sts/sum.py in <listcomp>(.0)
    524 
    525     # List the model parameters in canonical order
--> 526     param_vals_list = [param_map[p.name] for p in self.parameters]
    527     observation_noise_scale = param_vals_list[0]
    528 

KeyError: 'monthly_seasonality_1/_drift_scale'

[davmre]

In graph-mode TF, you can fix this by building the inference model in a
separate graph from the training model (either call
tf.reset_default_graph() in colab, or wrap the training and inference code
in separate with tf.Graph(): contexts).

On Thu, Sep 12, 2019 at 10:47 AM Sean Han [e-mail] wrote:

In the
Structural_Time_Series_Modeling_Case_Studies_Atmospheric_CO2_and_Electricity_Demand.ipynb

If I run the original code, which is:

demand_forecast_dist = tfp.sts.forecast(
model=demand_model,
observed_time_series=demand_training_data,
parameter_samples=q_samples_demand_,
num_steps_forecast=num_forecast_steps)

It works fine, if I create a new model, it will raise an error:

demand_model_new = build_model(demand)
demand_forecast_dist = tfp.sts.forecast(
model=demand_model_new,
observed_time_series=demand_training_data,
parameter_samples=q_samples_demand_,
num_steps_forecast=num_forecast_steps)

The error is:

in
10 observed_time_series=demand_training_data,
11 parameter_samples=q_samples_demand_,---> 12 num_steps_forecast=24 * 7 * 2)
/anaconda3/lib/python3.7/site-packages/tensorflow_probability/python/sts/forecast.py in forecast(model, observed_time_series, parameter_samples, num_steps_forecast)
304 tf.shape(input=observed_time_series))[-2]
305 observed_data_ssm = model.make_state_space_model(--> 306 num_timesteps=num_observed_steps, param_vals=parameter_samples)
307 (_, _, _, predictive_means, predictive_covs, _, _
308 ) = observed_data_ssm.forward_filter(observed_time_series, mask=mask)
/anaconda3/lib/python3.7/site-packages/tensorflow_probability/python/sts/structural_time_series.py in make_state_space_model(self, num_timesteps, param_vals, initial_state_prior, initial_step)
157 param_map=self._canonicalize_param_vals_as_map(param_vals),
158 initial_state_prior=initial_state_prior,--> 159 initial_step=initial_step)
160
161 def prior_sample(self,

/anaconda3/lib/python3.7/site-packages/tensorflow_probability/python/sts/sum.py in _make_state_space_model(self, num_timesteps, param_map, initial_step, initial_state_prior, constant_offset)
524
525 # List the model parameters in canonical order--> 526 param_vals_list = [param_map[p.name] for p in self.parameters]
527 observation_noise_scale = param_vals_list[0]
528
/anaconda3/lib/python3.7/site-packages/tensorflow_probability/python/sts/sum.py in (.0)
524
525 # List the model parameters in canonical order--> 526 param_vals_list = [param_map[p.name] for p in self.parameters]
527 observation_noise_scale = param_vals_list[0]
528
KeyError: 'monthly_seasonality_1/_drift_scale'

—
You are receiving this because you are subscribed to this thread.
Reply to this email directly, view it on GitHub
https://github.com/tensorflow/probability/issues/508?email_source=notifications&email_token=AAHSFCV6WSG6YT4UH57SUYDQJJ6EHA5CNFSM4IHDRAGKYY3PNVWWK3TUL52HS4DFVREXG43VMVBW63LNMVXHJKTDN5WW2ZLOORPWSZGOD6SWKFA#issuecomment-530933012,
or mute the thread
https://github.com/notifications/unsubscribe-auth/AAHSFCWJ2MABASW7G6NZQOLQJJ6EHANCNFSM4IHDRAGA
.

[vz415]

I wanted to hop onto this thread and ask a rookie question that may be relevant enough to avoid having to open a new issue, but let me know if I should open a new issue. I'm trying to forecast by using a combination of semilocal linear trend and linear regression. For the linear regression, I want to change the prior distribution from normal to something like an exponential weighted gaussian, so as to shift the curve to something experts think is more likely than just a plain normal distribution. How do I create a custom distribution that is a normal distribution summed with random exponential distribution?

I browsed through the tutorial on distributions and this didn't seem obvious to me at first glance. If I missed something or if anyone can help, I'd appreciate it.

[brianwa84]

Unless you have an issue with the sts models requiring your prior to be
gaussian, I'd say to file a new issue.
At a glance, the sum of two distributions is hard, but it seems like a
special case exists for this:
https://en.wikipedia.org/wiki/Exponentially_modified_Gaussian_distribution
You could also take a look at
https://github.com/tensorflow/probability/pull/445 , which could get you
partway there.

Brian Patton | Software Engineer | [email protected]

On Sun, Sep 29, 2019 at 5:08 PM vz415 [e-mail] wrote:

I wanted to hop onto this thread and ask a rookie question that may be
relevant enough to avoid having to open a new issue, but let me know if I
should open a new issue. I'm trying to forecast by using a combination of
semilocal linear trend and linear regression. For the linear regression, I
want to change the prior distribution from normal to something like an
exponential weighted gaussian, so as to shift the curve to something
experts think is more likely than just a plain normal distribution. How do
I create a custom distribution that is a normal distribution summed with
random exponential distribution?

I browsed through the tutorial on distributions and this didn't seem
obvious to me at first glance. If I missed something or if anyone can help,
I'd appreciate it.

—
You are receiving this because you are subscribed to this thread.
Reply to this email directly, view it on GitHub
https://github.com/tensorflow/probability/issues/508?email_source=notifications&email_token=AFJFSI45LQEDSJSTF4BE3GLQMEKM5A5CNFSM4IHDRAGKYY3PNVWWK3TUL52HS4DFVREXG43VMVBW63LNMVXHJKTDN5WW2ZLOORPWSZGOD736WPY#issuecomment-536341311,
or mute the thread
https://github.com/notifications/unsubscribe-auth/AFJFSIZMFR42FOZH6QUVYBTQMEKM5ANCNFSM4IHDRAGA
.

[vz415]

Thanks, Brian. I'll play around with customizing it thise weekend. #445 does look like it addresses my issue, but I'll see if I can play around with priors to see if you can enable multiplication/customization via a new PR or new issue. No promises on this happening soon - a little busy!

[cmp1]

Hello all,

I have quite a sparse label that I would like to predict a long horizon for, so have been interested in building this model as an experiment.

My question here however is whether it possible to perform multi-step forecasts with this configuration? i.e. predict the time series out to 365 days from today (t+365)?

Thank you!

[yug95]

hey i have a doubt how to include external regressor variable while forecasting ?

[vz415]

@cmp1 Yes, that's possible but your variance might be quite large depending on the size of your historical dataset.

@yug95 Check out this page to create a model that includes external regressor variables in your forecast. https://www.tensorflow.org/probability/api_docs/python/tfp/sts/LinearRegression

[cmp1]

@vz415 Thank you for this. Do you have any examples of where this has been performed? taking the lagged approach for say t+252 would be unwieldy unless we were using monthly data perhaps.

[vz415]

@cmp1 I tihnk the provided tensorflow probability article is a good start but I'm not sure if it's applicable to your data or problem, since it sounds like you want to predict much farther out as compared to historical data. Maybe open another issue or comment somewhere on stack exchange? As for converting to monthly - that could be an option but you lose a lot of information and it just depends on if you're ok with a simplified model. Hope this helps.

A broader question for this issue/channel - is there any interest in including a RMSE and/or more validation functions to validate forecasts on previous data? I wrote my own little script and would be happy to share/expand on it should it be of interest.