8schools.md

+++ title= "Eight Schools Model" date= 2018-07-05T10:19:11+05:30 draft= false tags=["gsoc", "gsoc18", "machine learning", "pymc"] categories=["gsoc18", "gsoc"] summary= """Non-Centered Eight schools model PyMC3 vs PyMC4""" math= false +++

We are finally at a state where we can demonstrate the use of the PyMC4 API side by side with PyMC3 and showcase the consistency in results by using non-centered eight schools model.

I will be comparing the PyMC3 and PyMC4 way of doing the same task. PyMC3 follows a context manager based approach whereas PyMC4 uses a functional approach to API design.

Imports

PyMC3

import pymc3 as pm
import numpy as np
import theano.tensor as tt # For Random Variable operation
        

PyMC4

import pymc4 as pm
import numpy as np
import tensorflow as tf # For Random Variable operation
from tensorflow_probability import edward2 as ed # For defining random variables
        

Observed data and configuration

PyMC3

J = 8 # No. of schools
y = np.array([28.,  8., -3.,  7., -1.,  1., 18., 12.])
sigma = np.array([15., 10., 16., 11.,  9., 11., 10., 18.])
        

PyMC4

J = 8 # No. of schools
y = np.array([28.,  8., -3.,  7., -1.,  1., 18., 12.])
sigma = np.array([15., 10., 16., 11.,  9., 11., 10., 18.])
        

Model Initialization

PyMC3

NonCentered_eight = pm.Model()
        

PyMC4

NonCentered_eight = pm.Model(num_schools=J, y=y, sigma=sigma )
# These attributes can be later used in model definition and even be modified using configure
        

Model Definition

PyMC3

with NonCentered_eight:
mu = pm.Normal('mu', mu=0, sd=5)
log_tau = pm.Normal('log_tau', mu=5, sd=1)
theta_prime = pm.Normal('theta_prime', mu=0, sd=1, shape=J)
theta = mu + tt.exp(log_tau) * theta_prime
obs = pm.Normal('obs', mu=theta, sd=sigma, observed=y)
trace = pm.sample()

theta = trace['mu'][:, np.newaxis] + np.exp(trace['log_tau'])[:, np.newaxis] * trace['theta_prime']
theta.mean(axis=0)

PyMC4

@model.define
def process(cfg):
mu = ed.Normal(loc=0., scale=5., name="mu")
log_tau = ed.Normal(
loc=5., scale=1., name="log_tau")
theta_prime = ed.Normal(
loc=tf.zeros(cfg.num_schools),
scale=tf.ones(cfg.num_schools),
name="theta_prime")
theta = mu + tf.exp(
log_tau) * theta_prime
y = ed.Normal(
loc=theta,
scale=np.float32(cfg.sigma),
name="y")

return y

theta = (
trace['mu'][:, np.newaxis] +
np.exp(trace['log_tau'])[:, np.newaxis] * trace['theta_prime'])
theta.mean(axis=0)

Observe Variables

PyMC3

# During model definition
        

PyMC4

NonCentered_eight.observe(y=y)
        

Sampling

PyMC3

with NonCentered_eight:
    trace = pm.sample()
        

PyMC4

trace = pm.sample(NonCentered_eight)
        

Calculating theta from theta_prime

PyMC3

theta = trace['mu'][:, np.newaxis] + np.exp(trace['log_tau'])[:, np.newaxis] * trace['theta_prime']
# Can be done inline using `theta = pm.Deterministic('theta', mu + tau * theta_tilde)`
        

PyMC4

theta = trace['mu'][:, np.newaxis] + np.exp(trace['log_tau'])[:, np.newaxis] * trace['theta_prime']