How can I use the value sampled by pymc3.Discreteuniform for slicing

Questions
1

I make some code by pymc3 for parameter inference of curve fitting about real experimental data.(Sorry for I use dummy data here.)
I have 2 questions.

  1. The data has a change point ob trend, I want to sampling curve fitting parameters and this changing point.
    I code pm.DiscreteUniform for estimate changing point and using sampled value for making array.
    But sampled value cannot use for slicing index. How should I refine my code? And would you please teach me good implementation.
    I find pymc2 can decorate function @pm.deteministic and pymc3 changed using pm.Deterministic.
    I found https://discourse.pymc.io/t/deterministic-rvs-with-slicing-converting-from-pymc-2-to-3/2395
    That can switching simply, how to write some complex logic deterministic function?
    And if you know better model please teach me.
import numpy as np
import pymc3 as pm
import matplotlib.pyplot as plt
np.random.seed(123)

def f1(x1, a1, b1):
    return a1 * x1 + b1


def f2(x2, a2, b2):
    return a2 * x2 + b2


def f(x, t, a1, b1, a2, b2):
    x1 = x[:t]
    x2 = x[t:] - x[t]
    y1 = f1(x1, a1, b1)
    y2 = f2(x2, a2, b2)
    return np.hstack((y1, y2))


def make_dummy_data(x):
    changing_point = 60
    a1 = 0.5
    b1 = 3
    a2 = -3
    b2 = a1 * changing_point - b1
    data = f(x, changing_point, a1, b1, a2, b2)
    data += np.random.normal(scale=10, size=len(x))
    return  data


x = np.arange(0, 100)
y = make_dummy_data(x)
plt.scatter(x, y)
plt.show()

with pm.Model() as model:
    t = pm.DiscreteUniform("t", lower=40, upper=70)
    a1 = pm.Normal("a1", mu=5)
    b1 = pm.Normal("b1", mu=3)
    a2 = pm.Normal("a2", mu=-3)
    b2 = pm.Deterministic("b2", a1 * t - b1)
    sigma = pm.HalfCauchy("sigma", 4)
    mu = pm.Deterministic("mu", f(x, t, a1, a2, b1, b2))
    y_pred = pm.Normal("y_pred", mu=mu, sd=sigma, observed=y)
    start = pm.find_MAP()
    step = pm.NUTS(scaling=start)
    trace_t = pm.sample(5000, step=step, start=start)

---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)
<ipython-input-56-0c37fa1c03ab> in <module>()
      6     b2 = pm.Normal("b2", mu=-33)
      7     sigma = pm.HalfCauchy("sigma", 4)
----> 8     mu = pm.Deterministic("mu", f(x, t, a1, a2, b1, b2))
      9     y_pred = pm.Normal("y_pred", mu=mu, sd=sigma, observed=y)
     10     start = pm.find_MAP()

<ipython-input-55-8087c83b76ce> in f(x, t, a1, b1, a2, b2)
      8 
      9 def f(x, t, a1, b1, a2, b2):
---> 10     x1 = x[:t]
     11     x2 = x[t:] - t
     12     y1 = f1(x1, a1, b1)

TypeError: slice indices must be integers or None or have an __index__ method

> <ipython-input-57-8a02a0e2d9dd>(11)f()
     10     import ipdb;ipdb.set_trace()
---> 11     x1 = x[:t]
     12     x2 = x[t:] - t

ipdb> type(t)
<class 'pymc3.model.FreeRV'>
  1. Actually I have several series of data.
x = np.arange(0, 100)
y = [make_dummy_data(x) for _ in range(10)]

Actual data x is the time.
How can I update my model corresponding to multiple data?
Can I include data as observed=y same way?
I want to estimate the variances of each time (x[0:len(x)]).
Would you teach me good model and how to code?

2

try casting x to a theano tensor using theano.shared

3

Really thank you! Sorry for the easy question.
I faced another bug but your advice helps me understanding about theano.tensor and I could solve problem.

1 Like