Sample method of pymc3

Questions
1

I’m new with pymc3 and beginner in probabilistic programming. I wrote a program for testing but it didn’t work. I have a problem in running of last line (sample method). after several hours it doesn’t still work and I can’t deal with it.

import numpy as np
import matplotlib.pyplot as plt
from pymc3 import Model, Normal, HalfNormal, find_MAP, NUTS, sample
from pymc3 import traceplot, summary
from scipy import optimize
import pandas as pd
from pymc3 import Exponential, StudentT, Deterministic
from pymc3.distributions.timeseries import GaussianRandomWalk
from theano.tensor import exp

np.random.seed(123)
alpha, sigma = 1, 1
beta = [1, 2.5]
size = 100
X1 = np.linspace(0, 1, size)
X2 = np.linspace(0,.2, size)
Y = alpha + beta[0]*X1 + beta[1]*X2 + np.random.randn(size)*sigma

returns = pd.read_csv('SP500.csv', index_col=0, parse_dates=True, sep='delimiter', header=None, engine='python')

basic_model = Model()
with basic_model:
    alpha = Normal('alpha', mu=0, sd=10)
    beta = Normal('beta', mu=0, sd=10, shape=2)
    sigma = HalfNormal('sigma', sd=1)
    mu = alpha + beta[0] * X1 + beta[1] * X2
    Y_obs = Normal('Y_obs', mu=mu, sd=sigma, observed=Y)

with Model() as sp500_model:
    nu = Exponential('nu', 1. / 10, testval=5.)
    sigma = Exponential('sigma', 1. / .02, testval=.1)
    s = GaussianRandomWalk('s', sigma - 2, shape=len(returns))
    volatility_process = Deterministic('volatility_process', exp(-2 * s))
    r = StudentT('r', nu, lam=1 / volatility_process)

    start = find_MAP(vars=[s], fmin=optimize.fmin_l_bfgs_b)
    step = NUTS(scaling=start)
    trace = sample(100, step, progressbar=False)

the run of program stopped in this point and I was waiting more than 5 hours but didn’t work.

console:

Only 100 samples in chain.
Multiprocess sampling (4 chains in 4 jobs)
NUTS: [r, s, sigma_log__, nu_log__]

please help me, what should I do?

my OS is Win 10 and my IDE is spyder of Anaconda3.

2

You did not assign any observed value in your sp500_model, is that intentional? If so it means you want to sample from the prior, which you should use pm.sample_prior_predictive, which generate sample from the prior distribution.
Otherwise, assign a variable for your observed value return, and sample with:

trace = sample(1000, tune=1000)

DO NOT use find_MAP(…) and step = NUTS(scaling=start), they usually result in bad initual statue of the transistional kernel of the sampler

3

thank you, I’ll try this and I hope my problem will be solved.
actually I’m learning pymc3 and I copied this code from a paper for testing and they didn’t assign any observed value (I don’t know this is nessesary or optional). for now, I want some ususal outputs. my teacher wanted a simple implemented program by pymc3.
you mean I don’t use find_MAP(…) and step = NUTS(scaling=start) and just trace = sample(1000, tune=1000) is enough?

4

I will suggest you to have a look at the examples in the doc.

Yes

5

thank you, I will do that