ICA独立成分分析—FastICA基于负熵最大
1. 概念
官方解釋:利用統計原理進行計算的方法,是一種線性變換。
?ICA分為基于信息論準則的迭代算法和基于統計學的代數方法兩大類,如FastICA算法,Infomax算法,最大似然估計算法等。
這里主要討論FastICA算法。
先來看ICA和PCA的區別:PCA是降維并提取不相關屬性,而ICA是降維并提取相互獨立的屬性(不相關不一定獨立,獨立一定不相關。不相關是指沒有線性關系,獨立是指沒有任何關系)。PCA是提取出最能表示原始事物的特征,而ICA是使每個分量最大化獨立,便于發現隱藏因素。
PCA的適用環境是數據為高斯分布時。而ICA不適用于高斯分布的數據。
ICA的兩條假設:
① 源信號之間互相獨立
② 每一個源信號為非高斯分布
【更新日志:2019-12-27】
為什么原信號為非高斯分布,戳這里
http://cis.legacy.ics.tkk.fi/aapo/papers/IJCNN99_tutorialweb/node9.html
混合信號源的三個特性:
①獨立性(Independence):源信號相互獨立,但是混合的信號卻不是,因為混合信號共享了源信號。
②正態化(Normality):根據中心極限定理,具有有限方差的獨立隨機變量的和的分布傾向于高斯分布;寬松點說就是,兩個獨立隨機變量的和比原始的獨立隨機變量中的任何一個更加接近高斯分布。這里我們考慮的每一個信號源都是獨立隨機變量。
③復雜度(Complexity):任何混合信號的時間復雜性都比構成它的最簡單的信號源復雜。
其實這幾點就告訴我們ICA的動機:如果組成混合信號的信號是獨立的,或者具有非高斯直方圖,或者較低復雜度,那么他們一定是信號源。
2. 白化
為了便于計算,需要先進行數據預處理:
對數據進行白化或球化處理去除各觀測信號之間的相關性,簡化后續獨立分量的提取過程,且算法的收斂性較好。
白化向量:若一零均值的隨機向量Z=(Z1,Z2,…ZM)T滿足E{ZZT}=I,其中I是單位矩陣。白化的本質是去相關,同PCA的目標是一樣的。稱源信號S(t)為白色的是因為對于零均值的獨立源信號S(t)=[S1(t),….,SN(t)]T,有E{Si Sj}=E{Si}E{Sj}=0(i≠j)且協方差矩陣是單位陣,cov(S)=I。
對于觀測信號X(t),找到線性變換W0,使X(t)投影到新的子空間后變成白化向量,即Z(t)=W0X(t),
其中W0為白化矩陣,Z為白化向量
利用主成分分析能得到
???????????????W0=Λ-1/2UT
其中U和Λ分別代表協方差矩陣CX的特征向量矩陣和特征值矩陣。
因此協方差矩陣:
E{ZZT}=E{Λ-1/2UTXXTUΛ-1/2}=Λ-1/2UTE{XXT}UΛ-1/2=Λ-1/2ΛΛ-1/2=I
將X(t)=AS(t)式代入Z(t)=W0X(t),且令W0A=A~,有
Z(t)=W0AS(t)=A~S(t)
多維情況下,混合矩陣A是N*N的,白化后的新的混合矩陣A~是正交矩陣,自由度降為N*(N-1)/2,所以說白化使得ICA問題的工作量幾乎減少了一半。
用PCA對觀測信號進行白化的預處理使原來所求的解混合矩陣退化成一個正交陣,減少了ICA的工作量,當觀測信號的個數大于源信號個數時,經過PCA降維也就是白化可以自動將觀測信號的數目降到與源信號數目維數相同。
?3.?FastICA算法
FastICA算法也叫固定點算法(Fixed-Point)算法,是一種快速尋優迭代算法,采用批處理的方式,每一步迭代由大量的樣本數據參與運算。
FastICA由基于峭度,基于似然最大,基于負熵最大等形式。這里介紹基于負熵最大的FastICA算法。
負熵的判決準則:由信息論理論可知,在所有等方差的隨機變量中,高斯變量的熵最大,所以可以利用熵來度量非高斯性,采用熵的修正形式,負熵。根據中心極限定理,若一個隨機變量X由許多相互獨立的隨機變量Si(i=1,2,3….,N)之和組成,只要Si具有有限的均值和方差,則不論其為何種分布,隨機變量X較Si更接近高斯分布。所以當高斯性度量達到最大的時候,說明完成各獨立成分的分離。
負熵的定義
Ng(Y)=H(YGauss)-H(Y)
其中YGauss是與Y具有相同方差的高斯隨機變量。H(.)為隨機變量的微分熵。
?
根據信息理論,在具有相同方差的隨機變量中,高斯分布的隨機變量具有最大的微分熵。當Y具有高斯分布時,Ng(Y)=0;Y的非高斯性越強,其微分熵越小,Ng(Y)的值越大,所以Ng(Y)可以作為隨機變量Y非高斯性的測度。采用負熵定義求解需要知道Y的概率密度分布函數,但是實際不可能,于是采用下面的近似公式:
Ng(Y)={E[g(Y)]-E[g(YGauss)]}2
其中E[.]為均值運算,g(.)為非線性函數,可取g1(y)=tanh(a1y)或g2(y)=y exp(-y2/2)或g3(y)=y3等非線性函數,這里1≤a1≤2,通常取a1=1
快速ICA的規則就是找到一個方向以便WTX(Y=WTX)具有最大的非高斯性,非高斯性用Ng(Y)={E[g(Y)]-E[g(YGauss)]}2給出的負熵的近似值來度量。WTX的方差約束為1,對于白化數據,等于約束W的范數為1.
FastICA的推導:
① WTX的負熵的最大近似值能通過對E{G(WTX)}進行優化取得。在E{( WTX)2}=||W||2=1的約束下,E{G(WTX)}的最優值能在滿足下式的點上獲得
E{Xg(WTX)}+βW=0
其中β=E{W0TXg(WTX)} 是一個恒定值,W0是優化后的W值。
②利用牛頓迭代法解①的方程。用F表示左邊的函數,得到F的雅克比矩陣JF(W)如下:
JF(W)=E{XXTg’(WTX)}-βI ??可以近似為第一項,即忽略βI
由于數據被球化,所以E{XXT}=I,所以E{XXTg’(WTX)}≈E{XXT}*E{g’(WTX)}=E{g’(WTX)}I。
從而雅克比矩陣變成了對角陣,并且比較容易求逆。因而得到下面的近似牛頓迭代公式:
?
這里的W*是W的新值,β= E{WTXg(WTX)},規格化能提高穩定性。
簡化后得到FastICA的迭代公式:
?
實踐中,FastICA算法中用的期望必須用他們的估計值代替。最好的估計是相應的樣本平均。理想情況下,所有的有效數據都應該參與計算,但是會降低運算速度,所以通常選取一部分樣本的平均來估計,樣本數目的多少對最后估計的精確度有很大影響。迭代中的樣本點應該分別選取,加入收斂不理想,可以增加樣本數量。
FastICA算法的步驟:
1.??????對觀測數據X進行中心化,使它的均值為0
2.??????對數據進行白化,X→Z
3.??????選擇需要估計的分量的個數m,設迭代次數p←1
4.??????選擇一個初始權矢量(隨機的)Wp。
5.??????令Wp=E{Zg(WTZ)}-E{g’(WTZ)}W,非線性函數g,可取g1(y)=tanh(a1y)或g2(y)=y exp(-y2/2)或g3(y)=y3等非線性函數
?? 6.
?
7.令Wp=Wp/||Wp||。
8. 假如Wp不收斂的話,返回第五步
9. 令p=p+1,如果p≤m,返回第四步
ICA.m
?
function Z = ICA( X )%去均值 [M,T]=size(X); %獲取輸入矩陣的行列數,行數為觀測數據的數目,列數為采樣點數 average=mean(X')'; %均值 for i=1:MX(i,:)=X(i,:)-average(i)*ones(1,T); end%白化/球化 Cx=cov(X',1); %計算協方差矩陣Cx [eigvector,eigvalue]=eig(Cx); %計算Cx的特征值和特征向量 W=eigvalue^(-1/2)*eigvector'; %白化矩陣 Z=W*X; %正交矩陣%迭代 Maxcount=10000; %最大迭代次數 Critical=0.00001; %判斷是否收斂 m=M; W=rand(m); for n=1:mWP=W(:,n); %初始權矢量(任意)%Y=WP'*Z;%G=Y.^3;%G為非線性函數,可取y^3等%GG=3*Y.^2; %G的導數count=0;LastWP=zeros(m,1);W(:,n)=W(:,n)/norm(W(:,n));while abs(WP-LastWP)&abs(WP+LastWP)>Criticalcount=count+1; %迭代次數LastWP=WP; %上次迭代的值%WP=1/T*Z*((LastWP'*Z).^3)'-3*LastWP;for i=1:mWP(i)=mean(Z(i,:).*(tanh((LastWP)'*Z)))-(mean(1-(tanh((LastWP))'*Z).^2)).*LastWP(i);endWPP=zeros(m,1);for j=1:n-1WPP=WPP+(WP'*W(:,j))*W(:,j);endWP=WP-WPP;WP=WP/(norm(WP));if count==Maxcountfprintf('未找到相應的信號');return;endendW(:,n)=WP; end Z=W'*Z; end?
ICATest.m
?
clear all; clc; N=200; n=1:N;%N為采樣本數 s1=2*sin(0.02*pi*n); %正弦信號 t=1:N; s2=2*square(100*t,50); %方波信號 a=linspace(1,-1,25); s3=2*[a,a,a,a,a,a,a,a];%鋸齒信號 s4=rand(1,N); %隨機噪聲 S=[s1;s2;s3;s4]; %信號組成4*N A=rand(4,4); X=A*S; %觀察信號%源信號波形圖 figure(1); subplot(4,1,1);plot(s1);axis([0 N -5,5]);title('源信號'); subplot(4,1,2);plot(s2);axis([0 N -5,5]); subplot(4,1,3);plot(s3);axis([0 N -5,5]); subplot(4,1,4);plot(s4);xlabel('Time/ms'); %觀察信號(混合信號)波形圖 figure(2); subplot(4,1,1);plot(X(1,:));title('觀察信號(混合信號)'); subplot(4,1,2);plot(X(2,:)); subplot(4,1,3);plot(X(3,:)); subplot(4,1,4);plot(X(4,:));Z=ICA(X);figure(3); subplot(4,1,1);plot(Z(1,:));title('分離信號'); subplot(4,1,2);plot(Z(2,:)); subplot(4,1,3);plot(Z(3,:)); subplot(4,1,4);plot(Z(4,:)); plot(Z(4,:)); xlabel('Time/ms');
結果
?
?
?
參考文獻:
http://www.cnblogs.com/tornadomeet/archive/2012/12/30/2839841.html
http://wenku.baidu.com/link?url=2gqaP9JugoT41p-9yvDRsZe6easW_NSNjhIakgt9lAR74GiYfZ2-lXBDu57DAYe8U8FZN5Q9uhBeO2icr32KfYuDER6IzruR9rh6smmY7gW
https://en.wikipedia.org/wiki/Independent_component_analysis
http://cs229.stanford.edu/notes/cs229-notes11.pdf
本文已經同步到微信公眾號中,公眾號與本博客將持續同步更新運動捕捉、機器學習、深度學習、計算機視覺算法,敬請關注
總結
以上是生活随笔為你收集整理的ICA独立成分分析—FastICA基于负熵最大的全部內容,希望文章能夠幫你解決所遇到的問題。
- 上一篇: 剑网3 PVP花间对战藏剑小技巧分享 花
- 下一篇: 热播职场剧惊现《逆水寒》情节 端游玩家看