當前位置：首頁 > 编程资源 > 编程问答 >内容正文

编程问答

martingale、markov chain、Monte Carlo、MCMC

發布時間：2023/12/2 编程问答 36 豆豆

生活随笔收集整理的這篇文章主要介紹了 martingale、markov chain、Monte Carlo、MCMC 小編覺得挺不錯的,現在分享給大家,幫大家做個參考.

文章結構如下：

1: MCMC

1.1 MCMC是什么
1.2 為什么需要MCMC
2：蒙特卡羅

2.1 引入
2.2 均勻分布，Box-Muller 變換
2.3 拒絕接受采樣（Acceptance-Rejection Sampling）
2.4 接受拒絕采樣的直觀解釋
2.5 接受拒絕采樣方法有效性證明
2.6 接受拒絕采樣方法python實現
2.7 蒙特卡羅方法小結
3：馬爾科夫鏈

3.1 馬爾科夫鏈概述
3.2 馬爾科夫鏈模型狀態轉移矩陣的性質
3.3 基于馬爾科夫鏈采樣
3.4 馬爾科夫鏈采樣小結
3.5 補充：PageRank：Google的民主表決式網頁排名技術
4：馬爾可夫鏈蒙特卡羅算法

4.1 馬爾科夫鏈的細致平穩條件(Detailed Balance Condition)
4.2 MCMC采樣
5： M-H采樣

5.1 M-H采樣算法
5.2 M-H采樣python實現
5.3 M-H采樣小結
6：Gibbs采樣

6.1 重新尋找合適的細致平穩條件
6.2 二維Gibbs采樣
6.3 多維Gibbs采樣
6.4 二維Gibbs采樣python實現
6.5 Gibbs采樣小結
7：參考文獻

一、基本知識

1. 條件概率

條件概率是指在某件事情已經發生的前提下，另一件事情在此基礎上發生的概率，舉例來說P(A丨B)表示B發生的基礎上，A也發生的概率，基本公式為：

2. 條件期望

在上述概率下的期望我們稱之為條件期望。
計算：

3. 隨機過程

在概率論概念中，隨機過程是隨機變量的集合。若一隨機系統的樣本點是隨機函數，則稱此函數為樣本函數，這一隨機系統全部樣本函數的集合是一個隨機過程。實際應用中，樣本函數的一般定義在時間域或者空間域。

顧名思義，它其實就是個過程，比如今天下雨，那么明天下不下雨呢？后天下不下雨呢？從今天下雨到明天不下雨再到后天下雨，這就是個過程。那么怎么預測N天后到底下不下雨呢？這其實是可以利用公式進行計算的，隨機過程就是這樣一個工具，把整個過程進行量化處理，用公式就可以推導出來N天后的天氣狀況，下雨的概率是多少，不下雨的概率是多少。說白了，隨機過程就是一些統計模型，利用這些統計模型可以對自然界的一些事物進行預測和處理，比如天氣預報，比如股票，比如市場分析，比如人工智能。

隨機過程是一族時間函數的集合，隨機過程的每個樣本函數（一個隨機信號）是一個確定的時間函數x(t)（t=0~T），隨機過程在一個確定的時刻t1是一個隨機變量X(t1)。
從上述對隨機過程描述中可以看出，從信號分析的角度來看，一個隨機過程就是針對一個實驗測量無數個的樣本，每個樣本就是一個隨機信號（要分清什么是隨機過程和隨機信號），而所有樣本測量的時間都是0~T（這看起來很麻煩和不可思議），然后，我們把這些樣本（每個隨機信號）放在一起擺放，見圖

其實，這就是一個隨機過程，然后，我們怎樣描述呢？我們從兩個角度來描述它（和上述定義對應）。
（1）從圖中縱向的角度（0~T），一個隨機過程是你這次針對于某個任務進行測量的所有樣本信號的集合。
（2）從橫向的角度（在0-T區間的任意時刻），你可以拿把“菜刀”在任意時刻橫著切一下，這時，你會看到所有你測量的樣本信號在這個時刻的幅度值都不一樣，它們是隨機的，為了表征這些幅度值的統計特性，我們用此時刻的隨機變量X(t)來表征，那么在0~T時間區間內有多少隨機變量呢？當然是無數個（相當于你切了無數刀），這些無數個隨機變量的集合就是隨機過程。
如果所有工作都像這樣去做的話，那么，所有隨機信號分析就沒有辦法完成了，所以，人們將問題進行簡化，盡量將隨機過程說成平穩的并且具有各態歷經的。

二、鞅（martingale）

在概率論中，鞅（martingale）是滿足下述條件的隨機過程：已知過去某一時刻 s 以及之前所有時刻的觀測值，若某一時刻 t 的觀測值的條件期望等于過去某一時刻 s 的觀測值，則稱這一隨機過程是鞅。而于博弈論中，鞅經常用來作為公平博弈的數學模型。

鞅的原名 martingale 原指一類于18世紀流行于法國的投注策略，稱為加倍賭注法[1]。這類策略中最簡單的一種策略是為博弈設計的。在博弈中，賭徒會擲硬幣，若硬幣正面向上，賭徒會贏得賭本，若硬幣反面向上，賭徒會輸掉賭本。這一策略使賭徒在輸錢后加倍賭金投注，為的是在初次贏錢時贏回之前輸掉的所有錢，同時又能另外贏得與最初賭本等值的收益。當賭徒的財產和可用時間同時接近無窮時，他擲硬幣后硬幣正面向上的概率會接近1，由此看來，加倍賭注法似乎是一種必然能贏錢的策略。然而，賭金的指數增長最終會導致使用這一策略的賭徒破產。

