ICA快速算法原理和程序

本文格式為Word版，下載可任意編輯——ICA快速算法原理和程序

試驗2：FastICA算法

一．算法原理：

獨立分量分析（ICA）的過程如下圖：在信源s(t)中各分量相互獨立的假設下，由觀測x(t)通過解混系統(tǒng)B把他們分開開來，使輸出y(t)迫近s(t)！

圖1-ICA的一般過程

ICA算法的研究可分為基于信息論準則的迭代估計方法和基于統(tǒng)計學的代數(shù)方法兩大類，從原理上來說，它們都是利用了源信號的獨立性和非高斯性?；谛畔⒄摰姆椒ㄑ芯恐?，各國學者從最大熵、最小互信息、最大似然和負熵最大化等角度提出了一系列估計算法。如FastICA算法,Infomax算法，最大似然估計算法等?；诮y(tǒng)計學的方法主要有二階累積量、四階累積量等高階累積量方法。本試驗探討FastICA算法。

1.數(shù)據(jù)的預處理

一般狀況下，所獲得的數(shù)據(jù)都具有相關性，所以尋常都要求對數(shù)據(jù)進行初步的白化或球化處理，由于白化處理可去除各觀測信號之間的相關性，從而簡化了后續(xù)獨立分量的提取過程，而且，尋常狀況下，數(shù)據(jù)進行白化處理與不對數(shù)據(jù)進行白化處理相比，算法的收斂性較好。

若一零均值的隨機向量Z??Z1,?,ZM?滿足EZZT?I，其中：I為單位矩陣，我

T??們稱這個向量為白化向量。白化的本質在于去相關，這同主分量分析的目標是一樣的。在ICA

中，對于為零均值的獨立源信號S?t???S1?t?,...,SN?t??T，有：

且協(xié)方差矩陣是單位陣cov?S??I，因此，源信號S?t?E?SiSj??E?Si?E?Sj??0,當i?j，

是白色的。對觀測信號X?t?，我們應當尋覓一個線性變換，使X?t?投影到新的子空間后變成白化向量，即：

Z?t??W0X?t?（2.1）

其中，W0為白化矩陣，Z為白化向量。

利用主分量分析，我們通過計算樣本向量得到一個變換

W0???1/2UT

其中U和?分別代表協(xié)方差矩陣CX的特征向量矩陣和特征值矩陣?？梢宰C明，線性變換W0滿足白化變換的要求。通過正交變換，可以保證UTU?UUT?I。因此，協(xié)方差矩陣：EZZT?E??1/2UTXXTU??1/2???1/2UTEXXTU??1/2???1/2???1/2?I（2.2）再將X?t??AS?t?式代入Z?t??W0X?t?，且令W0A?A,有

Z?t??W0AS?t??AS?t?（2.3）

由于線性變換A連接的是兩個白色隨機矢量Z?t?和S?t?，可以得出A一定是一個正交

??????~~~~變換。假使把上式中的Z?t?看作新的觀測信號，那么可以說，白化使原來的混合矩陣A簡化成一個新的正交矩陣A。證明也是簡單的：EZZT?EASSTAT?AESSTAT?AAT?I（2.4）其實正交變換相當于對多維矢量所在的坐標系進行一個旋轉。

在多維狀況下，混合矩陣A是N?N的，白化后新的混合矩陣A由于是正交矩陣，其

~???~~?~??~~~~自由度降為N??N?1?/2，所以說白化使得ICA問題的工作量幾乎減少了一半。

白化這種常規(guī)的方法作為ICA的預處理可以有效地降低問題的繁雜度，而且算法簡單，

用傳統(tǒng)的PCA就可完成。用PCA對觀測信號進行白化的預處理使得原來所求的解混合矩陣退化成一個正交陣，減少了ICA的工作量。此外，PCA本身具有降維功能，當觀測信號的個數(shù)大于源信號個數(shù)時，經(jīng)過白化可以自動將觀測信號數(shù)目降到與源信號維數(shù)一致。

2.FastICA算法

FastICA算法，又稱固定點(Fixed-Point)算法，是由芬蘭赫爾辛基大學Hyv?rinen等人提出來的。是一種快速尋優(yōu)迭代算法，與普通的神經(jīng)網(wǎng)絡算法不同的是這種算法采用了批處理的方式，即在每一步迭代中有大量的樣本數(shù)據(jù)參與運算。但是從分布式并行處理的觀點看該算法仍可稱之為是一種神經(jīng)網(wǎng)絡算法。FastICA算法有基于峭度、基于似然最大、基于負熵最大等形式，這里，我們介紹基于負熵最大的FastICA算法。它以負熵最大作為一個搜尋

方向，可以實現(xiàn)順序地提取獨立源，充分表達了投影追蹤（ProjectionPursuit）這種傳統(tǒng)線性變換的思想。此外，該算法采用了定點迭代的優(yōu)化算法，使得收斂更加快速、穩(wěn)健。

由于FastICA算法以負熵最大作為一個搜尋方向，因此先探討一下負熵判決準則。由信息論理論可知：在所有等方差的隨機變量中，高斯變量的熵最大，因而我們可以利用熵來度量非高斯性，常用熵的修正形式，即負熵。根據(jù)中心極限定理，若一隨機變量X由大量相互獨立的隨機變量Si?i?1,2,3,...N?之和組成，只要Si具有有限的均值和方差，則不管其為何種分布，隨機變量X較Si更接近高斯分布。換言之，Si較X的非高斯性更強。因此，在分開過程中，可通過對分開結果的非高斯性度量來表示分開結果間的相互獨立性，當非高斯性度量達到最大時，則說明已完成對各獨立分量的分開。

