2.3 模型優化概率方法

機器學習領域中有很多模型優化方法通過概率工具實現問題求解。從直觀上看，算法應該追求盡可能地穩定，引入概率方法有時會破壞算法的穩定性。但機器學習領域的很多模型優化問題都需要處理巨大的求解空間或可能存在的不確定性信息，此時難以通過基本優化方法實現求解，概率方法則可以通過適當的概率工具有效地解決這些問題。因此，模型優化的概率方法是機器學習理論研究和應用開發領域非常重要的基礎知識。本節主要介紹三種常用的模型優化概率方法，即隨機梯度法、最大期望法和蒙特卡洛法。

2.3.1 隨機梯度法

隨機梯度法是對梯度下降法的一種改進。梯度下降方法在機器學習模型進行優化時，其搜索方向的每次更新都需要計算所有訓練樣本。例如，對于模型f（X；β），其中β為模型參數向量，即β={β₁，β₂，…，β_t}，若用訓練樣本集合S={（X₁，y₁），（X₂，y₂），…，（X_n，y_n）}對該模型進行訓練并將經驗風險作為優化目標，則優化的目標函數為

使用梯度下降法對上述目標函數進行優化時，第k次迭代的搜索方向計算公式為

其中，β_k為第k次迭代的起始點。

由于訓練樣本集S包含n個樣本，故在計算P_k時需分別計算這n個樣本的損失函數對參數向量β的梯度并求出它們的均值。當n較大時，計算梯度需要耗費大量時間。由于提高模型泛化性能需要盡可能多的訓練樣本，故不能通過減少訓練樣本的方式來簡化計算。

事實上，當訓練樣本集S較大時，可用隨機梯度法實現對目標函數的優化。隨機梯度法的基本思想是隨機選擇少量訓練樣本來計算梯度，并將該梯度作為在全部訓練樣本上梯度的近似代替值用于梯度下降的迭代計算。隨機梯度法有很多具體的實現算法，隨機梯度下降法是其中最基本也是最常用的一種。下面具體介紹隨機梯度下降法。

隨機梯度下降法每次迭代的方向計算只與單個樣本有關。對于機器學習模型f（X；β）及其訓練樣本集S={（X₁，y₁），（X₂，y₂），…，（X_n，y_n）}，該方法每次隨機選取S中單個樣本代替全部樣本對模型參數進行一次更新，然后換另一個未參與訓練的樣本進行下次更新，當S中的全部樣本都參與更新計算之后，隨機調整S中所有樣本排列次序后重復以上過程。

具體地說，假設隨機梯度下降法在進行第k次迭代時隨機選擇的樣本為（X_i，y_i），則此次參數更新的搜索方向可通過該樣本的損失函數計算得到。令β_k為模型的當前參數向量，則用函數F′（β）=L（f_βk（X_i），y_i）實現此次迭代搜索方向的更新。需要注意的是，函數F′（β）為隨機選擇到的樣本（X_i，y_i）的損失函數，僅用于計算參數搜索的更新方向，模型優化的目標函數則保持不變。方向更新的具體計算公式為

顯然，使用上述公式僅需計算一次梯度便可實現方向更新。用代替P_k，并進行第k次參數更新的結果為

由以上分析可知：梯度下降法在n元訓練樣本集S上對模型進行優化時，每次參數更新需要執行n次梯度計算。使用隨機梯度下降法只需執行一次梯度計算，可以有效減小計算量，使得模型優化過程在大樣本場合變得可行。隨機梯度下降法的具體過程如下：

（1）初始化t=0，k=0；

（2）對S中的全部訓練樣本進行隨機排序，得到新的樣本集合S^t；

（3）若S^t中的樣本均參與過訓練過程，則令t=t+1，并重復步驟（2），否則隨機選擇S^t中未參與訓練的單個樣本，并根據該樣本數據計算參數更新方向P_k；

（4）根據式（2-36）更新參數向量得第k次參數更新的結果β_k+1；

（5）判斷迭代是否滿足算法終止條件，若滿足則終止算法；否則令k=k+1，重復步驟（3）、（4）。

事實上，目前的計算機大多采用多核架構，具有一定的并行計算能力，若每次只使用單個樣本確定參數搜索方向的更新，則會造成計算能力的浪費。因此，小批量隨機梯度下降法也是一種常用的隨機梯度法。

