R语言贝叶斯METROPOLIS-HASTINGS GIBBS 吉布斯采样器估计变点指数分布分析泊松过程车站等待时间...

本文主要是介绍R语言贝叶斯METROPOLIS-HASTINGS GIBBS 吉布斯采样器估计变点指数分布分析泊松过程车站等待时间...,希望对大家解决编程问题提供一定的参考价值,需要的开发者们随着小编来一起学习吧!

原文链接:http://tecdat.cn/?p=26578 

指数分布是泊松过程中事件之间时间的概率分布,因此它用于预测到下一个事件的等待时间,例如,您需要在公共汽车站等待的时间,直到下一班车到了点击文末“阅读原文”获取完整代码数据)。

相关视频

在本文中,我们将使用指数分布,假设它的参数 λ ,即事件之间的平均时间,在某个时间点 k 发生了变化,即:

1cffcf991100a4e82bd08311890131c4.png

我们的主要目标是使用 Gibbs 采样器在给定来自该分布的 n 个观测样本的情况下估计参数 λ、α 和 k。

吉布斯Gibbs 采样器

Gibbs 采样器是 Metropolis-Hastings 采样器的一个特例,通常在目标是多元分布时使用。使用这种方法,链是通过从目标分布的边缘分布中采样生成的,因此每个候选点都被接受。

Gibbs 采样器生成马尔可夫链如下:

  • 让 90e39b2c70b1edb9d8ae5db30ad3b4dd.png 是 Rd 中的随机向量,在时间 t=0 初始化 X(0)。

  • 对于每次迭代 t=1,2,3,...重复:

  • 设置 x1=X1(t-1)。

  • 对于每个 j=1,...,d:

  • 生成 X∗j(t) 从 d13f65a90c1184928941f6fbc7a2ca28.png, 其中 33b9b4b5a6b3146d79a1bc15f89da6d3.png 是给定 X(-j) 的 Xj的单变量条件密度。

  • 更新 a0dd6dbde59f4537bb93f99e0f465282.png.

  • 当每个候选点都被接受时,设置 a59c8d6c7047f211e02707f042a1ab8f.png.

  • 增加 t。

贝叶斯公式

变点问题的一个简单公式假设 f和 g 已知密度:

1601e9b57696db82b91b483e6c176e5c.png

其中 k 未知且 k=1,2,...,n。让 Yi为公交车到达公交车站之间经过的时间(以分钟为单位)。假设变化点发生在第 k分钟,即:

78a05e911039848bfb49091f250ad075.png

当 Y=(Y1,Y2,...,Yn) 时,似然 L(Y|k)由下式给出:

b87f3574979366ed5dc56085b27fe4af.png

假设具有独立先验的贝叶斯模型由下式给出:

656bd5e545ae437450fc34dc9d3e5049.png

数据和参数的联合分布为:

d7d810cecfe7f150b8e2a4b9ffafed75.png

其中,

b7633f06d43cdf3642e5a90d4b79ba9e.png

正如我之前提到的,Gibbs 采样器的实现需要从目标分布的边缘分布中采样,因此我们需要找到 λ、α 和 k 的完整条件分布。你怎么能这样做?简单来说,您必须从上面介绍的连接分布中选择仅依赖于感兴趣参数的项并忽略其余项。

相关视频

λ 的完整条件分布由下式给出:

919f9b902045e91ef249ee225cb22b5a.png

α 的完整条件分布由下式给出:

9479bfe7e8a0281e27d8037b982ad99b.png

k 的完整条件分布由下式给出:

7423e733db4b26b6ad9e0c2c5ef25bd6.png

计算方法

在这里,您将学习如何使用使用 R 的 Gibbs 采样器来估计参数 λ、α 和 k。

数据

首先,我们从具有变化点的下一个指数分布生成数据:

c632d5fb79838aeb168d90b9b2beabd6.png

set.seed(98712)
y <- c(rexp(25, rate = 2), rexp(35, rate = 10))

考虑到公交车站的情况,一开始公交车平均每2分钟一班,但从时间i=26开始,公交车开始平均每10分钟一班到公交车站。


点击标题查阅往期内容

80a884acb8261501485dd66467982cad.jpeg

R语言Gibbs抽样的贝叶斯简单线性回归仿真分析

outside_default.png

左右滑动查看更多

outside_default.png

01

3d9ca32719c644df60992e2e31e3de2a.png

02

1876f9f0d65700a0446114e2bb295ffe.png

03

226805c3fe01a6c2f2ee1be5dcf84cb7.png

04

d638d2dd6a2c0fb8db1aba5446877164.png

Gibbs采样器的实现

首先,我们需要初始化 k、λ 和 α。