負熵的定義：

Ng?Y??H?YGaus?s?H?Y?（2.5）式中，YGauss是一與Y具有一致方差的高斯隨機變量，H???為隨機變量的微分熵

H?Y???pY???lgpY???d?（2.6）

?根據(jù)信息理論，在具有一致方差的隨機變量中，高斯分布的隨機變量具有最大的微分熵。當Y具有高斯分布時，Ng?Y??0；Y的非高斯性越強，其微分熵越小，Ng?Y?值越大，所以Ng?Y?可以作為隨機變量Y非高斯性的測度。由于根據(jù)式（3.6）計算微分熵需要知道Y的概率密度分布函數(shù)，這顯然不切實際，于是采用如下近似公式：Ng?Y???E?g?Y???E?g?YGauss???2（2.7）

其中，E???為均值運算；g???為非線性函數(shù)，可取g1?y??tanh(a1y)，或這里，尋常我們取a1?1。1?a1?2，g2?y??yexp?y2/2或g3?y??y3等非線性函數(shù)，

快速ICA學習規(guī)則是找一個方向以便WXY?WX具有最大的非高斯性。這里，

T非高斯性用式（3.7）給出的負熵Ng(WX)的近似值來度量,WX的方差約束為1，對于

TT??T?T?白化數(shù)據(jù)而言，這等于約束W的范數(shù)為1。FastICA算法的推導如下。首先，WX的負熵的最大近似值能通過對EGWX??T??進行優(yōu)化來獲得。根據(jù)????Kuhn-Tucker條件，在

EWTX

?????W22?1的約束下，EGWTX的最優(yōu)值能在滿足下式的點上獲得。

EXgWTX??W?0（2.8）這里，?是一個恒定值，

??????EW0TXgW0TX,W0是優(yōu)化后的W值。下面我們利用牛

????頓迭代法解方程（3.8）。用F表示式（3.8）左邊的函數(shù)，可得F的雅可比矩陣JF?W?如下：JF?W??EXXTg'WTX??I（2.9）為了簡化矩陣的求逆，可以近似為（3.9）式的第一項。由于數(shù)據(jù)被球化，EXXT?I,所以，EXXTg'WTX???????????E?XX??E?g'?WX???E?g'?WX??I。因而雅可比矩陣變成了

TTT對角陣，并且能比較簡單地求逆。因而可以得到下面的近似牛頓迭代公式：

?W??W?EXgWTX??W/Eg'WTX??W?W/W??????????????（2.10）

這里，W是W的新值，??EWTXgWTX，規(guī)格化能提高解的穩(wěn)定性。簡化后就可以得到FastICA算法的迭代公式：

????W??EXgWTX?Eg'WTXWW?W/W??????????(2.11)

實踐中，F(xiàn)astICA算法中用的期望必需用它們的估計值代替。當然最好的估計是相應的樣本平均。理想狀況下，所有的有效數(shù)據(jù)都應當參與計算，但這會降低計算速度。所以尋常用一部分樣本的平均來估計，樣本數(shù)目的多少對最終估計的確切度有很大影響。迭代中的樣本點應當分別選取，假使收斂不理想的話，可以增加樣本的數(shù)量。

3.FastICA算法的基本步驟：

1.對觀測數(shù)據(jù)X進行中心化，使它的均值為0；2.對數(shù)據(jù)進行白化，X?Z。

3.選擇需要估計的分量的個數(shù)m，設迭代次數(shù)p?14.選擇一個初始權矢量（隨機的）Wp。

5.令Wp?EZgWpZ?Eg'WpZW,非線性函數(shù)g的選取見前文。6.Wp?Wp???p?1j?1T?????T????WTpWjWj。

7.令Wp?Wp/Wp。

8.假使Wp不收斂的話，返回第5步。9．令p?p?1，假使p?m，返回第4步。

二．MATLAB源程序及說明：

%下程序為ICA的調用函數(shù)，輸入為觀測的信號，輸出為解混后的信號functionZ=ICA(X)

%去均值

[M,T]=size(X);%獲取輸入矩陣的行/列數(shù)，行數(shù)為觀測數(shù)據(jù)的數(shù)目，列數(shù)為采樣點數(shù)

average=mean(X')';%均值fori=1:M

X(i,:)=X(i,:)-average(i)*ones(1,T);end

%白化/球化

Cx=cov(X',1);%計算協(xié)方差矩陣Cx

[eigvector,eigvalue]=eig(Cx);%計算Cx的特征值和特征向量W=eigvalue^(-1/2)*eigvector';%白化矩陣Z=W*X;%正交矩陣

%迭代

Maxcount=10000;%最大迭代次數(shù)Critical=0.00001;%判斷是否收斂

m=M;%需要估計的分量的個數(shù)W=rand(m);forn=1:m

WP=W(:,n);%初始權矢量（任意）%Y=WP'*Z;

%G=Y.^3;%G為非線性函數(shù)，可取y^3等%GG=3*Y.^2;%G的導數(shù)count=0;

LastWP=zeros(m,1);

W(:,n)=W(:,n)/norm(W(