鞅的概念首先是由保羅·皮埃爾·萊維于 1934 年提出的，但他只提出了離散時間的版本，而且沒有給予命名。直到 1939 年，約翰?維爾將此概念推廣到連續時間的情況，并且首次提出 martingale 這個名稱。約瑟夫·利奧·杜布（Joseph Leo Doob）等人在鞅的相關理論的初期發展做出重大貢獻，而完成這些工作的部分動機是為了表明成功的投注策略不可能存在。此外，伊藤清在分析應用方面作出了重要的貢獻。從 1970 年代開始，鞅論就在純粹數學和應用數學的很多領域中有廣泛的應用，特別是在數學物理和金融數學中。

以下我們討論T個離散階段下，鞅的相關定義和性質。

我們將一個隨機試驗的所有基本結果稱為樣本點，記為：
稱其構成的集合文件樣本空間：
其對應的概率空間記為：
記每個階段的狀態為S0，S1, …, ST，我們定義以下集合：
Ft可以理解為t時刻前所有狀態的信息，即整個隨機過程到達當前狀態的所有可能路徑。我們再定義:

這樣一來F就包含了所有時刻，所有狀態，所有路徑的信息。以下的性質顯而易見：

因為現t2時刻的信息肯定是要多于t1時刻的。
我們繼續有以下定義：

接下來我們給出條件期望的基本性質：

2.1 Martingale的定義

我們開始討論鞅的定義，鞅有以下三種等價的定義：

我們以下證明三個定義是等價的：
以上便是鞅的基本定義，為了更好的理解鞅，我們可以舉個例子。
假設賭徒在t時刻有 $M_t$ 的財富，如果賭局是公平的，那么他未來財富的期望應當是與當前財富相同，同時不受之前賭博策略的影響（因為我們求的是條件期望）。事實上martingale理論最早便是用于賭博策略的研究，我們現在討論的鞅是公平賭局，也存在像上鞅和下鞅這樣不公平賭局。至于為什么中文會翻譯成鞅呢，因為在法文里，martingale有兩層意思，一是“倍賭策略”（每輸一局就賭資加倍，只要贏一局就扳回所有損失），二是“馬韁繩”，將其翻譯成鞅（馬的韁繩）是用了它的第二層含義。

2.2 Martingale, supermartingale, submartingale

鞅：E[Xt | Fs] = Xs （t>s）；
上鞅：E[Xt | Fs] < Xs （t>s）；
下鞅：E[Xt | Fs] > Xs （t>s）。

通俗理解：Fs表示你在s時刻所掌握的所有信息，那么E[Xt | Fs]就表示你在s時刻所掌握的所有信息對X在未來t時刻預期的判斷，如何和X在s時刻是一致的說明游戲是公平的；上鞅說明對過程X是看衰的，下鞅則是看好。
鞅是公平賭博，下鞅是有利賭博，上鞅是不利賭博。

鞅直觀意義是一條（有固定趨勢的）線。
考慮一個函數f(t)，橫軸是時間軸，縱軸是函數值。這個函數可以是平的（就是 $f (t) = c$ ），也可以單調增，也可以單調減。因此從圖像上看是一條有固定趨勢的線。
如果僅僅是線，那么未免太乏味無趣。概率就是把固定取值的變量 $X = c$ 取代為可以同時取多個值、每個值都有一定權重的隨機變量 $X$ ，相當于引入了多個可能狀態。隨機過程進一步引入一個時間軸，考慮隨機變量的序列 $X_n$ （或者連續版本的 $X_t$ ），因此如果畫在剛剛那個圖里， $X_n(w)$ 可以看做某個大體上是線、局部有些抖動的點列，而 $X_n$ 就是這些點列的集合。
那么離散鞅是一種點列集（隨機過程），這個點列集大體趨勢是 $f (t) = c$ ，只是會有一些抖動，而且是考慮一個抖動圍繞這條線的點列集整體。類似的，上鞅、下鞅就是單調減、單調增函數的點列集。

2.3 Stopping time

關于隨機過程 X1，X2，X3，… 的停時是隨機變量 τ，這一隨機變量具有如下性質：對于每一個時間，事件 τ = t 的發生與否僅取決于 X1，X2，X3，…，Xt 的取值。從定義中可以感受到的直覺是在任一特定時刻 t，我們都可以知道在這一時刻隨機過程是否到了停時。現實生活中停時的例子如賭徒離開賭桌的時刻，這一時刻可能是賭徒以前贏得錢財的函數（例如，僅當他沒有錢時，他才可能離開賭桌），但是他不可能根據還未完成的博弈的結果來選擇離開還是留下。

上述停時定義滿足強條件，下面給出一個弱條件的停時定義：若事件 τ = t 的發生與否統計獨立于 Xt+1，Xt+2，… 但并不是完全決定于時刻 t 以及之前的過程歷史，則隨機變量 τ 是停時。雖然這是一個弱條件，但在需要用到停時的證明中的一些情況也算是足夠強的條件。

鞅的一個基本性質是若 $X_{t})_{{t>0}}$ 是下\上鞅且 $τ\tau$ 是停時，由 $Xtτ:=Xmin?{τ,t}X_{t}^{\tau }:=X_{{\min\{\tau ,t\}}}$ 定義的對應停止過程 $(Xtτ)t>0(X_{t}^{\tau })_{{t>0}}$ 也是下\上鞅。

停時鞅的概念引出了一系列定理，例如可選停止定理（又稱可選抽樣定理）：在特定條件下，停時的鞅的期望值等于其初始值。利用這一定理，我們可以證明對于一個壽命有限且房產有限的賭徒，成功的投注策略不可能存在。

