%matplotlib inline
import matplotlib.pyplot as plt
plt.rc('font', size=20) 

import numpy as np
import pymc
print(pymc.__version__)

true_m = 1.0
true_c = 0.0
true_tau = 1.0/(1.5*1.5)

def generate():
    # priors
    x = pymc.Uniform("x", lower=0, upper=10)
    tau = pymc.Gamma("tau", alpha=0.1, beta=0.1, value = true_tau, observed=True)

    m = pymc.Normal("m", mu=0.0, tau=1.0/100, value=true_m, observed=True)
    c = pymc.Normal("c", mu=0.0, tau=1.0/100, value=true_c, observed=True)

    @pymc.deterministic(plot=False)
    def mu(m=m, x=x, c=c):
        return m*x + c
    
    y = pymc.Normal("y", mu=mu, tau=tau)
    
    return locals()

# generate MCMC object
m1 = pymc.MCMC(generate())

m1.sample(iter=1000, burn=100)

tempy = m1.y.trace()
tempx = m1.x.trace()
plt.plot(tempx,tempy,'o')

tempx.shape

tempy.shape

N=10
xdata = tempx[:N]
ydata = tempy[:N]
plt.plot(xdata,ydata,'o')
plt.xlabel('x')
plt.ylabel('y')

print(xdata.size)
print(ydata.size)

def sampling():#(true_m, true_c, true_tau)
    # priors
    x = pymc.Uniform("x", lower=0, upper=10, value=xdata, observed=True)
    tau = pymc.Gamma("tau", alpha=0.1, beta=0.1, value = true_tau, observed=True)

    m = pymc.Normal("m", mu=0.0, tau=1.0/100)
    c = pymc.Normal("c", mu=0.0, tau=1.0/100)

    @pymc.deterministic(plot=False)
    def mu(m=m, x=x, c=c):
        return m*x + c
    
    y = pymc.Normal("y", mu=mu, tau=tau, value=ydata, observed=True)
    
    # posterior prediction. This will calculate values of x and y 
    # which are under the influence of the parameters, but are seperated 
    # from the data
    x_pred = pymc.Uniform("x_pred", lower=0, upper=10)
    y_pred = pymc.Normal("y_pred", mu=mu, tau=tau) # but the 'mu' IS tied to data 
    # if you really wanted though you could set up a set of uniformly spaced x vales
    # and calc mu for that, and get predited y values.  But the way how it is above 
    # will give a distribution of expected y values for each given provided x value.

    return locals()

# generate MCMC object
m2 = pymc.MCMC(sampling())

m2.sample(iter=1000, burn=500)

pymc.Matplot.plot(m2.m)

pymc.Matplot.plot(m2.c)

plt.plot( m2.c.trace() , m2.m.trace() ,'.')
plt.xlabel('c')
plt.ylabel('m')

tempx = m2.x_pred.trace()
tempy = m2.y_pred.trace()
print(tempx.shape)
print(tempy.shape)



plt.plot(m2.x_pred.trace(),
         m2.y_pred.trace(),'o')

pymc.Matplot.plot(m2.x_pred)

pymc.Matplot.plot(m2.y_pred)

m2.y_pred.shape

m2.x_pred.shape

m2.y_pred.parent_names