如何在 Python 中用 `pymc.MCMC` 绘制概率分布

Question

我知道我可以使用：

S = pymc.MCMC(model1)
from pymc import Matplot as mcplt
mcplt.plot(S)

这将给我一个包含三个图的图形，但我想要的只是直方图的一个图。然后我想对直方图进行归一化，然后绘制一条平滑的分布曲线而不是直方图的条形图。任何人都可以帮助我编写代码以便我可以获得分布的最终图吗？

Answer 1

如果您安装了 matplotlib 用于绘图，并且 scipy 用于进行核密度估计（many KDE functions exist), then you could do something similar to the following (based on this example，其中 'late_mean' 是采样的名称之一在那种情况下的参数）：

from pymc.examples import disaster_model
from pymc import MCMC
import numpy as np

M = MCMC(disaster_model) # you could substitute your own model

# perform sampling of model
M.sample(iter=10000, burn=1000, thin=10)

# get numpy array containing the MCMC chain of the parameter you want: 'late_mean' in this case
chain = M.trace('late_mean')[:]

# import matplotlib plotting functions
from matplotlib import pyplot as pl

# plot histogram (using 15 bins, but you can choose whatever you want) - density=True returns a normalised histogram
pl.hist(chain, bins=15, histtype='stepfilled', density=True)
ax = pl.gca() # get current axis

# import scipy gaussian KDE function
from scipy.stats import gaussian_kde

# create KDE of samples
kde = gaussian_kde(chain)

# calculate KDE at a range of points (in this case based on the current plot, but you could choose a range based on the chain)
vals = np.linspace(ax.get_xlim()[0], ax.get_xlim()[1], 100)

# overplot KDE
pl.plot(vals, kde(vals), 'b')

pl.xlabel('late mean')
pl.ylabel('PDF')

# show the plot
pl.show()

# save the plot
pl.savefig('hist.png', dpi=200)

如何在 Python 中用 `pymc.MCMC` 绘制概率分布

How to plot a probability distribution with `pymc.MCMC` in Python

histogram

mcmc

python-3.x

pymc