2.3.1 Optional stopping theorem

In probability theory, the optional stopping theorem (or Doob’s optional sampling theorem) says that, under certain conditions, the expected value of a martingale at a stopping time is equal to its initial expected value. Since martingales can be used to model the wealth of a gambler participating in a fair game, the optional stopping theorem says that, on average, nothing can be gained by stopping play based on the information obtainable so far (i.e., without looking into the future). Certain conditions are necessary for this result to hold true. In particular, the theorem applies to doubling strategies.

The optional stopping theorem is an important tool of mathematical finance in the context of the fundamental theorem of asset pricing.

2.3.1.1 discrete-time version of the theorem

2.3.1.2 Applications

2.3.1.3 Proof

三、markov chain

3.1 Markov Property

在學習馬爾可夫鏈之前，我們首先了解一下什么是馬爾可夫性質（或者叫馬爾可夫特性，Markov Property）。
假設你有一個系統，這個系統具有 $M$ 種可能的狀態，并且在狀態之間正在移動。真實世界中，我們常常見到這種例子，例如天氣從炎熱到溫和，或者是股票市場從熊市到牛市再到停滯的狀態（stagnant states）。

3.1.1 馬爾可夫性的定義

關于馬爾可夫性，我們給出了如下的Definition：

從上述的式子可以看出，t+1時刻的狀態包含了1,…,t時刻狀態的全部歷史信息，并且當我們知道t時刻的狀態后，我們只關注于環境的信息，而不用管之前所有狀態的信息，這就是馬爾可夫性，當論文中說某一狀態或其他信息符合馬爾可夫性時，我們也應當聯想到這個性質。

3.1.2 State Transition Matrix狀態傳輸矩陣

如果一個過程具有馬爾可夫性質，它也被成為馬爾可夫過程（Markov Process）。

3.2 Markov Processes馬爾可夫過程

在概率論及統計學中，馬爾可夫過程（英語：Markov process）是一個具備了馬爾可夫性質的隨機過程，因為俄國數學家安德雷·馬爾可夫得名。馬爾可夫過程是不具備記憶特質的（memorylessness）。換言之，馬爾可夫過程的條件概率僅僅與系統的當前狀態相關，而與它的過去歷史或未來狀態，都是獨立、不相關的。

3.3 markov chain馬爾可夫鏈

具備離散狀態的馬爾可夫過程，通常被稱為馬爾可夫鏈。

3.3.1 定義與解釋

3.3.1 定義舉例1

3.3.2 定義舉例2

Markov Chain 體現的是狀態空間的轉換關系，下一個狀態只決定于當前的狀態.通俗點說就是：你下一時刻在哪只跟你現在在哪有關，于你以前去過哪沒有關系. 如下圖：

這個狀態圖的轉換關系可以用一個轉換矩陣 T 來表示：

舉一個例子，如果當前狀態為 u(x) = (0.5, 0.2, 0.3), 那么下一個矩陣的狀態就是 u(x)T = (0.18, 0.64, 0.18), 依照這個轉換矩陣一直轉換下去，最后的系統就趨近于一個穩定狀態 (0.22, 0.41, 0.37) 。而事實證明無論你從那個點出發，經過很長的 Markov Chain 之后都會匯集到這一點。

滿足什么條件下經過很長的 Markov Chain 后系統會趨近一個穩定狀態呢，大概的條件如下：

Irreducibility. 即圖是聯通的，各個狀態之間都有過去的辦法，舉個不聯通的例子，比如爬蟲爬不到內部局域網的網頁
Aperiodicity. 非周期性, 即圖中遍歷不會陷入到一個死圈里，進去了再也出不來，有些網站為了防機器人，會專門設置這種陷阱
Detailed Balance，這是保證系統有穩態的一個重要條件，詳細說明見下面。
假設 p(x) 是最后的穩態，那么 detailed balance 可以用公式表示為：

什么意思呢？假設上面狀態圖 x1 有 0.22 元， x2 有 0.41 元，x3 有 0.37 元，那么 0.22×1 表示 x1 需要給 x2 錢，以此類推，手動計算，可以發現下一個狀態每個人手中的錢都沒有變。值得說明的是，這里體現了一個很重要的特性，那就是從一個高概率狀態 xi 向一個低概率狀態 x(i-1) 轉移的概率等于從這個低概率狀態向高概率狀態轉移的概率（reversible，至于要不要轉移又是另外一回事)
ergodicity, 遍歷性，是指不管事物現在出于什么狀態，在較長時間內，馬爾科夫過程逐漸趨于穩定狀況，而且穩定狀態與初始狀況無關。

3.3.3 通俗解釋

先說說我們村智商為0的王二狗，人傻不拉幾的，見人就傻笑，每天中午12點的標配，仨狀態：吃，玩，睡。這就是傳說中的狀態分布。

你想知道他n天后中午12點的狀態么？是在吃，還是在玩，還是在睡？這些狀態發生的概率分別都是多少？
先看個假設，他每個狀態的轉移都是有概率的，比如今天玩，明天睡的概率是幾，今天玩，明天也玩的概率是幾幾，看圖更清楚一點。

這個矩陣就是轉移概率矩陣P，并且它是保持不變的，就是說第一天到第二天的轉移概率矩陣跟第二天到第三天的轉移概率矩陣是一樣的。（這個叫時齊，有興趣的同學自行百度）。
有了這個矩陣，再加上已知的第一天的狀態分布，就可以計算出第N天的狀態分布了。

