R语言贝叶斯MCMC:GLM逻辑回归、Rstan线性回归、Metropolis Hastings与Gibbs采样算法实例

简介: R语言贝叶斯MCMC:GLM逻辑回归、Rstan线性回归、Metropolis Hastings与Gibbs采样算法实例

什么是频率学派?

在频率学派中,观察样本是随机的,而参数是固定的、未知的数量。

概率被解释为一个随机过程的许多观测的预期频率。

有一种想法是 "真实的",例如,在预测鱼的生活环境时,盐度和温度之间的相互作用有一个回归系数?

什么是贝叶斯学派?

在贝叶斯方法中,概率被解释为对信念的主观衡量。

所有的变量--因变量、参数和假设都是随机变量。我们用数据来确定一个估计的确定性(可信度)。

这种盐度X温度的相互作用反映的不是绝对的,而是我们对鱼的生活环境所了解的东西(本质上是草率的)。

目标

频率学派

保证正确的误差概率,同时考虑到抽样、样本大小和模型。

  • 缺点:需要对置信区间、第一类和第二类错误进行复杂的解释。
  • 优点:更具有内在的 "客观性 "和逻辑上的一致性。

贝叶斯学派

分析更多的信息能在多大程度上提高我们对一个系统的认识。

  • 缺点:这都是关于信仰的问题! ...有重大影响。
  • 优点: 更直观的解释和实施,例如,这是这个假设的概率,这是这个参数等于这个值的概率。可能更接近于人类自然地解释世界的方式。

实际应用中:为什么用贝叶斯

  • 具有有限数据的复杂模型,例如层次模型,其中

  • 实际的先验知识非常少

贝叶斯法则:

一些典型的贝叶斯速记法。


注意:

  • 贝叶斯的最大问题在于确定先验分布。先验应该是什么?它有什么影响?

目标:

计算参数的后验分布:π(θ|X)。

点估计是后验的平均值。

一个可信的区间是

你可以把它解释为一个参数在这个区间内的概率 。

计算

皮埃尔-西蒙-拉普拉斯(1749-1827)(见:Sharon Bertsch McGrayne: The Theory That Would Not Die)


  • 有些问题是可分析的,例如二项式似然-贝塔先验。
  • 但如果你有很多参数,这是不可能完成的操作
  • 如果你有几个参数,而且是奇数分布,你可以用数值乘以/整合先验和似然(又称网格近似)。
  • 尽管该理论可以追溯到1700年,甚至它对推理的解释也可以追溯到19世纪初,但它一直难以更广泛地实施,直到马尔科夫链蒙特卡洛技术的发展。

MCMC

MCMC的思想是对参数值θi进行 "抽样"。

回顾一下,马尔科夫链是一个随机过程,它只取决于它的前一个状态,而且(如果是遍历的),会生成一个平稳的分布。

技巧 "是找到渐进地接近正确分布的抽样规则(MCMC算法)。

有几种这样的(相关)算法。

  • Metropolis-Hastings抽样
  • Gibbs 抽样
  • No U-Turn Sampling (NUTS)
  • Reversible Jump

一个不断发展的文献和工作体系!

Metropolis-Hastings 算法

  1. 开始:
  2. 跳到一个新的候选位置:
  3. 计算后验:
  4. 如果
  5. 如果
  6. 转到第2步

Metropolis-Hastings: 硬币例子

你抛出了5个正面。你对θ的最初 "猜测 "是

MCMC:

