fixed-sampling-step

fixed-sampling-step.py

import pymc3 as pm
import numpy as np
import matplotlib.pyplot as plt
import arviz as az

from pymc3.step_methods.arraystep import BlockedStep

b = [2,1.5]
sigma = 2

n = 200
x = np.linspace(start=-20,stop = 20, num = n)
y = b[0]*x+b[1]
y_obs = y + sigma*np.random.randn(n)


def my_sampler(mu=0,sigma=2):
    b = np.random.normal(mu,sigma,2)
    return b

class MySamplingStep(BlockedStep):
    def __init__(self, var, mu, sigma):
            self.vars = [var]
            self.var = var
            self.mu = mu
            self.sigma = sigma

    def step(self, point: dict):
        new = point.copy()
        new[self.var.name] = my_sampler(self.mu,self.sigma)

        return new

mu = 0
sigma = 2

with pm.Model() as my_linear_model:
    bm = pm.Normal("bm", mu=0, sigma=2,shape=2)
    step_bm = MySamplingStep(bm, mu, sigma)
    noise = pm.Gamma("noise", alpha=2, beta=1)

    y_observed = pm.Normal("y_observed",mu=bm[0]*x+bm[1],sigma=noise,observed=y_obs)

    posterior2 = pm.sample(step = [step_bm], chains=1, draws=3000, tune=1000,return_inferencedata=True)