How to extend an existing model

Questions
1

If I have a model such as:

y = np.random.randn(10)

with pm.Model() as model_1:
    noise = pm.Gamma('noise', alpha=2, beta=1)
    y_observed = pm.Normal('y_observed', mu=0, sigma=noise, observed=y)

and after exploring this model I want to change the model, to:

with pm.Model() as model_2:
    noise = pm.Gamma('noise', alpha=2, beta=1)
    y_observed = pm.Exponential('y_observed', mu=noise, observed=y)

but I don’t want to rewrite the whole model again (which in reality would be larger, and 99 lines may be the same and one line changed or added). How do I do this? Doing:

with model_1 as model_2:
    y_observed = pm.Exponential('y_observed', mu=noise, observed=y)

results in an error because y_observed already exists. I just want to update its contents. I believe I can add variable by putting the model definition into a custom class and instantiating and building on that, but I don’t know how to update a model by changing variables. This includes Deterministics as well.

2 Likes
2

You can create the def function that will choose either of two models depending on the argument (e.g. you put the if-statement inside that function)

def F(parameter):
  with pm.Model() as model:
    ...
  if (parameter):
    with model:
      ...
  else:
    with model:
      ...
  ...
  return model
3

Thank you, this sounds like a good way to specify different configurations of a model without duplicating code. But is there any way to replace parts of a model (not data) that already exist?

1 Like
4

That’s a really interesting question! I like using a “functional” pattern, though you could do this with graph rewriting in Theano if you have a really specific need. Here’s an example of my approach:

def hierarchical_normal(name, model):
    with model:
        mean = pm.Normal(f'{name}_mean', 0, 1)
        sd = pm.HalfNormal(f'{name}_sd', 1)
    return model

def make_priors(model):
    model = hierarchical_normal('a', model)
    model = hierarchical_normal('b', model)
    return model

def model_version_one():
    model = pm.Model()
    model = make_priors(model)
    with model:
        c = pm.Normal('c', model['a'] + model['b'], observed=data)
    return model
1 Like