簡介
稀疏貝葉斯學習(Sparse Bayesian Learning,SBL)是稀疏信號重構的方法之一，其性能相當於重加權的$\ell_1$范數恢復方法，並且不需要設置正則化參數，在目標定位，生物醫學信號提取等方面被廣泛應用。但是其涉及復雜的數學知識包括高斯函數、最大似然估計、向量求導、貝葉斯估計、EM算法等讓很多人望而卻步。筆者在學習此部分內容也曾花費大量時間，為解決小伙伴們的煩惱，本系列文章將詳細解讀稀疏貝葉斯學習的基本原理及其對應的數學推導，大致分為幾塊，包括證據和后驗概率的計算、EM算法部分推導等。下面先對證據和后驗概率的計算推導進行敘述。以下需要用到的數學基礎包括高斯函數的基本性質、向量的求導。
模型
先考慮對一個向量的觀測，假設有觀測矩陣 $\bm{\Phi}\in C^{N\times M}$，對未知變量$\bm{\omega}\in C^{M\times1}$進行觀測，記為
\[\bm{t}=\bm{\Phi}\bm{\omega}+\bm{\epsilon}\qquad(1)\]
式中$t\in C^{N\times1}$，觀測矩陣也稱之為過完備基，這里假定$\bm{\omega}$是稀疏變量，即$\bm{\omega}$的大部分元素都為0，$\epsilon$ 為觀測噪聲。SBL要解決的問題是根據已知的$\bm{t}$和${\bm{\Phi}}$估計出$\bm{\omega}$,其實就是稀疏信號的重構。
首先解釋下貝葉斯公式：
\[p(\omega|t)=\frac{p(t|\omega)p(\omega)}{p(t)}\qquad\]
$p(\omega)$稱之為先驗概率，表示在觀測之前的概率，$p(\omega|t)$稱之為后驗概率，是觀測之后的概率，$p(t|\omega)$是似然概率，在求最大似然估計的時候就是使用的該概率形式，$p(t)$表示證據。很多情況下，我們要估計$\bm{\omega}$可由$argmax_\omega p(\omega|x)$求得，但上述后驗概率不易求得。因證據$p(x)$與$\bm{\omega}$無關，上述后驗概率最大化可由貝葉斯公式轉化為似然概率和先驗概率的乘積的最大化求得，即$argmax_\omega p(x|\omega)p(\omega)$。
證據推導
SBL采用了神經網絡里常用的自動相關決策理論(Automatic Relevance Determination)來獲取稀疏解。首先假定$\bm{\epsilon}$符合均值為0，方差為$\sigma^2\bm{I}_N$的高斯分布，則可得出$\bm{t}$符合均值為$\bm{\Phi}\bm{\omega}$，方差為$\sigma^2\bm{I}_N$的高斯分布,即
\[p(\bm{t}|\bm{\omega})=(2\pi\sigma^2)^{-N/2}exp[-\frac{1}{2\sigma^2}(\bm{t}-\bm{\Phi\omega})^H(\bm{t}-\bm{\Phi\omega})]\qquad(2)\]
根據ARD，其假定$\bm{\omega}$由超參數$\bm{\gamma}$產生，假定其$\omega_i$由$\gamma_i$控制，並符合均值為0，方差為$\gamma_i$的高斯分布，即
\[p(\bm{\omega};\bm{\gamma})=(2\pi)^{\frac{-M}{2}}\left|\bm{\Gamma}\right|^{-\frac{1}{2}}e^{-\frac{1}{2}\bm{\omega}^H\bm{\Gamma^{-1}\omega}}\qquad(3)\]
式中$\bm{\Gamma}=diag(\bm{\gamma})$。
利用全概率公式即可得第二類似然函數為
\[p(\bm{t};\bm{\gamma})=\int _{\bm{\omega}}{p(\bm{t}|\bm{\omega})p(\bm{\omega};\bm{\gamma})d\bm{\omega}}\]
將(2)和(3)代入到(4)中，可得
\[p(\bm{t};\bm{\gamma})=\int_{\bm{\omega}}(2\pi\sigma^2)^{-N/2}(2\pi)^{-M/2}\left|\bm{\Gamma}\right|^{-\frac{1}{2}}exp[-\frac{1}{2\sigma^2}(\bm{t}-\bm{\Phi\omega})^H(\bm{t}-\bm{\Phi\omega})-\frac{1}{2}\bm{\omega}^H\bm{\Gamma^{-1}\omega}]d\bm{\omega}\]
其實該式可以看成兩個高斯函數進行卷積，根據高斯函數性質知，兩個高斯函數卷積的結果仍為高斯函數。所以只需要求得卷積后的高斯函數的均值和期望，就相當於求出上式的積分了。
取其指數，令
\[L=-\frac{1}{2\sigma^2}(\bm{t}-\bm{\Phi\omega})^H(\bm{t}-\bm{\Phi\omega})-\frac{1}{2}\bm{\omega}^H\bm{\Gamma^{-1}\omega} \qquad(5)\]
進一步，可以得到
\[L=-\frac{1}{2\sigma^2}[\bm{\omega}^H(\bm{\Phi}^H\bm{\Phi}+\sigma^2\bm{\Gamma}^{-1})\bm{\omega}-\bm{t}^H\bm{\Phi}\bm{\omega}-\bm{\omega}^H\bm{\Phi}^H\bm{t}+\bm{t}^H\bm{t}]\]
$L$是關於$\bm{\omega}$的二次項。這里求解上述積分要用到高斯函數的以下性質：
\[\int_{\bm{\omega}}e^{-(\bm{A\omega}+\bm{b})^2}d\bm{\omega}=C\]
式中$\bm{A}$是矩陣，$\bm{b}$是向量，其維數應滿足上式的乘法規則。$C$是常數，具體是多少，我們可以不關注，感興趣的話可以自己推導或查閱相關文獻。我們需要關注的是似然函數對$\bm{\omega}$積分后$\bm{t}$項和$\sigma$項。現在的問題是我們需要將$L$表達成$-(\bm{A\omega}+\bm{b})^2+f(t,\sigma^2)$的樣式，並求得$f(t,\sigma^2)$。顯然，我們將滿足$\bm{A\omega}+\bm{b}=\bm{0}$的$\bm{\omega}$代入其中，即得到$f(t,\sigma^2)$。先求$\bm{\omega}$，下面通過求導完成。
\[\frac{dL}{d\bm{\omega}}=\frac{1}{\sigma^2}[(\bm{\Phi}^H\bm{\Phi}+\sigma^2\bm{\Gamma}^{-1})\bm{\omega}-\bm{\Phi}^H\bm{t}]\]
令$\frac{dL}{d\bm{\omega}}=0$可得
\[\bm{\omega}=(\bm{\Phi}^H\bm{\Phi}+\sigma^2\bm{\Gamma}^{-1})^{-1}\bm{\Phi}^H\bm{t}\qquad(6)\]
將(6)代入(5)中，得到
\[L=-\frac{1}{2\sigma^2}\bm{t}^H[\bm{I}-\bm{\Phi}(\bm{\Phi}^H\bm{\Phi}+\sigma^2\bm{\Gamma}^{-1})^{-1}\bm{\Phi}^H]\bm{t}\]
因此全概率公式積分后得
\[p(\bm{t};\bm{\gamma})=Cexp\{-\frac{1}{2\sigma^2}\bm{t}^H[\bm{I}-\bm{\Phi}(\bm{\Phi}^H\bm{\Phi}+\sigma^2\bm{\Gamma}^{-1})^{-1}\bm{\Phi}^H]\bm{t}\}\]
現在可以看出$p(\bm{t};\bm{\gamma})$是一個高斯分布，其均值為$\bm{0}$,協方差矩陣$\Sigma_t$滿足$\Sigma_t^{-1}=\frac{1}{\sigma^2}[\bm{I}-\bm{\Phi}(\bm{\Phi}^H\bm{\Phi}+\sigma^2\bm{\Gamma}^{-1})^{-1}\bm{\Phi}^H].$
$\Sigma_t$可由矩陣求逆公式得到，如下：
\[\Sigma_t=\sigma^2 \bm{I}+\bm{\Phi\Gamma\Phi}^H\]
到此，我們完成了證據或者叫第二類似然函數的概率分布的推導。
后驗概率推導
下面我們繼續完成后驗概率的推導，根據貝葉斯公式，有
\[p(\bm{\omega}|\bm{t};\bm{\gamma})=\frac{p(\bm{t}|\bm{\omega})p(\bm{\omega};\bm{\gamma})}{p(\bm{t};\bm{\gamma})} \qquad(7)\]
其實利用前面的結果，該式大部分都求得差不多了。證據（分母部分）已求得。分子部分是兩個高斯概率密度函數的乘積，其結果仍為高斯分布。再與分母部分相除，最終還是為高斯分布。將前面求得的結果分別代入到(7), 忽略常數部分，得
\[p(\bm{\omega}|\bm{t};\bm{\gamma})=exp\{-\frac{1}{2\sigma^2}[\bm{\omega}^H(\bm{\Phi}^H\bm{\Phi}+\sigma^2\bm{\Gamma}^{-1})\bm{\omega}-\bm{t}^H\bm{\Phi}\bm{\omega}-\bm{\omega}^H\bm{\Phi}^H\bm{t}+\bm{t}^H\bm{t}]+\frac{1}{2}\bm{t}^H\bm{\Sigma}^{-1}\bm{t}\}\]
其均值為指數部分對$\bm{\omega}$的一階導數零點，協方差矩陣的逆為指數部分對$\bm{\omega}$的二階導數。
故
\[\Sigma_{\omega}^{-1}=\frac{1}{\sigma^2}\bm{\Phi}^H\bm{\Phi}+\bm{\Gamma}^{-1}\]
\[\bm{\mu_{\omega}}=(\bm{\Phi}^H\bm{\Phi}+\sigma^2\bm{\Gamma}^{-1})^{-1}\bm{\Phi}^H\bm{t}\label{omega}\]
一般情況下，$M$往往遠大於$N$,所以求$\bm{\Sigma_{\omega}}$的逆的復雜度遠遠高於$\bm{\Sigma_t}$的逆的復雜度，所以運用矩陣和求逆公式將$\bm{\Sigma_{\omega}^{-1}}$轉化為求$\bm{\Sigma_t^{-1}}.$結果如下：
\[\Sigma_{\omega}=\bm{\Gamma}-\bm{\Gamma\Phi}^H\bm{\Sigma}_t^{-1}\bm{\Phi\Gamma}\]
\[\bm{\mu_{\omega}}=\bm{\Gamma\Phi}^H\bm{\Sigma}_t^{-1}\bm{t}\]
至此，關於稀疏貝葉斯算法中的證據和后驗概率的推導解釋完畢，對於多測量模式下（Multiple Mearsure Vector）的推導可以直接拓展過來，這里不進行詳述，可以參考相關文獻。
Reference
[1] D. P. Wipf and B. D. Rao, "Sparse Bayesian learning for basis selection," IEEE Transactions on Signal Processing, vol. 52, no. 8, pp. 2153-2164, 2004.
[2] D. P. Wipf and B. D. Rao, "An empirical Bayesian strategy for solving the simultaneous sparse approximation problem," IEEE Transactions on Signal Processing, vol. 55, no. 7, pp. 3704-3716, 2007.
附矩陣求逆公式
\[(A+UBV)^{-1}=A^{-1}-A^{-1}UB(I+VA^{-1}UB)^{-1}VA^{-1} \]