Skip to content
Chris Fonnesbeck edited this page Apr 12, 2016 · 1 revision
"""
Simple hierarchical linear model of price as a function of age.
"""

from pymc import *
from numpy import array

# Data
age = array([13, 14, 14,12, 9, 15, 10, 14, 9, 14, 13, 12, 9, 10, 15, 11, 15, 11, 7, 13, 13, 10, 9, 6, 11, 15, 13, 10, 9, 9, 15, 14, 14, 10, 14, 11, 13, 14, 10])
price = array([2950, 2300, 3900, 2800, 5000, 2999, 3950, 2995, 4500, 2800, 1990, 3500, 5100, 3900, 2900, 4950, 2000, 3400, 8999, 4000, 2950, 3250, 3950, 4600, 4500, 1600, 3900, 4200, 6500, 3500, 2999, 2600, 3250, 2500, 2400, 3990, 4600, 450,4700])/1000.

"""Original WinBUGS model:

model
	{
		for( i in 1 : N ) {
    		mu[i] <- a+b*age[i]
			price[i] ~ dnorm(mu[i], tau)		
		}	
		tau ~ dgamma(0.001, 0.001) 
		sigma <- 1 / sqrt(tau)
		a ~ dnorm(0, 1.0E-12)
		b ~ dnorm(0, 1.0E-12)
	}"""
	
# Constant priors for parameters
a = Normal('a', 0, 0.0001)
b = Normal('b', 0, 0.0001)

# Precision of normal distribution of prices
tau = Gamma('tau', alpha=0.1, beta=0.1)

@deterministic
def mu(x=age, a=a, b=b):
    # Linear age-price model
    return a + b*x

    
# Sampling distribution of prices
p = Normal('p', mu, tau, value=price, observed=True)
Clone this wiki locally