S1 是4月1號中午12點的的狀態分布矩陣 [0.6, 0.2, 0.2]，里面的數字分別代表吃的概率，玩的概率，睡的概率。
那么
4月2號的狀態分布矩陣 S2 = S1 * P (倆矩陣相乘)。
4月3號的狀態分布矩陣 S3 = S2 * P (看見沒，跟S1無關，只跟S2有關)。
4月4號的狀態分布矩陣 S4 = S3 * P (看見沒，跟S1，S2無關，只跟S3有關)。
…
4月n號的狀態分布矩陣 Sn = Sn-1 * P (看見沒，只跟它前面一個狀態Sn-1有關)。

總結：馬爾可夫鏈就是這樣一個任性的過程，它將來的狀態分布只取決于現在，跟過去無關！

就把下面這幅圖想象成是一個馬爾可夫鏈吧。實際上就是一個隨機變量隨時間按照Markov性進行變化的過程。

附：S2 的計算過程 (沒興趣的同學自行略過)

3.3.2 n步轉移概率

除了一步轉移以外，馬爾科夫鏈還有n步轉移，即通過n次達到目標狀態

其實是Markov Chain 的基本性質“無后效性”，即事物將來的狀態及其出現的概率的大小，只取決于該事物現在所處的狀態，而與以前時間的狀態無關

3.3.2.1 C-K方程

查普曼-柯爾莫格洛夫方程(Chapman-Kolmogorov equation,C-K equation)給出了計算 [公式] 步轉移概率的一個方法：

3.3.3 狀態分類

①可達的：i到j的概率為正

②常返的：從i出發后還能再回到i

③非常返的：從i出發后不一定能回到j

常返和非常返的區別！：常返指出發后可以無限次回到自身；非常返并不是無法回到自身，而是只能有限次回到自身。

常返類：相互可達的狀態組。

馬爾科夫鏈的分解：

一個馬爾科夫鏈可以分解成一個或者多個常返類，加上可能的非常返類。
一個常返態從它所屬的類里任何一個狀態出發都是可達的，但從其他類里的常返狀態出法是不可達的。
從任何一個常返狀態出發都不可到達非常返狀態
從一個非常返狀態出發，至少有一個常返態是可達的

常返類的周期：

若一個常返類能被分成幾個組，這幾個組之間的下一個狀態轉移必定不在自己組。則這個常返類被稱為有周期的。

命題得證。

很顯然，暫態也是一個類性質。而利用上述性質可以得到：有限馬爾可夫鏈的所有狀態不可能都是暫態，有限不可約馬爾可夫鏈的所有狀態都是常返態。

3.3.4 穩態

3.3.4.1 穩態概率舉例1

3.3.4.2 穩態概率舉例2

舉個具體的例子。社會學家把人按其經濟狀況分為3類：下層，中層，上層，我們用1,2,3表示這三個階層。社會學家發現決定一個人的收入階層最重要的因素就是其父母的收入階層。如果一個人的收入屬于下層類別，則它的孩子屬于下層收入的概率為0.65，屬于中層收入的概率為0.28，屬于上層收入的概率為0.07。從父代到子代，收入階層轉移概率如下

我們用P表示這個轉移矩陣，則

假設第1代人的階層比例為
則前10代人的階層分布如下

我們可以看到，在相同的轉移矩陣作用下，狀態變化最終會趨于平穩。對于第n代人的階層分布，我們有
從表達式上我們可以看到，π是一維向量，P是兩維矩陣，P進行足夠多次自乘后，值趨于穩定。

3.3.5 馬氏鏈平穩分布

在轉移矩陣P作用下達到的平穩狀態，我們稱之為馬氏鏈平穩分布。對于這個特性，有如下精彩定理

我在這里直觀的解釋一下上面定理

條件

（1）非周期馬氏鏈：馬氏鏈轉移要收斂，就一定不能是周期性的。不做特別處理，我們處理的問題基本上都是非周期性的，在此不做多余解釋。
（2）存在概率轉移矩陣P，任意兩個狀態是連通的：這里的連通可以不是直接相連，只要能夠通過有限次轉移到達即可。比如對于a, b, c狀態，存在a->b, b->c，則我們認為a到c是可達的。

結論

（1）不論初始狀態是什么，經過足夠多次概率轉移后，會存在一個穩定的狀態π。
（2）概率轉移矩陣自乘足夠多次后，每行值相等。即

3.3.5.1 馬爾科夫鏈平穩狀態定理的物理解釋

我們再用一個更加簡單的例子來闡明這個定理的物理含義。假設城市化進程中，農村人轉移為城市人的概率為0.5，城市人轉移為農村人的概率為0.1。

假設一開始有100個農村人，0個城市人，每代轉移人數如下

可以看到，城市化進程中馬爾科夫平穩狀態就是農村人轉移為城市人的速度等于城市人轉移為農村人的速度。對于上述轉移矩陣P，平穩分布為農村人17%，城市人83%。如果我們可以得到當前中國城市化轉移矩陣P，我們就可以算出中國最終城市化率大概為多少（這里不考慮P的變化）。同時如果我們知道了中國城市化人口比例，我們就能知道城市化進程還可以持續多少代人。

3.3.5.2 長期頻率解釋

3.3.6 Birth-Death process 生滅過程

3.3.6.1 平穩分布

3.3.6.2 性質

3.3.6.3 例子

3.3.6.4 Probability balance equations 概率平衡方程

3.3.7 吸收概率和吸收期望時間

3.3.7.1 平均時間

平均吸收時間：i的平均吸收時間為從狀態i開始直到到達吸收態所需要的期望步數

平均吸收時間計算：解平均吸收時間方程組

3.3.8 為什么馬爾可夫鏈這么重要

