liesel.experimental.pymc

A ModelInterface for PyMC models.

To use this module, the PyMC package must be installed. To do so, please install Liesel with the optional dependency PyMC:

$ pip install liesel[pymc]

Example: A linear model

This model is also used in the tests:

RANDOM_SEED = 123
rng = np.random.RandomState(RANDOM_SEED)

# set parameter values
num_obs = 100
sigma = 1.0
beta = [1, 1, 2]

# simulate covariates
x1 = rng.randn(num_obs)
x2 = 0.5 * rng.randn(num_obs)

# simulate outcome variable
y = beta[0] + beta[1] * x1 + beta[2] * x2 + sigma * rng.normal(size=num_obs)

basic_model = pm.Model()
with basic_model:
    # priors
    beta = pm.Normal("beta", mu=0, sigma=10, shape=3)
    sigma = pm.HalfNormal("sigma", sigma=1)
    # sigma is automatically transformed to real (log)
    # the new variable is called sigma_log__

    # predicted value
    mu = beta[0] + beta[1] * x1 + beta[2] * x2

    # track the predicted value of the first obs
    pm.Deterministic("mu[0]", mu[0])

    # distribution of response (likelihood)
    pm.Normal("y_obs", mu=mu, sigma=sigma, observed=y)

interface = PyMCInterface(basic_model, additional_vars=["sigma", "mu[0]"])
state = interface.get_initial_state()
builder = gs.EngineBuilder(1, 2)
builder.set_initial_values(state)
builder.set_model(interface)
builder.set_duration(1000, 2000)

builder.add_kernel(gs.NUTSKernel(["beta"]))
builder.add_kernel(gs.NUTSKernel(["sigma_log__"]))

builder.positions_included = ["sigma", "mu[0]"]

engine = builder.build()

engine.sample_all_epochs()
results = engine.get_results()
sum = gs.Summary(results)
sum

Transformations of RVs can be avoided by setting transform = None as a distribution argument.

Classes

PyMCInterface(model[, additional_vars])

An implementation of the Goose ModelInterface to be used with a PyMC model.