HSI
From: Bayesian Models for Astrophysical Data, Cambridge Univ. Press
(c) 2017, Joseph M. Hilbe, Rafael S. de Souza and Emille E. O. Ishida
you are kindly asked to include the complete citation if you used this material in a publication
​
Code 3.4 Bayesian normal linear model in Python.
==================================================
import numpy as np
from pymc3 import Model, sample, summary, traceplot
from pymc3.glm import glm
import pylab as plt
import pandas
from scipy.stats import uniform, norm
​
# Data
np.random.seed(1056) # set seed to replicate example
nobs = 250 # number of obs in model
x1 = uniform.rvs(size=nobs) # random uniform variable
beta0 = 2.0 # intercept
beta1 = 3.0 # angular coefficient
xb = beta0 + beta1 * x1 # linear predictor
y = norm.rvs(loc=xb, scale=1.0, size=nobs) # create y as adjusted
# Fit
df = pandas.DataFrame({'x1': x1, 'y': y}) # rewrite data
with Model() as model_glm:
glm('y ~ x1', df)
trace = sample(5000)
# Output
summary(trace)
# Show graphical output
traceplot(trace)
plt.show()
==================================================
Python 2.7 - Output on screen:
​
Intercept:
Mean SD MC Error 95% HPD interval
-------------------------------------------------------------------
2.040 0.133 0.003 [1.792, 2.284]
Posterior quantiles:
2.5 25 50 75 97.5
|--------------|==============|==============|--------------|
1.794 1.955 2.039 2.128 2.287
x1:
Mean SD MC Error 95% HPD interval
-------------------------------------------------------------------
3.134 0.226 0.005 [2.707, 3.574]
Posterior quantiles:
2.5 25 50 75 97.5
|--------------|==============|==============|--------------|
2.690 2.989 3.134 3.285 3.569
sd_log:
Mean SD MC Error 95% HPD interval
-------------------------------------------------------------------
-0.054 0.045 0.001 [-0.140, 0.037]
Posterior quantiles:
2.5 25 50 75 97.5
|--------------|==============|==============|--------------|
-0.140 -0.085 -0.056 -0.025 0.038
sd:
Mean SD MC Error 95% HPD interval
-------------------------------------------------------------------
0.948 0.043 0.001 [0.865, 1.033]
Posterior quantiles:
2.5 25 50 75 97.5
|--------------|==============|==============|--------------|
0.869 0.918 0.946 0.975 1.038
​