因為它的靜態分布（或叫平穩性分布，stationary distribution）性質。
我們用一下Wiki上給出的股票市場的例子（見下圖），這是一個具有靜態空間的馬爾可夫鏈（a continuous-time Markov chain with state-space），包含了牛市、熊市、停滯市場三種狀態（Bull market, Bear market, Stagnant market）。

它各狀態之間由轉移概率組成的轉移矩陣（ transition rate matrix）如下，記為 Q：

你可以嘗試從其他狀態出發，但最終的靜態分布一定是相同的。那么得出這個靜態分布有什么用處呢？

靜止狀態分布的重要性在于它可以讓你在隨機的時間里，對于一個系統的任意一個狀態確定其概率。（The stationary state distribution is important because it lets you define the probability for every state of a system at a random time.）

對于上述例子，你可以說整個系統，62.5%的概率系統將處于牛市，31.25%的概率將處于熊市，而6.25%的概率將處于停滯態。

直觀上說，你可以想想成在一個狀態鏈上進行隨機游走。你在一個狀態出，然后你決定你前往下一個狀態的方法是在給定當前狀態下觀察下一個狀態的概率分布。

3.3.9 小結

馬爾科夫鏈是一個數學對象，包含一系列狀態以及狀態之間的轉移概率，如果每個狀態轉移到其他狀態的概率只與當前狀態相關，那么這個狀態鏈就稱為馬爾科夫鏈。有這樣一個馬爾科夫鏈之后，我們可以任取一個初始點，然后根據狀態轉移概率進行隨機游走。假設能夠找到這樣一個馬爾科夫鏈，其狀態轉移概率正比于我們想要采樣的分布（如貝葉斯分析中的后驗分布），采樣過程就變成了簡單地在該狀態鏈上移動的過程。

3.4 連續時間馬爾科夫鏈CTMC

3.4.1 定義

3.4.1.1 定義1:一般定義

和起始時間t無關的話，我們稱這是時間齊次的馬爾科夫鏈。這個轉移矩陣和離散時間不同的是，離散時間給出的是一步轉移概率，但是連續馬爾科夫鏈的轉移概率給出的是和時間相關的。

3.4.1.2 定義2

我們考慮連續時間馬爾科夫鏈從一個狀態 i開始，到狀態發生變化，比如變成j所經過的時間，由于馬爾科夫鏈的馬爾科夫性，這個時間是具有無記憶性的，所以這個時間是服從指數分布的。這和離散時間馬爾科夫鏈是密切相關的，離散時間馬爾科夫鏈中的時間是離散時間，因為由無記憶性，所以是服從幾何分布的。

這樣我們就可以這樣定義連續時間馬爾科夫鏈。馬爾科夫鏈是這樣的一個過程。

(i)在轉移到不同的狀態 $i$ 前，它處于這個狀態的時間是速率為 $v_i$ 的指數分布。
(ii)當離開狀態 $i$ 時，以某種概率 $P_{ij}$ 進入下一個狀態 $j$ ,當然 $P_{ij}$ 滿足

對比于半馬爾科夫鏈，我們可以發現，連續時間馬爾科夫鏈是一種特殊的半馬爾科夫鏈，在一個狀態所待的時間是只不過是一個具體的分布–指數分布，而半馬爾科夫鏈只是說所待的時間是任意的一個隨機時間。

3.4.1.3 定義3

3.4.1.4 舉例

3.4.2 兩個定義(1和2)之間的關系

3.4.3 CTMC的極限概率

3.4.4 CTMC極限分布的應用

3.4.5 時間可逆性

3.4.6 截止的馬爾科夫鏈

四、martingale與markov chain的關系

二者都是隨機過程中的一種過程。

Martingale的詞本意是指馬車上馬夫控制馬前進的韁繩（如果我記得沒錯的話），所以從詞源來看刻畫了一種過程前進（未來）與現在出發點關系的含義。具體來說韁繩的作用是使得馬車的前進方向與現在所沖的方向一致，所以在概率上來解釋就是未來發生的所有路徑可能性的期望值與現在的出發點一致。

從這個意義上來說Matingale的核心是說明了一個過程在向前演化的過程中的穩定性性質。但它并沒有說明這個過程是如何到達這一時間點的（是否由上一個時間點所在的位置決定，matingale并沒有說明）。再用馬車的例子來說，Martingale告訴了你馬車在未來是怎么向前走的，中間會有左右的波動（比如馬、車夫走神了，路上有坑馬要繞開，etc.），但整體來說馬是沿著一條直線向前走的。

而馬爾科夫過程的核心在于點明了過程的演化是無記憶性的。還拿馬車來舉栗子，假設車夫喝醉了，他沒有意識并在一個很大的廣場，馬車下一刻前進的方向并不需要是一條直線（經過車夫與馬的直線，這種情況下韁繩是繃直的，是martingale），或者說韁繩由于車夫沒繃緊是松垮的，這種情況下馬車在下一刻可以去任何一個方向，整體上來說前進方向也不必須有什么穩定性規律可循，但整個過程唯一的共性是馬邁出前腿的時候，能夠到達的所有可能范圍，是由它的后腿（你現在所在的位置）決定的（但馬可能扭一扭屁股，身子彎曲一下，所以不必須走直線，不必需走直線，不必需走直線！），而并沒有由上一步馬所在的位置決定，這也就是所謂的無記憶性。

所以從這兩個角度來理解，兩個名詞

有共性：

都從一個過程的全生命角度描畫了一個過程的演進性質，

有重疊：
當還是馬車例子的時候,martingale也是一個markov（因為雖然走直線，但下一刻的位置還是只由現在決定，只是馬身子不能扭曲，不能改變方向），這個例子在概率上最熟悉的模型就是brownian motion了；而反過來，馬車未來位置由現在決定，又走直線，所以此時markov process 也是一個martingale （例子還是brownian motion）；

