马尔可夫链蒙特卡罗法

Posted 2020-12-21 botak

马尔可夫链蒙特卡罗法

蒙特卡罗法

思想：假设概率分布的定义已知，然后通过随机抽样获得概率分布的随机样本，通过得到的随机样本对概率分布的特征进行分析。
for example：从随机抽出的样本中计算出样本均值，从而得到总体的期望。
蒙特卡罗方法的核心：随机抽样
主要有直接抽样，接受-拒绝抽样，重要性抽样

随机抽样

接受拒绝抽样
input：抽样的目标概率分布的概率密度函数(p(x))
output：概率分布的随机样本(x_1,x_2,...,x_n)
parameters：样本数n
建议分布：(q(x))，概率分布：(p(x))

[u<= frac{p(x^* )}{cq(x^ * )} ]

数学期望估计

样本均值：(hat{f_n})

[hat{f_n} = frac{1}{n}sum_{i=1}^nf(x_i) ]

根据大数定律可知，当样本容量增大时，样本均值以概率1收敛于数学期望：

[hat{f_n} -> E_{p(x)}[f(x)]~,~n->infty ]

[E_{p(x)}[f(x)] = frac{1}{n}sum_{i=1}^nf(x_i) ]

积分计e算

[int_{mathcal{X}}h(x)dx ]

将(h(x))分解成为一个函数(f(x))和一个密度函数(p(x))

[int_{mathcal{X}}h(x)dx = int_{mathcal{X}}f(x)p(x)dx = E_{p(x)}[f(x)] ]

summery

一般的蒙特卡罗法中的抽样分布是独立的，而马尔可夫链蒙特卡罗法中的抽样样本不是独立的，样本序列形成马尔可夫链

马尔可夫链

基本定义

马尔可夫性：

[P(X_t|X_0,X_1,X_2,...,X_{t-1}) = P(X_t|X_{t-1})~,~t=1,2,...,n ]

马尔可夫链（具有马尔可夫性的随机序列）也称之为马尔可夫过程：

[X = lbrace X_0,X_1,X_2,...,X_t,... brace ]

马尔可夫链的转移概率分布(决定了马尔可夫链的特性)：

[P(X_t|X_{t-1})~,~t=1,2,3,...,n ]

时间齐次的马尔可夫链（转移概率分布(P(X_t|X_{t-1}))与t无关）：

[P(X_{t+s}|X_{t-1+s}) = P(X_t|X_{t-1})~,~t=1,2,3,...,n ]

另外，还有n阶马尔可夫链，n阶马尔可夫链可以转化成为1阶马尔可夫链。

离散状态马尔可夫链

• 转移概率矩阵和状态分布
  (p_{ij} = (X_t = i | X_{t-1} = j)~,~i=1,2,...)满足(p_{ij}>=0,sum_ip_{ij}=1),且满足这两个条件的矩阵称之为随机矩阵,马尔可夫链(X=lbrace X_0,X_1,...,X_t,... brace)在时刻t的概率分布称之为t的状态分布，其中，(pi_i(t))指的是时刻t状态为i的概率(P(X_t = i))记做,：
  [pi(t) = [pi_1(t),...,pi_n(t)]^T ]

• 平稳分布（P：状态转移矩阵(P = (p_{ij}))）
  [pi = Ppi ]

直观上，如果马尔可夫链的平稳分布存在，那么以该平稳分布作为初始分布，而向未来随机状态转移，之后的任何一个状态也是平稳状态。
平稳分布的充分必要条件：

[x_i = sum_{j}p_{ij}x_j ]

[x_i>=0 ]

[sum_ix_i=1 ]

马尔可夫链X在时刻t的状态分布，可以由在时刻（t-1）的状态分布，以及转移概率决定：

[pi(t) = Ppi(t-1) ]

(pi(t) = P^tpi(0)),这里的(P^t)称为t步转移概率矩阵。

连续状态马尔可夫链

转移核(P(x,A))定义为：

[P(x,A) = int_Ap(x,y)dy ]

[P(X_t = A|X_t-1 = x) = P(x,A) ]

马尔可夫链的额性质

• 不可约
  (P(X_t = i|X_0=j)>0)
  从任意状态出发，当经过充分长的时间后，可以达到任意一个状态。
• 非周期
  ({t:P(X_t=i|X_)=i)>0} ext{的最大公约数为1})
  一个非周期的马尔可夫链，不存在一个状态，从这一个状态出发，再返回这个状态所经历的时间呈现一定的周期性。
• 正常返
  (lim_{t->+infty}p_{ij}^t>0)
  一个正常返的马尔可夫链，其中任意一个状态，从其他任何一个状态出发，时间趋于无穷时，首次转移到这个状态的概率不为0
• 遍历定理
  不可约的非周期的正常返的马尔可夫链，当时间趋于无穷时，马尔可夫链的状态总会趋于平稳分布
• 可逆马尔可夫链
  (p_{ij}pi_j = p_{ji}pi_i)
  可逆的马尔可夫链，无论是面向未来还是面向过去，任意一个时刻的状态分布均是平稳分布。
• 细致平衡方程
  (Ppi = pi)

马尔可夫链蒙特卡罗法

基本想法

燃烧期(0~m)
(hat{f(x)} = frac{1}{n-m}sum_{i=n-m+1}^nf(x_i))

基本步骤

step 1:在随机变量(x)的状态空间里构造一个满足条件的马尔可夫链，使得其平稳分布为目标分布(p(x))
step 2:从状态空间的某一点(x_0)出发，用构造的马尔可夫链进行随机游走，产生样本序列(x_0,x_1,x_2,..)
step 3:应用马尔可夫链的遍历定理，确定正整数m和n，得到样本集合({x_{m+1},...,x_n})得到函数的均值(hat{f(x)} = frac{1}{n-m}sum_{i=n-m+1}^nf(x_i))
important questions:
one:如何定义马尔可夫链，保证马尔可夫链蒙特卡罗法的条件成立
two:如何确定m，保证样本的无偏性
three:如何确定n，保证遍历均值的精度

Metropolis-Hastings算法

基本原理

• 马尔可夫链
  Metropolis-Hastings算法采用转移核为(p(x,x^{prime}))
  [p(x,x^{prime}) = q(x,x^{prime})alpha(x,x^{prime}) ]

• 建议分布(q(x,x^{prime}))和接受分布(alpha(x,x^{prime}))
  [alpha(x,x^{prime}) = minlbrace1, frac{p(x^{prime}q(x^{prime},x))}{p(x)q(x,x^{prime})} brace ]
  转移核：
  [p(x,x^{prime}) = ?? ]
  构造出来了：
  [p(x)p(x,x^{prime}) = p(x^{prime})p(x^{prime},x) ]

• 满条件分布
  []

Metropolis-Hastings算法

input:抽样的目标分布的密度函数(p(x)),函数(f(x))
output:(p(x))的随机样本(x_{m+1},..,x_n)，函数样本的均值
parameters:m，n

单分量Metropolis-Hastings法

Metropolis-Hastings算法需要对多变量分布进行抽样，有时候对多元变量分布的抽样是苦难的，可以对多元变量的每一个变量的条件分布依次进行抽样，这就是单分量Metropolis-Hastings法

吉布斯抽样

基本做法是，从联合概率分布定义满条件概率分布，依次对满条件分布进行抽样，得到样本的序列。