如何在matlab中使用metropolis-hasting算法

如何在matlab中使用metropolis-hasting算法,第1张

MH算法也是一种基于模拟的MCMC技术,一个很重要的应用是从给定的概率分布中抽样。主要原理是构造了一个精妙的Markov链,使得该链的稳态 是你给定的概率密度。它的好处,不用多说,自然是可以对付数学形式复杂的概率密度。有人说,单维的MH算法配上Gibbs Sampler几乎是“无敌”了。

今天试验的过程中发现,MH算法想用好也还不简单,里面的转移参数设定就不是很好弄。即使用最简单的高斯漂移项,方差的确定也是个头疼的问题,需要不同问题不同对待,多试验几次。当然你也可以始终选择“理想”参数。

还是拿上次的混合高斯分布来做模拟,模拟次数为500000次的时候,概率分布逼近的程度如下图。虽然几个明显的"峰"已经出来了,但是数值上还是 有很大差异的。估计是我的漂移方差没有选好。感觉还是inverse sampling好用,迭代次数不用很多,就可以达到相当的逼近程度。

试了一下MH算法,

R Code:

p=function(x,u1,sig1,u2,sig2){

(1/3)(1/(sqrt(2pi)15)exp(-05(x-70)^2/15^2)+1/(sqrt(2pi)11)exp(-05(x+80)^2/11^2)+1/(sqrt(2pi)sig1)exp(-05(x-u1)^2/sig1^2)+1/(sqrt(2pi)sig2)exp(-05(x-u2)^2/sig2^2))

}

MH=function(x0,n){

x=NULL

x[1] = x0

for (i in 1:n){

x_can= x[i]+rnorm(1,0,325)

d= p(x_can,10,30,-10,10)/p(x[i],10,30,-10,10)

alpha= min(1,d)

u=runif(1,0,1)

if (u<alpha){

x[i+1]=x_can}

else{

x[i+1]=x[i]

}

if (round(i/100)==i/100) print(i)

}

x

}

z=MH(10,99999)

z=z[-10000]

a=seq(-100,100,02)

plot(density(z),col=1,main='Estimated Density',ylim=c(0,002),lty=1)

points(a, p(a,10,30,-10,10),pch='',col=2,lty=2)

legend(60,002,c("True","Sim (MH)"),col=c(1,2),lty=c(1,2))

A = [

1 002 004 002 002;

004 1 008 004 004;

003846 003846 1 003846 003846;

008 008 016 1 008;

005 005 01 005 1

];

b = [403241; 29156763; 292425; 2217813; 2931695];

x=A\b

得到结果:

x =

3813195

2421583

2552374

1125977

2308420

matlab高斯拟合的初始值问题:最近数学建模用到了高斯混合模型GMM,高斯分布是一种在自然界大量的存在的、最为常见的分布形式。我们的目标用高斯混合模型对数据集进行拟合聚类,得到一个概率的分布,最后做出来是很好的效果。

fs = 100; %采样频率

fc = 10; %正弦波中心频率

t = 0:1/fs:1; %时间

x = sin(2pifct); %正弦波

y = awgn(x,10); %加入awgn噪声,10指的是信噪比为10dB

以上就是关于如何在matlab中使用metropolis-hasting算法全部的内容,包括:如何在matlab中使用metropolis-hasting算法、急求,帮忙用matlab或c/c++高斯迭代法程序求x1,x2,x3,x4,x5、matlab高斯拟合的初始值问题等相关内容解答,如果想了解更多相关内容,可以关注我们,你们的支持是我们更新的动力!

欢迎分享,转载请注明来源:内存溢出

原文地址: http://outofmemory.cn/zz/10635589.html

(0)
打赏 微信扫一扫 微信扫一扫 支付宝扫一扫 支付宝扫一扫
上一篇 2023-05-10
下一篇 2023-05-10

发表评论

登录后才能评论

评论列表(0条)

保存