但更重要的是兩個過程本質上不是在講一回事：比如還是馬車車夫，喝醉了但走在一個三維空間，這時候它是一個markov process,但是由于方向不確定，此時已經不是martingale而變成了一個local martingale；而反過來，假設有一個錯幀宇宙，空間共享但時間差一天，這時候同一個馬車走在不同的宇宙里（但行走軌跡獨立），韁繩拉直，此時兩架馬車都走之前，兩架馬車組成的系統是一個martingale，但是由于下一時刻前進的方向與宇宙1中的此時有關，也與宇宙2中的昨天有關，所以兩架馬車組成的系統就不再是一個markov了。

總結一下，brownian motion (wiener process)既是markov process 又是 martingale; 而markov process 與martingale是相交而非包含與反包含關系。只能說你中有我我中有你，但你不屬于我我也不屬于你

五、Monte Carlo

5.1 蒙特卡羅方法引入

蒙特卡洛（Monte Carlo）方法來自于摩納哥的蒙特卡洛賭場，許多紙牌類游戲需要計算其勝利的概率。用它作為名字大概是因為蒙特卡羅方法是一種隨機模擬的方法，這很像賭博場里面的扔骰子的過程。最早的蒙特卡羅方法都是為了求解一些不太好求解的求和或者積分問題。比如積分：
如果我們很難求解出f(x)的原函數，那么這個積分比較難求解。當然我們可以通過蒙特卡羅方法來模擬求解近似值。(我們可以將蒙特卡洛理解為簡單的模擬，通過模擬的情景來計算其發生的概率。)
如何模擬呢？假設我們函數圖像如下圖:

5.2 概率分布采樣

上一節我們講到蒙特卡羅方法的關鍵是得到x的概率分布。如果求出了x的概率分布，我們可以基于概率分布去采樣基于這個概率分布的n個x的樣本集，帶入蒙特卡羅求和的式子即可求解。但是還有一個關鍵的問題需要解決，即如何基于概率分布去采樣基于這個概率分布的n個x的樣本集。

5.3 接受-拒絕采樣

5.4 蒙特卡羅方法小結

使用接受-拒絕采樣，我們可以解決一些概率分布不是常見的分布的時候，得到其采樣集并用蒙特卡羅方法求和的目的。但是接受-拒絕采樣也只能部分滿足我們的需求，在很多時候我們還是很難得到我們的概率分布的樣本集。比如：

從上面可以看出，要想將蒙特卡羅方法作為一個通用的采樣模擬求和的方法，必須解決如何方便得到各種復雜概率分布的對應的采樣樣本集的問題。而馬爾科夫鏈就是幫助找到這些復雜概率分布的對應的采樣樣本集的白衣騎士。

5.5 應用場景

通常情況下，許多概率的計算非常復雜甚至是無法計算的，但是我們可以通過計算機對期待發生的場景進行大量的模擬，從而計算出其發生的概率。簡單的公式表達為：

學習過概率的同學應當知道這就是典型的頻率學派的觀點。

對于蒙特卡洛方法，經典的例子就是計算 $π\pi$ 值。在（-1，1）之間隨機取兩個數，如果在單位圓內，則記一次，在圓外則不計入次數。

import random import numpy as np from math import sqrt, pinum1 = 1000 # simulation times freq1 = 0 x1 = [] y1 = [] for i in range(1,num1+1):x_i, y_i = random.uniform(-1, 1), random.uniform(-1, 1)square_distance = x_i**2+y_i**2x1.append(x_i)y1.append(y_i)if square_distance <= 1: # 圓的方程x^2+y^2=1freq1 += 1simulation_pi1 = 4*freq1/num1;num2 = 10000 # simulation times freq2 = 0 x2 = [] y2 = [] for i in range(num2):x_i, y_i = random.uniform(-1, 1), random.uniform(-1, 1)square_distance = x_i**2+y_i**2x2.append(x_i)y2.append(y_i)if square_distance <= 1: # 圓的方程x^2+y^2=1freq2 += 1simulation_pi2 = 4*freq2/num2;% matplotlib inline import matplotlib.pyplot as pltfig, ax = plt.subplots(figsize=(8,8))theta = np.linspace(0,2*pi,500) x,y = np.cos(theta), np.sin(theta) ax.plot(x, y, color='red', linewidth=1.0)ax.scatter(x1,y1,color='blue',alpha=0.75) ax.set_xlim(-1,1); ax.set_ylim(-1,1);

結果為：

真實的pi值： 3.141592653589793
蒙特卡洛方法下得出的pi值（模擬1000次）： 3.188
相對誤差（模擬1000次）： 1.4771917153924754 %
蒙特卡洛方法下得出的pi值（模擬10000次）： 3.1292
相對誤差（模擬10000次）： 0.39447041536821975 %

2. 小結

蒙特卡洛就是一種模擬事情發生的方法。

六、MCMC

蒙特卡洛馬爾可夫鏈 Markov Chain Monte Carlo簡稱MCMC，是一個抽樣方法，用于解決難以直接抽樣的分布的隨機抽樣模擬問題。用來在概率空間，通過隨機采樣估算興趣參數的后驗分布。

簡而言之，就是用馬爾科夫鏈方法來抽樣。

抽樣算法的主要任務是找到符合給定分布p(x)的一系列樣本。對于簡單的分布，可以通過基本的抽樣算法進行抽樣。大多數分布都是不容易直接抽樣的，馬爾可夫鏈蒙特卡羅算法解決了不能通過簡單抽樣算法進行抽樣的問題，是一種重要的實用性很強的抽樣算法。