p.old <- prior *likelihood 
while(length(thetas) <= n){
  theta.new <- theta + rnorm(1,0,0.05)
  p.new <- prior *likelihood 
  if(p.new > p.old | runif(1) < p.new/p.old){
    theta <- theta.new
    p.old <- p.new
  }

画图:

hist(thetas\[-(1:100)\] )
curve(6*x^5 )

采样链:调整、细化、多链

  • 那个 "朝向 "平稳的初始过渡被称为 "预烧期",必须加以修整。
  • 怎么做?用眼睛看
  • 采样过程(显然)是自相关的。
  • 如何做?通常是用眼看,用acf()作为指导。
  • 为了保证你收敛到正确的分布,你通常会从不同的位置获得多条链(例如4条)。
  • 有效样本量

MCMC 诊断法

R软件包帮助分析MCMC链。一个例子是线性回归的贝叶斯拟合(α,β,σ

plot(line)


预烧部分:

plot(line\[\[1\]\], start=10)

MCMC诊断法

查看后验分布(同时评估收敛性)。

density(line)


参数之间的关联性,以及链内的自相关关系

levelplot(line\[\[2\]\])
acfplot(line)

统计摘要

运行MCMC的工具(在R内部)

逻辑Logistic回归:婴儿出生体重低

logitmcmc(low~age+as.factor(race)+smoke )


plot(mcmc)

MCMC与GLM逻辑回归的比较

MCMC与GLM逻辑回归的比较

对于这个应用,没有很好的理由使用贝叶斯建模,除非--你是 "贝叶斯主义者"。你有关于回归系数的真正先验信息(这基本上是不太可能的)。

一个主要的缺点是 先验分布棘手的调整参数。

但是,MCMC可以拟合的一些更复杂的模型(例如,层次的logit MCMChlogit)。

Metropolis-Hastings

Metropolis-Hastings很好,很简单,很普遍。但是对循环次数很敏感。而且可能太慢,因为它最终会拒绝大量的循环。

Gibbs 采样


在Gibbs吉布斯抽样中,你不是用适当的概率接受/拒绝,而是用适当的条件概率在参数空间中行进。并从该分布中抽取一次。

然后你从新的条件分布中抽取下一个参数。

比Metropolis-Hastings快得多。有效样本量要高得多!

BUGS(OpenBUGS,WinBUGS)是使用吉布斯采样器的贝叶斯推理。

JAGS是 "吉布斯采样器"

其他采样器

汉密尔顿蒙特卡洛(HMC)--是一种梯度的Metropolis-Hastings,因此速度更快,对参数之间的关联性更好。

No-U Turn Sampler(NUTS)--由于不需要固定的长度,它的速度更快。这是STAN使用的方法(见http://arxiv.org/pdf/1111.4246v1.pdf)。


(Hoffman and Gelman 2011)

其他工具

你可能想创建你自己的模型,使用贝叶斯MC进行拟合,而不是依赖现有的模型。为此,有几个工具可以选择。

  • BUGS / WinBUGS / OpenBUGS (Bayesian inference Using Gibbs Sampling) - 贝叶斯抽样工具的鼻祖(自1989年起)。WinBUGS是专有的。OpenBUGS的支持率很低。
  • JAGS(Just Another Gibbs Sampler)接受一个用类似于R语言的语法编写的模型字符串,并使用吉布斯抽样从这个模型中编译和生成MCMC样本。可以在R中使用rjags包。
  • Stan(以Stanislaw Ulam命名)是一个类似于JAGS的相当新的程序--速度更快,更强大,发展迅速。从伪R/C语法生成C++代码。安装:http://mc-stan.org/rstan.html**
  • Laplace’s Demon 所有的贝叶斯工具都在R中:http://www.bayesian-inference.com/software

STAN


要用STAN拟合一个模型,步骤是:

  1. 为模型生成一个STAN语法伪代码(在JAGS和BUGS中相同
  2. 运行一个R命令,用C++语言编译该模型
  3. 使用生成的函数来拟合你的数据

STAN示例--线性回归

STAN代码是R(例如,具有分布函数)和C(即你必须声明你的变量)之间的一种混合。每个模型定义都有三个块。

_1_.数据块:

int n; //
  vector\[n\] y; // Y 向量

这指定了你要输入的原始数据。在本例中,只有Y和X,它们都是长度为n的(数字)向量,是一个不能小于0的整数。

_2_. 参数块

real beta1;  // slope

这些列出了你要估计的参数:截距、斜率和方差。


_3_. 模型块

sigma ~ inv_gamma(0.001, 0.001); 
    yhat\[i\] <- beta0 + beta1 * (x\[i\] - mean(x));}
    y ~ normal(yhat, sigma);

注意:

  • 你可以矢量化,但循环也同样快
  • 有许多分布(和 "平均值 "等函数)可用
请经常参阅手册! https://github.com/stan-dev/stan/releases/download/v2.9.0/stan-reference-2.9.0.pdf

2. 在R中编译模型

你把你的模型保存在一个单独的文件中, 然后用stan_model()命令编译这个模型。

这个命令是把你描述的模型,用C++编码和编译一个NUTS采样器。相信我,自己编写C++代码是一件非常非常痛苦的事情(如果没有很多经验的话),而且它保证比R中的同等代码快得多。

注意:这一步可能会很慢。

3. 在R中运行该模型

这里的关键函数是sampling()。还要注意的是,为了给你的模型提供数据,它必须是列表的形式

模拟一些数据。

X <- runif(100,0,20)
Y <- rnorm(100, beta0+beta1*X, sigma)

进行取样!

sampling(stan, Data)

这里有大量的输出,因为它计算了


print(fit, digits = 2)

MCMC诊断法

为了应用coda系列的诊断工具,你需要从STAN拟合对象中提取链,并将其重新创建为mcmc.list。

extract(stan.fit
alply(chains, 2, mcmc)


相关文章
|
4月前
|
机器学习/深度学习 数据采集 算法
Python实现Naive Bayes贝叶斯分类模型(GaussianNB、MultinomialNB算法)项目实战
Python实现Naive Bayes贝叶斯分类模型(GaussianNB、MultinomialNB算法)项目实战
104 1
|
22天前
|
机器学习/深度学习 算法 数据安全/隐私保护
基于贝叶斯优化CNN-LSTM网络的数据分类识别算法matlab仿真
本项目展示了基于贝叶斯优化(BO)的CNN-LSTM网络在数据分类中的应用。通过MATLAB 2022a实现,优化前后效果对比明显。核心代码附带中文注释和操作视频,涵盖BO、CNN、LSTM理论,特别是BO优化CNN-LSTM网络的batchsize和学习率,显著提升模型性能。
|
27天前
|
机器学习/深度学习 算法 数据安全/隐私保护
基于贝叶斯优化卷积神经网络(Bayes-CNN)的多因子数据分类识别算法matlab仿真
本项目展示了贝叶斯优化在CNN中的应用,包括优化过程、训练与识别效果对比,以及标准CNN的识别结果。使用Matlab2022a开发,提供完整代码及视频教程。贝叶斯优化通过构建代理模型指导超参数优化,显著提升模型性能,适用于复杂数据分类任务。
|
2月前
|
存储 算法 测试技术
预见未来?Python线性回归算法:数据中的秘密预言家
【9月更文挑战第11天】在数据的海洋中,线性回归算法犹如智慧的预言家,助我们揭示未知。本案例通过收集房屋面积、距市中心距离等数据,利用Python的pandas和scikit-learn库构建房价预测模型。经过训练与测试,模型展现出较好的预测能力,均方根误差(RMSE)低,帮助房地产投资者做出更明智决策。尽管现实关系复杂多变,线性回归仍提供了有效工具,引领我们在数据世界中自信前行。
50 5
|
3月前
|
机器学习/深度学习 人工智能 算法
【人工智能】线性回归模型:数据结构、算法详解与人工智能应用,附代码实现
线性回归是一种预测性建模技术,它研究的是因变量(目标)和自变量(特征)之间的关系。这种关系可以表示为一个线性方程,其中因变量是自变量的线性组合。
70 2
|
3月前
|
存储 算法 定位技术
预见未来?Python线性回归算法:数据中的秘密预言家
【8月更文挑战第3天】站在数据的海洋边,线性回归算法犹如智慧的预言家,揭示着房价的秘密。作为房地产投资者,面对复杂的市场,我们可通过收集房屋面积、位置等数据并利用Python的pandas及scikit-learn库,建立线性回归模型预测房价。通过评估模型的均方根误差(RMSE),我们可以更精准地判断投资时机,让数据引领我们走向成功的彼岸。
24 1
|
3月前
|
机器学习/深度学习 算法 数据可视化
Python数据分析高手修炼手册:线性回归算法,让你的数据说话更有力
【8月更文挑战第1天】在数据驱动时代,掌握数据分析技能至关重要。线性回归是最基础且强大的工具之一,能从复杂数据中提炼简单有效的模型。本文探索Python中线性回归的应用并通过实战示例加深理解。线性回归建立变量间线性关系模型:Y = β0 + β1*X + ε。使用scikit-learn库进行实战:首先安装必要库,然后加载数据、训练模型并评估性能。示例展示了如何使用`LinearRegression`模型进行房价预测,包括数据可视化。掌握线性回归,让数据“说话”更有力。
42 2
|
4月前
|
数据采集 机器学习/深度学习 算法
Python实现多元线性回归模型(statsmodels OLS算法)项目实战
Python实现多元线性回归模型(statsmodels OLS算法)项目实战
652 2
|
4月前
|
机器学习/深度学习 数据采集 算法
Python实现贝叶斯岭回归模型(BayesianRidge算法)并使用K折交叉验证进行模型评估项目实战
Python实现贝叶斯岭回归模型(BayesianRidge算法)并使用K折交叉验证进行模型评估项目实战
130 1
|
4月前
|
机器学习/深度学习 人工智能 分布式计算
算法金 | 最难的来了:超参数网格搜索、贝叶斯优化、遗传算法、模型特异化、Hyperopt、Optuna、多目标优化、异步并行优化
机器学习中的超参数调优是提升模型性能的关键步骤,包括网格搜索、随机搜索、贝叶斯优化和遗传算法等方法。网格搜索通过穷举所有可能的超参数组合找到最优,但计算成本高;随机搜索则在预设范围内随机采样,降低计算成本;贝叶斯优化使用代理模型智能选择超参数,效率高且适应性强;遗传算法模拟生物进化,全局搜索能力强。此外,还有多目标优化、异步并行优化等高级技术,以及Hyperopt、Optuna等优化库来提升调优效率。实践中,应结合模型类型、数据规模和计算资源选择合适的调优策略。
175 0
算法金 | 最难的来了:超参数网格搜索、贝叶斯优化、遗传算法、模型特异化、Hyperopt、Optuna、多目标优化、异步并行优化