n <- length(y) # 样本的观察值的数量
lci <- 10000 # 链的大小
aba <- alpha <- k <- numeric(lcan)
k\[1\] <- sample(1:n,

现在,对于算法的每次迭代,我们需要生成 λ(t)、α(t) 和 k(t),如下所示(记住如果 k+1>n 没有变化点):

b485bac07e6ecf4b180a6b1d532889b2.png

for (i in 2:lcan){kt <- k\[i-1\]# 生成lambdalambda\[i\] <- rgamma# 生成α# 产生k   for (j in 1:n) {L\[j\] <- ((lambda\[i\] / alpha\[i# 删除链条上的前9000个值
bunIn <- 9000

结果

在本节中,我们将介绍 Gibbs 采样器生成的链及其参数 λ、α 和 k 的分布。参数的真实值用红线表示。

f70db7276e50ac0089f359fbe6aa49c8.png

a5e2c46a564bf03106d5cdd3cd1dcf53.png

a615a0765702c2cf046eae7cc35a4812.png

下表显示了参数的实际值和使用 Gibbs 采样器获得的估计值的平均值:

res <- c(mean(k\[-(1:bun)\]), mean(lmba\[-(1:burn)\]), mean(apa\[-(1:buI)\]))
resfil

1b218e5b11c300914bd9744ab0204dc7.png

结论

从结果中,我们可以得出结论,使用 R 中的 Gibbs 采样器获得的具有变点的指数分布对参数 k、λ 和 α 的估计值的平均值接近于参数的实际值,但是我们期望更好估计。这可能是由于选择了链的初始值或选择了 λ 和 α的先验分布。


0d7d3d6314626062a58658de3a4d2671.png

点击文末“阅读原文”

获取全文完整资料。

本文选自《R语言贝叶斯METROPOLIS-HASTINGS GIBBS 吉布斯采样器估计变点指数分布分析泊松过程车站等待时间》。

点击标题查阅往期内容

R语言马尔可夫MCMC中的METROPOLIS HASTINGS,MH算法抽样(采样)法可视化实例

python贝叶斯随机过程:马尔可夫链Markov-Chain,MC和Metropolis-Hastings,MH采样算法可视化

Python贝叶斯推断Metropolis-Hastings(M-H)MCMC采样算法的实现

Metropolis Hastings采样和贝叶斯泊松回归Poisson模型

Matlab用BUGS马尔可夫区制转换Markov switching随机波动率模型、序列蒙特卡罗SMC、M H采样分析时间序列

R语言RSTAN MCMC:NUTS采样算法用LASSO 构建贝叶斯线性回归模型分析职业声望数据

R语言BUGS序列蒙特卡罗SMC、马尔可夫转换随机波动率SV模型、粒子滤波、Metropolis Hasting采样时间序列分析

R语言Metropolis Hastings采样和贝叶斯泊松回归Poisson模型

R语言贝叶斯MCMC:用rstan建立线性回归模型分析汽车数据和可视化诊断

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

R语言贝叶斯Poisson泊松-正态分布模型分析职业足球比赛进球数

R语言用Rcpp加速Metropolis-Hastings抽样估计贝叶斯逻辑回归模型的参数

R语言逻辑回归、Naive Bayes贝叶斯、决策树、随机森林算法预测心脏病

R语言中贝叶斯网络(BN)、动态贝叶斯网络、线性模型分析错颌畸形数据

R语言中的block Gibbs吉布斯采样贝叶斯多元线性回归

Python贝叶斯回归分析住房负担能力数据集

R语言实现贝叶斯分位数回归、lasso和自适应lasso贝叶斯分位数回归分析

Python用PyMC3实现贝叶斯线性回归模型

R语言用WinBUGS 软件对学术能力测验建立层次(分层)贝叶斯模型

R语言Gibbs抽样的贝叶斯简单线性回归仿真分析

R语言和STAN,JAGS:用RSTAN,RJAG建立贝叶斯多元线性回归预测选举数据

R语言基于copula的贝叶斯分层混合模型的诊断准确性研究

R语言贝叶斯线性回归和多元线性回归构建工资预测模型

R语言贝叶斯推断与MCMC:实现Metropolis-Hastings 采样算法示例

R语言stan进行基于贝叶斯推断的回归模型

R语言中RStan贝叶斯层次模型分析示例

R语言使用Metropolis-Hastings采样算法自适应贝叶斯估计与可视化

R语言随机搜索变量选择SSVS估计贝叶斯向量自回归(BVAR)模型

WinBUGS对多元随机波动率模型:贝叶斯估计与模型比较

R语言实现MCMC中的Metropolis–Hastings算法与吉布斯采样

R语言贝叶斯推断与MCMC:实现Metropolis-Hastings 采样算法示例

R语言使用Metropolis-Hastings采样算法自适应贝叶斯估计与可视化

视频:R语言中的Stan概率编程MCMC采样的贝叶斯模型

R语言MCMC:Metropolis-Hastings采样用于回归的贝叶斯估计

7dc750d88919ffffbea3b5e334d7e7e8.png

97cc385b6c9ba87ec02d9abd39628c51.jpeg

459e12c9de6eaa66e42c063c4a867876.png

这篇关于R语言贝叶斯METROPOLIS-HASTINGS GIBBS 吉布斯采样器估计变点指数分布分析泊松过程车站等待时间...的文章就介绍到这儿,希望我们推荐的文章对编程师们有所帮助!



http://www.chinasem.cn/article/710158

相关文章

最新版IDEA配置 Tomcat的详细过程

《最新版IDEA配置Tomcat的详细过程》本文介绍如何在IDEA中配置Tomcat服务器,并创建Web项目,首先检查Tomcat是否安装完成,然后在IDEA中创建Web项目并添加Web结构,接着,... 目录配置tomcat第一步,先给项目添加Web结构查看端口号配置tomcat    先检查自己的to

Redis主从复制实现原理分析

《Redis主从复制实现原理分析》Redis主从复制通过Sync和CommandPropagate阶段实现数据同步,2.8版本后引入Psync指令,根据复制偏移量进行全量或部分同步,优化了数据传输效率... 目录Redis主DodMIK从复制实现原理实现原理Psync: 2.8版本后总结Redis主从复制实

锐捷和腾达哪个好? 两个品牌路由器对比分析

《锐捷和腾达哪个好?两个品牌路由器对比分析》在选择路由器时,Tenda和锐捷都是备受关注的品牌,各自有独特的产品特点和市场定位,选择哪个品牌的路由器更合适,实际上取决于你的具体需求和使用场景,我们从... 在选购路由器时,锐捷和腾达都是市场上备受关注的品牌,但它们的定位和特点却有所不同。锐捷更偏向企业级和专

SpringBoot集成SOL链的详细过程

《SpringBoot集成SOL链的详细过程》Solanaj是一个用于与Solana区块链交互的Java库,它为Java开发者提供了一套功能丰富的API,使得在Java环境中可以轻松构建与Solana... 目录一、什么是solanaj?二、Pom依赖三、主要类3.1 RpcClient3.2 Public

Android数据库Room的实际使用过程总结

《Android数据库Room的实际使用过程总结》这篇文章主要给大家介绍了关于Android数据库Room的实际使用过程,详细介绍了如何创建实体类、数据访问对象(DAO)和数据库抽象类,需要的朋友可以... 目录前言一、Room的基本使用1.项目配置2.创建实体类(Entity)3.创建数据访问对象(DAO

C语言线程池的常见实现方式详解

《C语言线程池的常见实现方式详解》本文介绍了如何使用C语言实现一个基本的线程池,线程池的实现包括工作线程、任务队列、任务调度、线程池的初始化、任务添加、销毁等步骤,感兴趣的朋友跟随小编一起看看吧... 目录1. 线程池的基本结构2. 线程池的实现步骤3. 线程池的核心数据结构4. 线程池的详细实现4.1 初

SpringBoot整合kaptcha验证码过程(复制粘贴即可用)

《SpringBoot整合kaptcha验证码过程(复制粘贴即可用)》本文介绍了如何在SpringBoot项目中整合Kaptcha验证码实现,通过配置和编写相应的Controller、工具类以及前端页... 目录SpringBoot整合kaptcha验证码程序目录参考有两种方式在springboot中使用k

Spring中Bean有关NullPointerException异常的原因分析

《Spring中Bean有关NullPointerException异常的原因分析》在Spring中使用@Autowired注解注入的bean不能在静态上下文中访问,否则会导致NullPointerE... 目录Spring中Bean有关NullPointerException异常的原因问题描述解决方案总结

python中的与时间相关的模块应用场景分析

《python中的与时间相关的模块应用场景分析》本文介绍了Python中与时间相关的几个重要模块:`time`、`datetime`、`calendar`、`timeit`、`pytz`和`dateu... 目录1. time 模块2. datetime 模块3. calendar 模块4. timeit

python-nmap实现python利用nmap进行扫描分析

《python-nmap实现python利用nmap进行扫描分析》Nmap是一个非常用的网络/端口扫描工具,如果想将nmap集成进你的工具里,可以使用python-nmap这个python库,它提供了... 目录前言python-nmap的基本使用PortScanner扫描PortScannerAsync异