馬爾可夫鏈蒙特卡羅算法（簡寫為MCMC）的核心思想是找到某個狀態空間的馬爾可夫鏈，使得該馬爾可夫鏈的穩定分布就是我們的目標分布p(x)。這樣我們在該狀態空間進行隨機游走的時候，每個狀態x的停留時間正比于目標概率p(x)。在用MCMC進行抽樣的時候，我們首先引進一個容易抽樣的參考分布q(x)，在每步抽樣的過程中從q(x)里面得到一個候選樣本y, 然后按照一定的原則決定是否接受該樣本，該原則的確定就是要保證我們得到的原本恰好服從p(x)分布。

在了解什么是MCMC的時候，我們需要考慮一個情況。舉例，我們已經知道一個分布（例如beta分布）的概率密度函數PDF，那么我們怎么樣從這個分布中提取樣本呢？

MCMC給我們提供了一種方式來從概率分布中進行采樣，而這種方法在貝葉斯統計中的后驗概率進行采樣十分方便。

1. 動因

MCMC通常用于解決高維度的積分和最優化問題，這兩種問題也是機器學習、物理、統計、計量等領域的基礎。比如：

2. MCMC的定義

Wikipedia的解釋如下:

Markov Chain Monte Carlo (MCMC) methods are a class of algorithms for sampling from a probability distribution based on constructing a Markov chain that has the desired distribution as its stationary distribution. The state of the chain after a number of steps is then used as a sample of the desired distribution. The quality of the sample improves as a function of the number of steps.
馬爾可夫鏈蒙特卡羅方法是一類以期望分布為平穩分布的馬爾可夫鏈為基礎，對概率分布進行抽樣的算法。經過一系列步驟之后，鏈的狀態將用作期望分布的樣本。樣品的質量隨著步驟數量的增加而提高。

MCMC方法提供了可以創建一個以Beta分布為其平穩分布（stationary distribution）的馬爾科夫鏈的算法，從而使取樣變得簡單，因為我們可以從一個均勻分布中取樣(這相對容易)。

如果我們基于一些算法，重復地從一個隨機的狀態開始，遍歷到下一個狀態，我們將會創建一個馬爾可夫鏈，一個以beta分布作為其平穩分布的馬爾可夫鏈。

這類算法，一個典型的代表就是Metropolis-Hastings Algorithm算法。

MCMC由梅特羅波利斯(Metropolis)于1949年基于馬爾可夫鏈的基本性質提出，下面介紹一下與馬爾可夫鏈相關的性質。

3. 相關性質

一、穩定分布

穩定分布是指當我們給定狀態空間的一個初始分布 $π0\pi_0$ 以后，按照轉移矩陣進行跳轉最終達到的穩定狀態。

即每個狀態的流出概率等于該狀態注入的概率，該方程稱為全局平衡方程。滿足該方程的分布稱為轉移矩陣A的穩定分布。

二、細致平衡

對于一個穩定分布，我們可以給其一個更強的條件限制，使得任意一個狀態 $i$ 滿足如下條件，

下面我們介紹基本的梅特羅波利斯算法(Metropolis algorithm)。假設p(x)是我們的目標分布函數，我們想得到一系列服從該分布的樣本。我們考慮樣本的產生過程構成一個馬爾可夫鏈，并且讓p(x)是該馬爾可夫鏈的穩定分布，那么該樣本序列就服從p(x)。現在p(x)是已知的，問題的關鍵在于如何構造這個產生過程，也即如何構造馬爾可夫鏈。

即目標概率p(x)滿足細致平衡方程，因此p(x)是轉移概率的穩定分布。

最后，回顧一下梅特羅波利斯抽樣的主要思想。我們首先構造了一個馬爾可夫鏈，接著證明了p(x)滿足其細致平衡方程，進而說明p(x)是該鏈的穩定分布。然后將抽樣的過程看成是在馬爾可夫鏈狀態空間進行跳轉的過程，跳轉的候選狀態由參考分布q(x)產生。最后得到一個的跳轉序列，該序列在每個狀態x的停留時間與p(x)成比，即服從p(x)分布。

4. Metropolis-Hastings (MH) Algorithm算法

梅特羅波利斯－黑斯廷斯算法（英語：Metropolis–Hastings algorithm）是統計學與統計物理中的一種馬爾科夫蒙特卡洛（MCMC）方法，用于在難以直接采樣時從某一概率分布中抽取隨機樣本序列。得到的序列可用于估計該概率分布或計算積分（如期望值）等。梅特羅波利斯－黑斯廷斯或其他MCMC算法一般用于從多變量（尤其是高維）分布中采樣。對于單變量分布而言，常會使用自適應判別采樣（adaptive rejection sampling）等其他能抽取獨立樣本的方法，而不會出現MCMC中樣本自相關的問題。

該算法的名稱源于美國物理學家尼古拉斯·梅特羅波利斯與加拿大統計學家W·K·黑斯廷斯

為了更形象地理解這個算法，我們用下面這個例子來類比。假設我們想知道某個湖的水容量以及這個湖中最深的點，湖水很渾濁以至于沒法通過肉眼來估計深度，而且這個湖相當大，網格近似法顯然不是個好辦法。為了找到一個采樣策略，我們請來了兩個好朋友小馬和小萌。經過討論之后想出了如下辦法，我們需要一個船（當然，也可以是竹筏）和一個很長的棍子，這比聲吶可便宜多了，而且我們已經把有限的錢都花在了船上。

（1）隨機選一個點，然后將船開過去。

（2）用棍子測量湖的深度。

（3）將船移到另一個地點并重新測量。

