實驗三、主成分分析

1、實驗三、主成分分析一、實驗?zāi)康恼莆栈陔x散K-L變換做特征提取的原理和特點，編寫基于離散K-L變換的特征提取算法程序。二、實驗內(nèi)容編寫基于類間離差矩陣Sb的DKLT特征提取程序，要求：（1）在程序注解中，列寫出程序的輸入、輸出、算法描述、變量說明； （2）對c類問題，求輸入樣本x的類間離差矩陣Sb；（3）求Sb的本征值12d及其對應(yīng)的本征矢量T=(t1,t2,td)，d<c； （4）取d=2，做變換y=Tx，觀察y的分布情況；（5）選一種聚類算法（例如c-均值算法），分別對x和y進(jìn)行聚類，對比分析其結(jié)果。三、實驗原理、方法和手段在實際問題的研究中，往往會涉及眾多有關(guān)的變量。但是，變量太多

2、不但會增加計算的復(fù)雜性，而且也會給合理地分析問題和解釋問題帶來困難。一般說來，雖然每個變量都提供了一定的信息，但其重要性有所不同，而在很多情況下，變量間有一定的相關(guān)性，從而使得這些變量所提供的信息在一定程度上有所重疊。因而人們希望對這些變量加以“改造”，用為數(shù)極少的互補(bǔ)相關(guān)的新變量來反映原變量所提供的絕大部分信息，通過對新變量的分析達(dá)到解決問題的目的。 總體主成分1）定義設(shè) X1，X2，Xp 為某實際問題所涉及的 p 個隨機(jī)變量。記 X=(X1，X2，,Xp)T，其協(xié)方差矩陣為=(ij)pp=E(X-E(X)(X-E(X)T,它是一個 p 階非負(fù)定矩陣。設(shè)TY1=l1X=l11X1+l12X2

3、+ +l1pXpTY2=l2X=l21X1+l22X2+ +l2pXp（1） Y=lTX=lX+lX+ +lXpp11p22pppp則有Var(Yi)=Var(liTX)=liTli,i=1,2,.,p,TCov(Yi,Yj)=Cov(liTX,lTjX)=lilj,j=1,2,.,p.（2）第 i 個主成分：一般地，在約束條件liTli=1及Cov(Yi,Yk)=liTlk=0,k=1,2,.,i-1.下，求 li 使 Var(Yi)達(dá)到最大，由此 li 所確定的Yi=liTX稱為 X1，X2，Xp 的第 i 個主成分。2） 總體主成分的計算設(shè) 是X=(X1,X2,.,Xp)T的協(xié)方差矩陣，

4、的特征值及相應(yīng)的正交單位化特征向量分別為12 p0及e1,e2,.,ep,則 X 的第 i 個主成分為Yi=eiTX=ei1X1+ei2X2+ +eipXp,i=1,2,.,p, （3） 此時TVar(Yi)=eiei=i,i=1,2,.,p,TCov(Yi,Yk)=eiek=0,ik.3）總體主成分的性質(zhì) 主成分的協(xié)方差矩陣及總方差記 Y=(Y1,Y2,.,Yp)T 為主成分向量，則 Y=PTX，其中P=(e1,e2,.,ep)，且Cov(Y)=Cov(PTX)=PTP=Diag(1,2,.,p),由此得主成分的總方差為Var(Y)=ii=1i=1ppi=tr(PP)=tr(PP)=tr()

5、=Var(Xi),TTi=1p即主成分分析是把 p 個原始變量 X1，X2，Xp 的總方差Var(X)ii=1p分解成 p 個互不相關(guān)變量 Y1，Y2，Yp的方差之和，即Var(Y)ii=1p而 Var(Yk)=k。 第 k 個主成分的貢獻(xiàn)率：ii=1p；i前m個主成分累計貢獻(xiàn)率：i=1i=1pmi，它表明前 m 個主成分Y1，Y2，Ym綜合提供 X1，iX2，Xp中信息的能力。 主成分 Yi 與變量 Xj 的相關(guān)系數(shù)由于 Y=PTX，故 X=PY，從而Xj=e1jY1+e2jY2+ +epjYp,Cov(Yi,Xj)=ieij.由此可得 Yi 與 Xj 的相關(guān)系數(shù)為Y,X=ijCov(Y,X

6、)=e=ij （4）4）標(biāo)準(zhǔn)化變量的主成分在實際問題中，不同的變量往往有不同的量綱，由于不同的量綱會引起各變量取值的分散程度差異較大，這時總體方差則主要受方差較大的變量的控制。為了消除由于量綱的不同可能帶來的影響，常采用變量標(biāo)準(zhǔn)化的方法，即令Xi*=,i=1,2,.,p, （5）其中 i=E(Xi),ii=Var(Xi). 這時*TX*=(X1,X2,.,X*)p的協(xié)方差矩陣便是X=(X1,X2,.,Xp)T的相關(guān)矩陣 =(ij)pp，其中ij=E(Xi*X*j)=利用 X 的相關(guān)矩陣 作主成分分析，有如下結(jié)論：Cov(X,X). （6）*T設(shè) X*=(X1為標(biāo)準(zhǔn)化的隨機(jī)向量，其協(xié)方差矩陣（即

7、 X 的相關(guān)矩陣）為 ,X2,.,X*)p ，則 X*的第 i 個主成分為X-*Yi*=(ei*)TX*=ei*1+ei*2+ +eip,i=1,2,.,p. （7）并且Var(Yi=1p*i)=Var(Xi*)=p, （8）*ii=1*ii=1*i1i2pp*T的特征值，e=(e,e,.,eip其中 12 *)為相應(yīng)于特征值 i*的正p0為交單位特征向量。i*第 i 個主成分的貢獻(xiàn)率：；p前 m 個主成分的累計貢獻(xiàn)率：i=1m*ip；*。 Yi*與Xi*的相關(guān)系數(shù)為Y*,X*=ijij（2）樣本主成分前面討論的是總體主成分，但在實際問題中，一般 （或）是未知的，需要通過樣本來估計。設(shè)xi=(

8、xi1,xi2,.,xip)T,i=1,2,.,n.為取自X=(X1,X2,.,Xp)T的一個容量為n的簡單隨機(jī)樣本，則樣本協(xié)方差矩陣及樣本相關(guān)矩陣分別為1nS=(sij)pp=(xk-)(xk-)T,n-1k=1（9） s,R=(rij)pp=其中1nT=(1,2,.,p),j=xij,j=1,2,.,p,ni=11nsij=(xki-i)(xkj-j),i,j=1,2,.,p.n-1k=1分別以 S 和 R 作為 和的估計，然后按總體主成分分析的方法作樣本主成分分析。四、實例編程實現(xiàn)某市為了全面分析機(jī)械類個企業(yè)的經(jīng)濟(jì)效益，選擇了8個不同的利潤指標(biāo)，14企業(yè)關(guān)于這8個指標(biāo)的統(tǒng)計數(shù)據(jù)如下表所示

9、，試進(jìn)行主成分分析。表1 14家企業(yè)的利潤指標(biāo)的統(tǒng)計數(shù)據(jù)分析：樣本均值向量為：=(27.97910.9509.1008.54311.06414.6141.55214.686)T,樣本協(xié)方差矩陣為：168.33360.35745.75741.21557.90671.6728.602101.62037.20716.82515.50523.53529.0294.78544.02324.84324.33536.47849.2783.62939.41024.42336.28349.1463.67538.718S=56.04675.4045.00259.723 103.0186.82174.5231.13

10、76.722102.707168.3360.35745.75841.21657.90671.6728.602101.6260.35737.20716.82515.50523.53529.0294.784644.02345.75816.82524.84324.33536.47849.2783.62939.4141.21615.50524.33524.42336.28349.1463.674738.718S=57.90623.53536.47836.28356.04675.4045.002259.723 71.67229.02949.27849.14675.404103.026.821574.52

11、38.6024.78463.6293.67475.00226.82151.1376.7217101.6244.02339.4138.71859.72374.5236.7217102.71由于S中主對角線元素差異較大，因此我們樣本相關(guān)矩陣R出發(fā)進(jìn)行主成分分析。樣本相關(guān)矩陣R為：1 0.76266 0.70758 0.64281 0.59617 0.54426 0.62178 0.772851 0.553410.51434 0.51538 0.468880.73562 0.7121410.98793 0.9776 0.974090.68282 0.78019 1 0.98071 0.97980.6

12、9735 0.77306R= 1 0.992350.62663 0.78718 0.6303 0.72449 1 1 0.62202 1矩陣R的特征值及相應(yīng)的特征向量分別為：R的特征值及貢獻(xiàn)率見下表前 前3個標(biāo)準(zhǔn)化樣本主成分中各標(biāo)準(zhǔn)化變量 xi*=對應(yīng)特征向量，由此得到3個標(biāo)準(zhǔn)化樣本主成分為*y1=0.32113x*+0.29516x+0.38912x+0.38472x+0.37955x+0.37087x+0.31996x+0.35546x12345678*y2=-0.4151x1-0.59766x2+0.22974x3+0.27869x4+0.31632x5+0.37151x6-0.2781

13、4x7-0.15684x8*y=-0.45123x+0.10303x-0.039895x+0.053874x-0.037292x+0.075186x+0.77059x-0.42478x123456783(i=1,2,.,8)前的系數(shù)即為注意到，y1近似是8個標(biāo)準(zhǔn)化變量xi*=(i=1,2,.,8)的等權(quán)重之和，是反映各企業(yè)總效應(yīng)的綜合指標(biāo)，y1的值越大，則企業(yè)的效益越好。由于y1的貢獻(xiàn)率高達(dá)76.708%，故若用y1的得分值對各企業(yè)進(jìn)行排序，能從整體上反映企業(yè)之間的效應(yīng)差別。將S中sii的值及中各i的值以及各企業(yè)關(guān)于xi的觀測值代入y1的表達(dá)式中，可求得各企業(yè)y1的得分及其按其得分由大到小的排

14、序結(jié)果。所以，第9Matlab程序：coeff,score,latent=princomp(X)注：該函數(shù)使用協(xié)方差陣作主成分分析。主成分分析程序a=;b=corrcoef(zscore(a)%計算相關(guān)系數(shù)矩陣D=tril(b)%得到三角矩陣d,v=eig(b)%計算特征值和特征向量y1=zscore(a)*d(:,7)%計算第一主成分?jǐn)?shù)值f1,i1=sort(y1);f2,i2=sort(i1);flipud(i1),flipud(f1),f2%第一主成分得分排序y2=zscore(a)*d(:,6)%計算第二主成分?jǐn)?shù)值f1,i1=sort(y2);f2,i2=sort(i1);flipud

15、(i1),flipud(f1),f2%第二主成分得分排序y3=zscore(a)*d(:,5)%計算第三主成分?jǐn)?shù)值f1,i1=sort(y3);f2,i2=sort(i1);flipud(i1),flipud(f1),f2%第三主成分得分排序y4=zscore(a)*d(:,4)%計算第三主成分?jǐn)?shù)值f1,i1=sort(y4);f2,i2=sort(i1);flipud(i1),flipud(f1),f2%第si主成分得分排序y5=zscore(a)*d(:,3)%計算第一主成分?jǐn)?shù)值f1,i1=sort(y1);f2,i2=sort(i1);flipud(i1),flipud(f1),f2%第一主成分得分排序y=y5*(0