简单易学的机器学习算法——马尔可夫链蒙特卡罗方法MCMC
来源:互联网 发布:杨颖开的淘宝店铺名字 编辑:程序博客网 时间:2024/06/05 10:39
- 一马尔可夫链
- 马尔可夫链
- 转移概率
- 马尔可夫链的平稳分布
- 二马尔可夫链蒙特卡罗方法
- 基本思想
- 细致平稳条件
- Metropolis采样算法
- 1Metropolis采样算法的基本原理
- 2Metropolis采样算法的流程
- 3Metropolis算法的解释
- 4实验
- 参考文献
对于一般的分布的采样,在很多的编程语言中都有实现,如最基本的满足均匀分布的随机数,但是对于复杂的分布,要想对其采样,却没有实现好的函数,在这里,可以使用马尔可夫链蒙特卡罗(Markov Chain Monte Carlo, MCMC)方法,其中Metropolis-Hastings采样和Gibbs采样是MCMC中使用较为广泛的两种形式。
MCMC的基础理论为马尔可夫过程,在MCMC算法中,为了在一个指定的分布上采样,根据马尔可夫过程,首先从任一状态出发,模拟马尔可夫过程,不断进行状态转移,最终收敛到平稳分布。
一、马尔可夫链
1、马尔可夫链
设
也就是说状态转移的概率只依赖于前一个状态。称这个变量为马尔可夫变量,其中,
马尔可夫链指的是在一段时间内随机变量
2、转移概率
马尔可夫链是通过对应的转移概率定义的,转移概率指的是随机变量从一个时刻到下一个时刻,从状态
记
假设状态的数目为
3、马尔可夫链的平稳分布
对于马尔可夫链,需要注意以下的两点:
- 1、周期性:即经过有限次的状态转移,又回到了自身;
- 2、不可约:即两个状态之间相互转移;
如果一个马尔可夫过程既没有周期性,又不可约,则称为各态遍历的。
对于一个各态遍历的马尔可夫过程,无论初始值
且这个平稳分布
其中,
二、马尔可夫链蒙特卡罗方法
1、基本思想
对于一个给定的概率分布
2、细致平稳条件
对于一个各态遍历的马尔可夫过程,若其转移矩阵为
则
3、Metropolis采样算法
Metropolis采样算法是最基本的基于MCMC的采样算法。
3.1、Metropolis采样算法的基本原理
假设需要从目标概率密度函数
其中,
在Metropolis采样算法的过程中,首先初始化状态值
这样的过程一直持续到采样过程的收敛,当收敛以后,样本
3.2、Metropolis采样算法的流程
基于以上的分析,可以总结出如下的Metropolis采样算法的流程:
- 初始化时间
t=1 - 设置
u 的值,并初始化初始状态θ(t)=u - 重复一下的过程:
- 令
t=t+1 - 从已知分布
q(θ∣θ(t−1)) 中生成一个候选状态θ(∗) - 计算接受的概率:
α=min(1,p(θ(∗))p(θ(t−1))) - 从均匀分布
Uniform(0,1) 生成一个随机值a - 如果
a⩽α ,接受新生成的值:θ(t)=θ(∗) ;否则:θ(t)=θ(t−1)
- 令
- 直到
t=T
3.3、Metropolis算法的解释
要证明Metropolis采样算法的正确性,最重要的是要证明构造的马尔可夫过程满足如上的细致平稳条件,即:
对于上面所述的过程,分布为
其中,
对于选择该已知的分布,在Metropolis采样算法中,要求该已知的分布必须是对称的,即
Qi,j=Qj,i ,即q(θ=θ(t)∣θ(t−1))=q(θ=θ(t−1)∣θ(t))
常用的符合对称的分布主要有:正态分布,柯西分布以及均匀分布等。
接下来,需要证明在Metropolis采样算法中构造的马尔可夫链满足细致平稳条件。
因此,通过以上的方法构造出来的马尔可夫链是满足细致平稳条件的。
3.4、实验
假设需要从柯西分布中采样数据,我们利用Metropolis采样算法来生成样本,其中,柯西分布的概率密度函数为:
那么,根据上述的Metropolis采样算法的流程,接受概率
代码如下:
'''Date:20160629@author: zhaozhiyong'''import randomfrom scipy.stats import normimport matplotlib.pyplot as pltdef cauchy(theta): y = 1.0 / (1.0 + theta ** 2) return yT = 5000sigma = 1thetamin = -30thetamax = 30theta = [0.0] * (T+1)theta[0] = random.uniform(thetamin, thetamax)t = 0while t < T: t = t + 1 theta_star = norm.rvs(loc=theta[t - 1], scale=sigma, size=1, random_state=None) #print theta_star alpha = min(1, (cauchy(theta_star[0]) / cauchy(theta[t - 1]))) u = random.uniform(0, 1) if u <= alpha: theta[t] = theta_star[0] else: theta[t] = theta[t - 1]ax1 = plt.subplot(211)ax2 = plt.subplot(212) plt.sca(ax1)plt.ylim(thetamin, thetamax)plt.plot(range(T+1), theta, 'g-')plt.sca(ax2)num_bins = 50plt.hist(theta, num_bins, normed=1, facecolor='red', alpha=0.5)plt.show()
- 1
- 2
- 3
- 4
- 5
- 6
- 7
- 8
- 9
- 10
- 11
- 12
- 13
- 14
- 15
- 16
- 17
- 18
- 19
- 20
- 21
- 22
- 23
- 24
- 25
- 26
- 27
- 28
- 29
- 30
- 31
- 32
- 33
- 34
- 35
- 36
- 37
- 38
- 39
- 40
- 41
- 1
- 2
- 3
- 4
- 5
- 6
- 7
- 8
- 9
- 10
- 11
- 12
- 13
- 14
- 15
- 16
- 17
- 18
- 19
- 20
- 21
- 22
- 23
- 24
- 25
- 26
- 27
- 28
- 29
- 30
- 31
- 32
- 33
- 34
- 35
- 36
- 37
- 38
- 39
- 40
- 41
实验的结果:
对于Metropolis采样算法,其要求选定的分布必须是对称的,为了弥补这样的一个缺陷,在下一篇中,介绍一下Metropolis-Hastings采样算法,其是Metropolis采样算法的推广形式。
参考文献
- 1、马尔可夫链蒙特卡罗算法
- 2、受限玻尔兹曼机(RBM)学习笔记(一)预备知识
- 3、LDA数学八卦
- 简单易学的机器学习算法——马尔可夫链蒙特卡罗方法MCMC
- 简单易学的机器学习算法——马尔可夫链蒙特卡罗方法MCMC
- 简单易学的机器学习算法——马尔可夫链蒙特卡罗方法MCMC
- 简单易学的机器学习算法——Apriori算法
- 简单易学的机器学习算法——Apriori算法
- 简单易学的机器学习算法——EM算法
- 简单易学的机器学习算法——朴素贝叶斯
- 简单易学的机器学习算法——Logistic回归
- 简单易学的机器学习算法——lasso
- 简单易学的机器学习算法——Softmax Regression
- 简单易学的机器学习算法——lasso
- 简单易学的机器学习算法——kMeans
- 简单易学的机器学习算法——AdaBoost
- 简单易学的机器学习算法——Gibbs采样
- 简单易学的机器学习算法——Softmax Regression
- 简单易学的机器学习算法——Softmax Regression
- 简单易学的机器学习算法——Softmax Regression
- 简单易学的机器学习算法——朴素贝叶斯
- SpringMVC学习笔记(1)
- 最大子序列和问题
- [Codeforces Gym 101243][2016acm/icpc Quarterfinal Central region of Russia]
- 请求删除任务,OSTaskDelReq()
- rocketmq
- 简单易学的机器学习算法——马尔可夫链蒙特卡罗方法MCMC
- 人脸识别人证比对技术SDK
- “听说你是个程序员,不能把时间消耗在阅读微信公众号上”
- 阿里云linux服务器上使用iptables设置安全策略的方法
- 自定义组合控件流程
- 0/1分数规划 【POJ2976】Dropping tests
- 智慧城市管理平台案例分享
- Java设计模式----1.代理模式
- POJ 2664 Prerequisites? G++