Updating fitted models fails if mcmc_samples > 0

[raffg]

The Updating fitted models section of the documentation describes how to warm start a model with parameters from a previously-trained model. This works fine of the models are trained with mcmc_samples=0 but if MCMC sampling is being used then the procedure as described in the documentation fails with IndexError: invalid index to scalar variable.

df = pd.read_csv('../examples/example_wp_log_peyton_manning.csv')
df1 = df.loc[df['ds'] < '2016-01-19', :]  # All data except the last day
m1 = Prophet(mcmc_samples=300).fit(df1) # A model fit to all data except the last day

%timeit m2 = Prophet(mcmc_samples=300).fit(df)  # Adding the last day, fitting from scratch
%timeit m2 = Prophet(mcmc_samples=300).fit(df, init=stan_init(m1))  # Adding the last day, warm-starting from m1

IndexError: invalid index to scalar variable.

The shapes of each parameter are different between MAP and MCMC: MAP, k: (1, 1) MCMC, k: (600,) MAP, m: (1, 1) MCMC, m: (600,) MAP, delta: (1, 25) MCMC, delta: (600, 25) MAP, sigma_obs: (1, 1) MCMC, sigma_obs: (600,) MAP, beta: (1, 26) MCMC, beta: (600, 26)

To solve this problem, I modified the stan_init function as I believe is appropriate:

def stan_init(m):
    res = {}
    if m.mcmc_samples == 0:
        for pname in ['k', 'm', 'sigma_obs']:
            res[pname] = m.params[pname][0][0]
        for pname in ['delta', 'beta']:
            res[pname] = m.params[pname][0]
    else:
        for pname in ['k', 'm', 'sigma_obs']:
            res[pname] = m.params[pname]
        for pname in ['delta', 'beta']:
            res[pname] = m.params[pname]       
    return res

But I am now getting ValueError: Invalid specification of initial values..

How can I get the warm start technique working with MCMC sampling?

A second thought I had was, is there any way to pre-train a model with MAP estimation, which is faster, and then use those parameters to warm start a model with MCMC sampling? I know the number of changepoints will most likely be different and this will cause problems, but is there any way to use some of that MAP information to speed things up with MCMC?

[bletham]

Sorry for the slow reply. The shape of the inputs will be the same whether you're doing MCMC or not. k, m, and sigma_obs need to be floats; delta and beta need to be 1-d arrays. So you'd need to turn it into a point estimate, like the posterior mean.

Or as you note you can fit a model with mcmc_samples=0 to get the MAP estimate, and then should be able to use the warm-starting code directly as given there while setting mcmc_samples to something >0 to do MCMC warm-started from the MAP estimate. And I agree that that makes sense, I will frequently warm-start MCMC from the MAP solution.

[raffg]

No worries, and thanks for getting back to this! I'm still not able to get it to work though. Do you mean I need to update the stan_init function with something like this?

def stan_init(m):
    res = {}
    for pname in ['k', 'm', 'sigma_obs']:
        res[pname] = np.mean(m.params[pname])  # when mcmc_samples=0, this is the same as m.params[pname][0][0]
    for pname in ['delta', 'beta']:
        res[pname] = np.mean(m.params[pname], axis=0)  # when mcmc_samples=0, this is the same as m.params[pname][0]
    return res

That results in the same shape for all parameters whether MAP or MCMC. But I still get an error with it:

df = pd.read_csv('../examples/example_wp_log_peyton_manning.csv')
df1 = df.loc[df['ds'] < '2016-01-19', :]  # All data except the last day

# MAP works fine with no errors
m1_map = Prophet().fit(df1) # A model fit to all data except the last day
m2_map = Prophet().fit(df, init=stan_init(m1_map))  # Adding the last day, warm-starting from model1

# MCMC results in a ValueError
m1_mcmc = Prophet(mcmc_samples=300).fit(df1)
m2_mcmc = Prophet(mcmc_samples=300).fit(df, init=stan_init(m1_mcmc))  # Adding the last day, warm-starting from model1

I get an error on that last line:

---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-27-91b8a1bad3f6> in <module>
      8 # MCMC results in a ValueError
      9 m1_mcmc = Prophet(mcmc_samples=300).fit(df1)
---> 10 m2_mcmc = Prophet(mcmc_samples=300).fit(df, init=stan_init(m1_mcmc))  # Adding the last day, warm-starting from model1

