【MCMC】pymc2库估计逻辑回归参数

最新推荐文章于 2022-11-26 22:36:40 发布

Quant0xff

最新推荐文章于 2022-11-26 22:36:40 发布

阅读量567

点赞数

分类专栏： # 机器学习 # 统计推断与概率计算 # Python 文章标签： pyMC2

本文链接：https://blog.csdn.net/qq_18822147/article/details/108773664

版权

机器学习同时被 3 个专栏收录

34 篇文章 3 订阅

订阅专栏

Python

27 篇文章 1 订阅

订阅专栏

统计推断与概率计算

13 篇文章 0 订阅

订阅专栏

该博客介绍了在小白鼠实验中如何利用MCMC（Markov Chain Monte Carlo）方法来估计逻辑回归模型的参数。通过设置观测数据、先验分布并定义似然函数，博主展示了如何进行参数推断，并最终绘制了参数的统计结果和风险随药物剂量变化的曲线。

摘要由CSDN通过智能技术生成

Estimating parameters of a logistic model

在小白鼠实验中，设有 $n = 5$ 只小白鼠，每只小白鼠的死亡概率为 $p$ ，可以建立如下logit模型
$\left\{ \begin{aligned} &y\sim Bin(n, p)\\ &logit(p)=\alpha+\beta x\\ &\alpha\sim\mathcal{N}(0, 5)\\ &\beta\sim\mathcal{N}(0, 10) \end{aligned} \right.$
设置样本观测点和 $\alpha$ 和 $\beta$ 的先验分布，定义似然函数进行MCMC推断

def invlogit(x):
    return pymc.exp(x)/(1+pymc.exp(x))

# 观测点
n = 5*np.ones(4)
x = np.array([-0.896, -0.296, -0.053, 0.727])
y_obs = np.array([0, 1, 3, 5])

# 定义先验
alpha = pymc.Normal('alpha', mu=0, tau=1.0/5**2)
beta = pymc.Normal('beta', mu=0, tau=1.0/10**2)

# 定义似然
p = pymc.InvLogit('p', alpha+beta*x)
y = pymc.Binomial('y_obs', n=n, p=p, value=y_obs, observed=True)

# 推断
m = pymc.Model([alpha, beta, y])
mc = pymc.MCMC(m)
mc.sample(iter=11000, burn=10000)

# beta参数的统计结果
print(beta.stats())

xp = np.linspace(-1, 1, 100)
a = alpha.stats()['mean']
b = beta.stats()['mean']
plt.plot(xp, invlogit(a+b*xp).value, linewidth=1, color='red')
plt.scatter(x, y_obs/5, s=50)
plt.xlabel('Log does of drug')
plt.ylabel('Risk of death')

plt.show()
pymc.Matplot.plot(mc)