小批量隨機梯度下降法的基本步驟與隨機梯度下降法類似，基本思想是首先從S中隨機選擇一小批訓練樣本代替全部樣本以實現對模型參數搜索方向的更新，然后在下次迭代搜索中再換另一小批不同樣本進行搜索方向的更新計算，當訓練樣本集S中的全部樣本都參與過一次參數的更新之后，將S中的全部樣本進行隨機排序后重新劃分成不同小批次再去更新搜索方向，直至目標函數的取值接近最小值。顯然，隨機梯度下降法可以看作是小批量隨機梯度下降法在每個小批量樣本數均為1時的一種特殊情況。

現結合實例說明小批量隨機梯度下降法的具體過程。假設訓練集S共有m=100000個樣本，將其中全部樣本隨機排序后得到集合S′，即有

首先，將S′劃分為若干小批量訓練樣本集，例如，將S′等分為1000個規模為100的小批量訓練樣本集S₁，S₂，…，S₁₀₀₀，即有

第1個小批量訓練集，…，，

第2個小批量訓練集，…，，……

第1000個小批量訓練集，…，。

然后，在每次迭代過程中隨機選擇一個當前未參與訓練的小批次訓練集進行搜索方向的更新計算。假設當前模型參數向量為β_k并選擇第i個批次進行方向更新，則有

據此計算更新方向為

依據上述更新方向對當前參數向量β_k進行更新

當所有劃分的批次均參與了訓練過程而未達到算法終止條件時，對S′進行隨機排序后重新劃分成若干小批次樣本集合，并重復上述搜索方向和參數更新過程，直至達到算法終止條件時結束算法。

隨機梯度法每次僅選擇單個或小批量樣本來確定搜索方向和實現對模型參數的更新，可有效減少花費在梯度計算上的時間。該方法的計算結果存在一定的隨機性，即不能保證每次更新的方向都是目標函數減小的方向，參數更新偶爾會使得目標函數取值增大。然而，取自相同總體的樣本服從相同的概率分布，并且所有訓練樣本均參與迭代計算，故在參數更新的整體過程中，隨機梯度下降法可以保證目標函數的取值不斷下降。隨機梯度法的參數更新過程如圖2-6所示。

圖2-6 隨機梯度法參數更新示意圖

另一方面，隨機梯度方法的參數更新方向存在著一定的隨機性，這也為避免陷入局部最優提供了可能。正常情況下梯度下降算法很難跳出局部最優，但是隨機梯度法打破了參數更新的穩定性，這有可能會使迭代過程中出現跳出當前局部最小值范圍的情況。

目前，很多機器學習任務的訓練樣本集規模都很大，由以上分析可知隨機梯度方法能夠在保證優化效果的前提下有效減少在大規模訓練樣本集上對模型進行優化的計算量。因此，隨機梯度法是目前最常用的模型優化方法之一。

2.3.2 最大期望法

有些機器學習模型含有隱含變量，通常難以對這些具有隱含變量的模型直接進行參數求解或估計，此時可通過迭代求解隱含變量（或其函數）數學期望最大值的方法實現對帶隱含變量模型的參數優化求解。這類優化方法通常稱為最大期望法或EM算法，主要包括兩個計算步驟，即計算數學期望E的步驟和計算函數最大值優化M的步驟。可用EM算法解決帶隱含變量模型參數的最大似然估計或者最大后驗估計問題。這里主要討論如何通過EM算法求解參數的最大似然估計，最大后驗估計的EM求解方法與此類似。

假設某學校男生和女生身高分別服從兩個參數不同的高斯（正態）分布，現從該學校所有學生當中隨機抽取100名學生，分別測量他們的身高數據但不記錄性別信息。根據這些數據對男生身高和女生身高所服從正態分布進行參數估計時，則存在一個隱含的性別數據。由于不知道性別信息，無法直接得知性別數據的分布，故無法直接求得已知樣本數據所滿足的似然函數。然而，在給定參數情況下，可以結合已知樣本計算出某個同學性別的概率分布情況，并且可以知道在給定參數情況下性別與身高的聯合概率分布。

更一般地，在一個包含隱含數據Z的模型參數估計問題中，X為可直接觀測數據，即已知樣本，β為模型參數向量。由于隱含數據Z的存在，通常無法直接得知已知樣本X取值下的似然函數L（β|X）=p（X|β），但可知參數給定情況下X和Z的聯合概率分布p（X，Z|β），以及在參數向量和可觀測數據給定情況下Z取值狀態的條件概率分布p（Z|X，β）。