（4）按如下方式比較兩點的測量結果。

如果新的地點比舊的地點水位深，那么在筆記本上記錄下新的測量值并重復過程（2）。
如果新的地點比舊的地點水位淺，那么我們有兩個選擇：接受或者拒絕。接受意味著記錄下新的測量值并重復過程（2）；拒絕意味著重新回到上一個點，再次記錄下上一個點的測量值。

如何決定接受還是拒絕新的測量值呢？這里的一個技巧便是使用Metropolis-Hastings準則，即接受新的測量值的概率正比于新舊兩點的測量值之比。

事實上，理論保證了在這種情形下，如果我們能采樣無數次，最終能得到完整的后驗。幸運地是，實際上對于很多問題而言，我們只需要相對較少地采樣就可以得到一個相當準確的近似。

5. 數學表達

6. 算法介紹

7. 代碼實現MCMC采樣beta分布

import random # Lets define our Beta Function to generate s for any particular state. # We don't care for the normalizing constant here.def beta_s(w,a,b): # beta distribution pdfreturn w**(a-1)*(1-w)**(b-1)# This Function returns True if the coin with probability P of heads comes heads when flipped. def random_coin(p):unif = random.uniform(0,1)if unif>=p:return Falseelse:return True# This Function runs the MCMC chain for Beta Distribution. def beta_mcmc(N_hops,a,b):states = []cur = random.uniform(0,1)for i in range(0,N_hops):states.append(cur)next = random.uniform(0,1)ap = min(beta_s(next,a,b)/beta_s(cur,a,b),1) # Calculate the acceptance probabilityif random_coin(ap):cur = nextreturn states[-1000:] # Returns the last 1000 states of the chain

在定義了相關函數后，我們來檢查一下我們采樣的結果和實際beta概率分布之間的差距

import numpy as np from matplotlib import pylab as pl import scipy.special as ss % matplotlib inline plt.rcParams['figure.figsize'] = (17.0, 4.0)# Actual Beta PDF. def beta(a, b, i):e1 = ss.gamma(a + b)e2 = ss.gamma(a)e3 = ss.gamma(b)e4 = i ** (a - 1)e5 = (1 - i) ** (b - 1)return (e1/(e2*e3)) * e4 * e5# Create a function to plot Actual Beta PDF with the Beta Sampled from MCMC Chain. def plot_beta(a, b):Ly = []Lx = []i_list = np.mgrid[0:1:100j]for i in i_list:Lx.append(i)Ly.append(beta(a, b, i))pl.plot(Lx, Ly, label="Real Distribution: a="+str(a)+", b="+str(b))plt.hist(beta_mcmc(100000,a,b),density=True,bins =25, histtype='step',label="Simulated_MCMC: a="+str(a)+", b="+str(b))pl.legend()pl.show()plot_beta(0.1, 0.1) plot_beta(1, 1) plot_beta(2, 3)

相關圖形如下所示：

從上面的采樣值可以看出，我們的采樣值還是很接近beta分布本身。雖然我們用了beta分布最為例子，但實際上其他的分布也是可以類似來進行采樣。

針對那些我們無法直接利用共軛分布得到的后驗分布，特別是高維的分布，MCMC方法更顯的尤為重要。

8. MH Algorithm的進一步說明

從貝葉斯的角度看，Metropolis-Hastings算法使得我們能夠從任意分布中以概率p(x)得到采樣值，只要我們能算出某個與p(x)成比例的值。這一點很有用，因為在類似貝葉斯統計的許多問題中，最難的部分是計算歸一化因子，也就是貝葉斯定理中的分母。Metropolis-Hastings算法的步驟如下：

最后，我們會得到一連串記錄值，有時候也稱采樣鏈或者跡。如果一切都正常進行，那么這些采樣值就是后驗的近似。在采樣鏈中出現次數最多的值就是對應后驗中最可能的值。該過程的一個優點是：后驗分析很簡單，我們把對后驗求積分的過程轉化成了對采樣鏈所構成的向量求和的過程。

強烈建議閱讀以下博文，作者用一個簡單的例子實現了metropolis方法，并將其用于求解后驗分布，文中用非常好看的圖展示了采樣的過程，同時簡單討論了最初選取的步長是如何影響結果的。
https://github.com/findmyway/Bayesian-Analysis-with-Python/blob/master/MCMC-sampling-for-dummies.ipynb

9. 一句話總結MCMC

MCMC generates samples from the posterior distribution by constructing a reversible Markov-chain that has as its equilibrium distribution the target posterior distribution.
https://twiecki.io/blog/2015/11/10/mcmc-sampling/

10. MCMC如何解決具體問題

尋找p(x)最大值（simulated annealing）
轉換核的混搭和循環
Gibbs 抽樣
Monte Carlo EM（期望-最大化）
輔助變量抽樣（Auxiliary variable samplers)
可逆跳躍MCMC（Reversible jump MCMC）

詳見

-----------------------------------------------
鏈接：
https://zhuanlan.zhihu.com/p/83314877
https://www.zhihu.com/question/26887292/answer/310904395
https://zh.wikipedia.org/wiki/%E9%9E%85_(%E6%A6%82%E7%8E%87%E8%AE%BA)
https://www.zhihu.com/question/52778636/answer/133133860
https://www.zhihu.com/question/31173033/answer/113202955
https://zhuanlan.zhihu.com/p/116725922
https://zhuanlan.zhihu.com/p/86995916
https://en.wikipedia.org/wiki/Metropolis%E2%80%93Hastings_algorithm

總結

以上是生活随笔為你收集整理的martingale、markov chain、Monte Carlo、MCMC的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。