~\miniconda3\lib\site-packages\fbprophet\forecaster.py in fit(self, df, **kwargs)
   1162                 self.params[par] = np.array([self.params[par]])
   1163         elif self.mcmc_samples > 0:
-> 1164             self.params = self.stan_backend.sampling(stan_init, dat, self.mcmc_samples, **kwargs)
   1165         else:
   1166             self.params = self.stan_backend.fit(stan_init, dat, **kwargs)

~\miniconda3\lib\site-packages\fbprophet\models.py in sampling(self, stan_init, stan_data, samples, **kwargs)
    224         )
    225         args.update(kwargs)
--> 226         self.stan_fit = self.model.sampling(**args)
    227         out = dict()
    228         for par in self.stan_fit.model_pars:

~\miniconda3\lib\site-packages\pystan\model.py in sampling(self, data, pars, chains, iter, warmup, thin, seed, init, sample_file, diagnostic_file, verbose, algorithm, control, n_jobs, **kwargs)
    794                                               diagnostic_file=diagnostic_file,
    795                                               algorithm=algorithm,
--> 796                                               control=control, **kwargs)
    797 
    798         # number of samples saved after thinning

~\miniconda3\lib\site-packages\pystan\misc.py in _config_argss(chains, iter, warmup, thin, init, seed, sample_file, diagnostic_file, algorithm, control, **kwargs)
    481         inits_specified = True
    482     if not inits_specified:
--> 483         raise ValueError("Invalid specification of initial values.")
    484 
    485     ## only one seed is needed by virtue of the RNG

ValueError: Invalid specification of initial values.

Do you know where I'm going wrong?

And a secondary question: if I build a model with MAP and use it to warm start an MCMC model, won't the number of changepoints most likely be different? And therefore, the length of m.params['delta'] will also be different, resulting in that same ValueError as above?

[bletham]

OK I figured out what the issue was. There were two issues - the first was the shape issue, where k, m, and sigma_obs needed to be floats, and delta and beta needed to be 1-d arrays. But once this was fixed there was an additional issue that for sampling, PyStan wants the initial specificiation to be a callable, and not just a dictionary (with MAP it can be either). In fbprophet the stan backend is converting the dict of intial conditions into a callable here: https://github.com/facebook/prophet/blob/4f34de036390bc0e66bdc9ffe2d89361ccbc3f07/python/fbprophet/models.py#L233 So we just have to pass it in as a callable too and then everything works. Here is working code. Note also the modification to the init function; this function works for both MAP and MCMC.

from fbprophet import Prophet
import pandas as pd
import time
import numpy as np

# Python
def stan_init_warmstart(m):
    res = {}
    for pname in ['k', 'm', 'sigma_obs']:
        res[pname] = np.mean(m.params[pname])
    for pname in ['delta', 'beta']:
        res[pname] = np.mean(m.params[pname], axis=0)
    return res

df = pd.read_csv('../examples/example_wp_log_peyton_manning.csv')
# Shorten so this example runs faster
df = df[df['ds'] > '2013-01-01'].copy()
df1 = df.loc[df['ds'] < '2016-01-19', :]  # All data except the last day
m1 = Prophet().fit(df1) # A model fit to all data except the last day

# MAP fit with no warm-start
t1 = time.time()
m1 = Prophet().fit(df)  # Adding the last day, fitting from scratch
print(time.time() - t1)

# MCMC with no warm-start
t1 = time.time()
m2 = Prophet(mcmc_samples=200).fit(df)
print(time.time() - t1)

# MAP fit warm-start
t1 = time.time()
m3 = Prophet().fit(df, init=stan_init_warmstart(m1))
print(time.time() - t1)

# MCMC warm-started from posterior mean of previous MCMC
t1 = time.time()
m4 = Prophet(mcmc_samples=200).fit(df, init=lambda: stan_init_warmstart(m1))
print(time.time() - t1)

# MCMC warm-started from MAP
m5 = Prophet().fit(df)

t1 = time.time()
m6 = Prophet(mcmc_samples=200).fit(df, init=lambda: stan_init_warmstart(m5))
print(time.time() - t1)

Note that warm-starting MCMC doesn't make it faster, because it still does the same specified number of samples; it just probably makes the chains better, so then you might be able to get away with fewer samples and less burn-in.

As for your secondary question on warm-starting MCMC with the MAP estimate: the number of changepoints is fixed, and is the kwarg n_changepoints that defaults to 25.

[raffg]

Ah, perfect! The lambda was the key. Thanks so much! It would have taken me months to figure that out.