現根據以上信息對模型參數向量β進行最大似然估計。最大似然估計通過最大化似然函數實現，然而已知樣本X取值狀態的似然函數L（β|X）=p（X|β）卻由于隱含數據Z的存在而難以直接求得，故難以直接使用似然函數對模型參數向量β進行最大似然估計。由于已知參數β給定條件下X和Z的聯合概率分布為p（X，Z|β），故將上述似然函數轉化為

考慮其對數似然，有

根據最大似然估計思想，只需求得對數似然lnL（β|X）的最大值即可得到模型參數β的最大似然估計值。但是由于上式存在隱含數據Z和積分的對數，直接求解其最大值較為困難，故用迭代逼近方法實現對數似然最大值的估計。為此，將對數似然做如下變形

其中，p（Z）為隱含數據Z的某一分布。

由于對數函數為凸函數，故如下不等式成立

由此可得對數似然的下界函數

取p（Z）=p（Z|X，β_t），則可得

略去下界函數B（β，β_t）中與待求參數向量β無關的項，得到如下Q函數

上式表示對隱含變量Z的函數L（β|X，Z）在概率分布p（Z|X，β_t）下的數學期望。由于B（β，β_t）≤lnL（β|X），故可通過迭代選取不同下界函數B（β，β_t）最大值的方式逐步逼近對數似然lnL（β|X）的最大值，迭代逼近的具體過程如圖2-7所示。

圖2-7 EM算法逼近似然函數最大值

由于Q（β，β_t）是通過省略B（β，β_t）中與待求參數向量β無關項得到的，用Q（β，β_t）迭代求解對數似然最大值與用B（β，β_t）求解對數似然最大值具有相同的效果，故EM算法通常直接使用Q函數進行優化計算。由以上分析可得EM算法的基本步驟如下。

（1）設置初始參數β₀和迭代停止條件；

（2）E步（期望步）：根據可直接觀測數據X和當前參數向量取值β_t計算Q（β，β_t）：

Q（β，β_t）=∫lnL（β|Z）p（Z|X，β_t）dZ

（3）M步（最大化步）：最大化Q（β，β_t）并根據Q（β，β_t）的最大值更新參數β_t的取值：

（4）判斷是否滿足迭代停止條件，若滿足則停止迭代；否則令t=t+1并返回步驟（2）。

EM算法可隨機選擇初始參數β₀，但應注意EM算法對初始參數β₀具有一定的敏感性。如果初值β₀選取不當，則可能會使迭代結果陷入局部最優。

通常將迭代停止條件設為相鄰兩次參數值變化不大或前后兩次參數更新使得Q（β，β_t）值變化很小時停止更新，即有

|β_t+1-β_t|<ε或|Q（β_t+1，β_t）-Q（β_t，β_t）|<ε

【例題2.9】假設X₁，X₂取自指數分布Y（θ）=θe^-xθ，且X₁，X₂不相關。若X₁=10而X₂缺失，試用EM算法實現參數θ的最大似然估計。

【解】由于X₁，X₂為取自指數分布Y（θ）=θe^-xθ的離散值，則有

由于X₁，X₂同分布，則有，從而有

為求出使得Q（θ，θ_t）最大化的參數θ，可將使得Q（θ，θ_t）最大化的參數作為第t+1次的估計值θ_t+1，即得迭代公式

當迭代收斂時，有θ_*=2θ-/（10θ_*+1），解得θ_*=0.1，即迭代算法求得參數估計值收斂于0.1，即參數θ的最大似然估計值為θ_*=0.1。□

2.3.3 蒙特卡洛法

蒙特卡洛法是一類以概率統計理論為基礎的隨機型數值算法。一般來說，當一個隨機算法滿足采樣次數越多，其輸出結果越近似最優解這一特性時，便可稱為蒙特卡洛法。該方法通常將待求解問題與某一概率模型聯系起來，利用從大量樣本中獲得的信息來完成對所求參數的概率估計，由此實現對實際問題的求解。

例如，可用蒙特卡洛法近似計算一個不規則湖面的面積。假設圍住湖面的長方形面積為M，湖面面積S為未知數。由于湖面形狀不規則無法直接求出S，為此向包含湖面的長方形區域內隨機撒布n個點，假設其中有k個點落在湖面當中，則可得到湖面大小S的估計值

顯然，撒布點數n越大，估計值就越接近于湖面的真實面積S。這是因為根據大數定律和中心極限定理，重復進行大量試驗時，事件A出現的頻率接近事件A發生的概率。

蒙特卡洛法最早用于近似求解難以精確計算的定積分，對于某函數φ（x）的積分

若能找到定義在（a，b）上的一個函數f（x）和概率密度函數p（x），滿足φ（x）=f（x）p（x），則可將上述積分I轉化為

如此一來，就將原積分I的求解問題轉化為了求解函數f（x）在p（x）分布上的數學期望問題。假設在分布p（x）上做隨機采樣得到的樣本集合為{x₁，x₂，…，x_n}，則可用這些樣本來對I進行估計，即有

由于

故I_p是I的無偏估計，即對于任意分布p（x），積分I的蒙特卡洛法估計值是無偏的。

有時在分布p（x）上進行采樣會出現一些問題，例如采樣偏差較大等，甚至在某些時候根本無法對其進行采樣。此時，可通過變換積分的分解形式找到另外一個易于采樣的分布q（x），將積分I轉化為如下形式

其中，q（x）為某一分布，可將f（x）[p（x）/q（x）]看作是某一函數。

令g（x）=f（x）p（x）/q（x），則有

此時可將積分I看作是函數g（x）在分布q（x）上的期望。假設在分布q（x）上進行隨機采樣獲得的樣本集合為{，，…，}，則可用該樣本集合對I進行估計，即有

在新分布q（x）上進行采樣的過程被稱為重要性采樣，p（x）/q（x）被稱為重要性權重。由于積分I任意分解形式通過蒙特卡洛法得到的估計值都是無偏的，故應適當選擇新分布q（x），使得估計值I_q的方差盡可能地達到最小。對于I_q的方差Var（I_q）

由于其后項與分布q（x）無關，故只需最小化期望E_q[f²（x）p²（x）/q²（x）]。

根據琴生不等式，有

故取q（x）為如下分布時可使得方差Var（I_q）達到最小

重要性采樣可以改進原采樣，如圖2-8a所示，若在分布p（x）上進行采樣，則很可能僅僅采樣到函數f（x）的極端值；如圖2-8b所示，若在分布q（x）上進行重要性采樣，則函數f（x）的取值較為豐富。

在實際應用當中，重要性采樣并非總能奏效。在很多情況下，重要性采樣q（x）難以取得理論上的最佳分布q^*（x），而只能取到一個方差較大的可行分布，難以滿足任務需求。此時用馬爾可夫鏈蒙特卡洛法（Markov Chain Monte Carlo，MCMC）近似計算待求參數概率分布的方式解決。

假設離散型隨機變量X的取值范圍是{X₁，X₂，…，X_n}，則稱該集合為X的狀態空間。馬爾可夫鏈是一個關于離散型隨機變量取值狀態的數列，從X的隨機初始狀態X_i開始，馬爾可夫鏈依據一個只與前一時序狀態相關的狀態轉移分布P（X^t+1|X^t）來確定下一時序的狀態，其中X^t和X^t+1分別表示隨機變量X在第t時序和第t+1時序的狀態。

圖2-8 不同采樣方式對比示意圖

a）在分布p（x）上進行采樣 b）在分布q（x）上進行采樣

MCMC法同時運行多條馬爾可夫鏈，這些馬爾可夫鏈在同一時序內具有相同的狀態空間，即服從同一個狀態概率分布。隨著馬爾可夫鏈的不斷更新，其狀態概率分布最終會收斂于某一分布。具體來說，MCMC法從一個初始概率分布p⁰中隨機選擇多條馬爾可夫鏈的初始狀態，p⁰的具體形式為

其中，表示在時序為0的初始狀態分布下X=X_i的概率。

由于這些馬爾可夫鏈的狀態空間是一個n維空間，故可將概率分布p⁰表示為一個n維向量p⁰=（，，…，）^T。此時，馬爾可夫鏈下一時序的分布可表示為

上式說明馬爾可夫鏈在第t+1時序的狀態分布只與前一時序的狀態分布p^t和狀態轉移分布P（X^t+1|X^t）相關。

接著考慮所有同時運行的馬爾可夫鏈狀態概率分布情況。若用矩陣A表示P（X^t+1|X^t），矩陣A中第i行第j列元素a_ij表示從狀態X_j轉移到狀態X_i的概率P（X_i|X_j），則這些馬爾可夫鏈在同一時序上狀態概率分布的變化可表示為

狀態概率分布會隨著時序t不斷變大而最終收斂，即有p^t+1=p^t。此時可用所求狀態概率分布作為真實狀態概率分布的一種近似分布。

事實上，由于貝葉斯學派認為模型的未知參數服從于某一概率分布，并將參數的概率分布看作一種狀態概率分布，故可通過MCMC法求解其近似分布的方式實現模型優化。

MCMC法的基本思想如圖2-8所示。對于參數的概率分布，MCMC法從隨機初始點開始進行采樣，隨著采樣的不斷進行，其樣本點的概率分布會逐步收斂。圖中圓點表示概率分布還未收斂時采樣到的點，叉點表示概率分布收斂后采樣到的點。由于最終收斂的概率分布是對參數概率分布的一種近似，故在其中進行大量采樣便可求出該近似概率分布，即求得參數概率分布的一種近似分布。

圖2-9 MCMC法的基本思想

顯然，MCMC法的關鍵在于如何確定狀態轉移概率分布P（X^t+1|X^t），如果P（X^t+1|X^t）的選擇合適，則狀態概率分布會收斂到待求參數的概率分布。從收斂性角度考慮，MCMC法的狀態轉移概率分布P（X^t+1|X^t）應保證多條馬爾可夫鏈同時滿足狀態分布收斂的條件。

事實上，當狀態概率分布p（X）與狀態轉移概率分布P（X^t+1|X^t）滿足細致平穩條件時，狀態轉移分布P（X^t+1|X^t）就可以使狀態概率分布收斂，此時的狀態分布p（X）稱為狀態轉移分布P（X^t+1|X^t）的平穩分布。

所謂細致平穩條件，是指在狀態分布p（X）中，利用狀態轉移分布對其狀態進行更新時，由X_i轉移到X_j的概率應與由X_j轉移到X_i的概率相同，即有p（X_i）P（X_j|X_i）=p（X_j）P（X_i|X_j）。

在實際應用中找到一個狀態轉移分布并使其滿足細致平穩條件有時是一件比較困難的事情，通常需要使用某些采樣方法構造一個滿足條件的狀態轉移分布。其中最簡單的采樣方法為MCMC采樣。具體地說，對于狀態分布p（X），可隨機選擇一個狀態轉移分布Q（X^t+1|X^t），它們通常不滿足細致平穩條件，即

但可選擇一對轉移算子α（X_i，X_j）和α（X_j，X_i）滿足下列等式

可用轉移算子α（X_i，X_j）和α（X_j，X_i）幫助確定是否接受采樣結果。要想使得上述等式成立，最簡單的方法是令

并取狀態轉移分布為P（X_i|X_j）=Q（X_j|X_i）α（X_i，X_j），則有狀態分布p（X）即為P（X_i|X_j）的平穩分布，此時就可通過MCMC方法獲得一個收斂的狀態概率分布。

MCMC采樣的基本步驟如下：

（1）隨機選定一個狀態轉移分布Q（X_j|X_i）并記狀態分布為p（X），設定收斂前狀態轉移次數上限τ和采樣個數n；

（2）從任意分布中采樣獲得初始狀態X₀；

（3）依據條件狀態分布Q（X_j|X_i）進行采樣，即采樣值為X^t+1=Q（X|X^t）；

（4）從均勻分布U[0，1]中采樣并記為u_t，若u_t<α（X^t，X^t+1）=p（X^t）Q（X^t+1|X^t），則接受X^t+1，并令t=t+1，否則，拒絕此次轉移并保持t值不變；

（5）當t>τ+n時終止算法。此時共取到τ+n個樣本，但前τ個樣本在分布未收斂時獲得，不能作為MCMC采樣結果，故只將第τ+1到τ+n個樣本作為采樣結果。

除了MCMC采樣之外，還有MH采樣、Gibbs采樣等獲取樣本方法。其中MH采樣是對MCMC采樣的改進，Gibbs采樣是對MH采樣的改進，這里不再贅述。