PCA人臉識別理論基礎(chǔ)附源碼

1、1PCA 與人臉識別及其理論基礎(chǔ)1.1 問題描述1對于一幅圖像可以看作一個(gè)由像素值組成的矩陣，也可以擴(kuò)展開，看成一個(gè)矢量，如一 幅 N*N 象素的圖像可以視為長度為 N2 的矢量，這樣就認(rèn)為這幅圖像是位于 N2 維空間中的 一個(gè)點(diǎn)，這種圖像的矢量表示就是原始的圖像空間，但是這個(gè)空間僅是可以表示或者檢測圖 像的許多個(gè)空間中的一個(gè)。不管子空間的具體形式如何，這種方法用于圖像識別的基本思想 都是一樣的，首先選擇一個(gè)合適的子空間，圖像將被投影到這個(gè)子空間上，然后利用對圖像 的這種投影間的某種度量來確定圖像間的相似度，最常見的就是各種距離度量。 K-L 變換1PCA 方法是由 Turk 和 Pentla

2、d 提出來的，它的基礎(chǔ)就是 Karhunen-Loeve 變換(簡稱 KL變換)，是一種常用的正交變換。下面我們首先對 K-L 變換作一個(gè)簡單介紹:假設(shè) X 為 n 維的隨機(jī)變量，X 可以用 n 個(gè)基向量的加權(quán)和來表示:nX = iii =1式中:i 是加權(quán)系數(shù)， i 是基向量，此式還可以用矩陣的形式表示:1 2 n 1 2 nX = ( , ,", )( , ,", )T= 取基向量為正交向量，即T 1i = jT則系數(shù)向量為： j = 0i j j = I = T X綜上所述，K-L 展開式的系數(shù)可用下列步驟求出:步驟一 求隨即向量 X 的自相關(guān)矩陣 R = E X T

3、 X ，由于沒有類別信息的樣本集的 µ 均值向量，常常沒有意義，所以也可以把數(shù)據(jù)的協(xié)方差矩陣 K_L 坐標(biāo)系的產(chǎn)生矩陣，這里 µ 是總體均值向量。= E ( x µ )( x µ )T 作為步驟二 求出自相關(guān)矩陣或協(xié)方差矩陣 R 的本征值 i 和本征向量 i ， = (1 ,i ,",n )步驟三 展開式系數(shù)即為 = T XK_L 變換的實(shí)質(zhì)是建立了一個(gè)新的坐標(biāo)系，將一個(gè)物體主軸沿特征矢量對齊的旋轉(zhuǎn)變換，這個(gè)變換解除了原有數(shù)據(jù)向量的各個(gè)分量之間相關(guān)性，從而有可能去掉那些帶 有較少信息的坐標(biāo)系以達(dá)到降低特征空間維數(shù)的目的。 利用 PCA 進(jìn)行人臉

4、識別完整的 PCA 人臉識別的應(yīng)用包括幾個(gè)步驟：人臉圖像預(yù)處理；讀入人臉庫，訓(xùn)練形成 特征子空間；把訓(xùn)練圖像和測試圖像投影到上一步驟中得到的子空間上；選擇一定 的距離函數(shù)進(jìn)行識別。下面詳細(xì)描述整個(gè)過程（源碼見faceRec.m）。1.讀入人臉庫 歸一化人臉庫后，將庫中的每人選擇一定數(shù)量的圖像構(gòu)成訓(xùn)練集，其余構(gòu)成測試集。設(shè)歸一化后的圖像是 n*m,按列相連就構(gòu)成 N=n*m 維矢量，可視為 N 維空間中的一個(gè)點(diǎn)，可 以通過 K-L 變換用一個(gè)低維子空間描述這個(gè)圖像。2.計(jì)算 K- L 變換的生成矩陣 所有訓(xùn)練樣本的協(xié)方差矩陣為（以下三個(gè)等價(jià)）：1.T TCA = (Mk =1xk ixk ) /

5、 M mx imxA2. C= ( Ai AT ) / M（1） M 3. CA= i =1( x m )( x m )T i x i xA = 1 ,2 , .,M ,i = xi mx , mx 是平均人臉, M 訓(xùn)練人臉數(shù)，協(xié)方差矩陣 CA 是一個(gè)N*N的矩陣, N 是 xi 的維數(shù)。A為了方便計(jì)算特征值和特征向量，一般選用第2個(gè)公式。根據(jù)K - L 變換原理，我們所求 的新坐標(biāo)系即由矩陣 Ai AT 的非零特征值所對應(yīng)的特征向量組成。直接求N*N大小矩陣 C 的 特征值和正交歸一特征向量是很困難的, 根據(jù)奇異值分解原理（見段落1.2.5和1.2.6），可以通過求解 AT i A 的特征

6、值和特征向量來獲得 AT i A 的特征值和特征向量，。N *r 在計(jì)算得到 CA 的所有非零特征值 0 , 1 ,", r 1 （從大到小排序，1 r < M ）及其對應(yīng)的 單位正交特征向量 u0 , u1 ,", ur 1 后，可以得到特征空間 U = u0 , u1 ,", ur 1 ，從而可以計(jì)算一張圖片X在特征空間上的投影系數(shù)（也可以理解為X在空間U中的坐標(biāo)）：3.識別Y = U T * X r *1（2）利用公式（2），首先把所有訓(xùn)練圖片進(jìn)行投影，然后對于測試圖片也進(jìn)行同樣的投影， 采用判別函數(shù)對投影系數(shù)進(jìn)行識別。1.2 PCA 的理論基礎(chǔ) 投影

7、2設(shè) d 維樣本 x1 ， x2 ， " , xn ，以及一個(gè) d 維基 w ，那么標(biāo)量：i iy = wT x是相當(dāng)于 xi 在基上的坐標(biāo)值。如果 w= 1 , yi 就是把 xi 向方向?yàn)?w 的直線進(jìn)行投影的結(jié)果，可以從圖 1 看到。推廣之，如果有一組基（m 個(gè)）組成的空間 W = w1 , w 2 ,", w m ，那么可T m*1以得到 xi 在空間 W 上的坐標(biāo)為： Y = W x 。ix證明： wT x =w又yw x cosx cos = y ,w =1i wT x = y圖 1 投影圖進(jìn)一步，表達(dá)式 w = m + ae 表示 w 是一條通過點(diǎn) m ，方向

8、為 e 的直線。 PCA 的作用及其統(tǒng)計(jì)特性3采用 PCA 對原始數(shù)據(jù)的處理，通常有三個(gè)方面的作用降維、相關(guān)性去除、概率估計(jì)。下 面分別進(jìn)行介紹：去除原始數(shù)據(jù)相關(guān)性從統(tǒng)計(jì)學(xué)上 講， E X E( X )Y E(Y ) 稱為隨 機(jī)變量 X 與 Y 協(xié)方差，記 為Cov( X , Y ) 。令 XY =Cov( X , Y )，稱為隨機(jī)變量 X 與 Y 的相關(guān)系數(shù)。 XY = 1 則 X 與D( X )D(Y )Y 是相關(guān)的， XY = 0 ，則 X 與 Y 是不相關(guān)的。AAT命題 1 對于矩陣 A 來說，如果是一個(gè)對角陣，那么 A 中的向量是非相關(guān)的。由 PCA 處理的人臉庫數(shù)據(jù)的非相關(guān)性可以從

9、兩點(diǎn)進(jìn)行說明。（1） 基底的非相關(guān)性T特征空間基 U = u0 , u1 ,", ur 1 是非相關(guān)的，即 UU= I 。（2） 投影系數(shù)的非相關(guān)性T由 SVD 可知 A = 1 ,2 ,.,M U V, 其中 i = xi mx , mx 是平均人臉。根據(jù)公式（2）可以把 A 映射到特征空間上，得到：關(guān)的，可由下面得到證明：B = U T * A ，其中 B 是非相Y 的協(xié)方差矩陣為： C= 1 BBT =1 U T AAT U =1 2（3）B MMM由命題 1 可知，B 是非相關(guān)的。統(tǒng)計(jì)參數(shù)（均值及方差）均值即 mx -平均人臉。命題 2 隨機(jī)變量方差越大，包含的信息越多，當(dāng)一

10、個(gè)變量方差為 0 時(shí)，該變量為常數(shù)，不 含任何信息。用 PCA 計(jì)算主分量，就是尋找一組向量，使得原始數(shù)據(jù) A = 1 ,2 ,.,M 在這組 向量上的投影值的方差盡可能大。最大方差對應(yīng)的向量就是第一主成份，以后遞推就是第二主成份，第三主成份。用 PCA 計(jì)算主分量就是求原始數(shù)據(jù) A = 1 ,2 ,.,M （其中 i = xi mx ）協(xié)方差T矩陣的特征向量 U = u0 , u1 ,", ur 1 ，由公式（3）可知，P = uiA = ( p1 , p2 ,", pm ) 是 A 在 ui上的投影值，其中 P 的方差就是 ui 對應(yīng)的特征值 i ，可以理解為：命題 3

11、 所有原始數(shù)據(jù)在主分量 ui 上的投影值方差為 i 。降維如果在原始空間表示一幅 n*m 大小的圖片 X，那么需要一個(gè) Nn*m 維矢量，但是 當(dāng)用公式(2)把它映射到特征空間后，只需要一個(gè) r*1 維的向量就可。另外，由命題 2 可知，可以根據(jù)方差的大小來判斷特征向量的重要性。由 ORL 圖片 庫的 200 個(gè)人臉計(jì)算得到的特征值呈圖 2 分布，可知特征向量重要性呈指數(shù)下降，據(jù)此 可以只選用前面幾個(gè)重要的特征向量來構(gòu)建特征空間。通過計(jì)算，前 71 個(gè)特征值占了 90.17%，因此 r 可以取 71 而非 200，從而達(dá)到進(jìn)一 步降維的作用。圖 2 特征值的分布 特征臉U = u0 , u1

12、,", ur 1 中的每一個(gè)單位向量都構(gòu)成一個(gè)特征臉，如圖 3 所示。由這些特征臉?biāo)鶑埑傻目臻g稱為特征臉子空間，需要注意對于正交基的選擇的不同考慮，對應(yīng)較大特 征值的特征向量(正交基)也稱主分量，用于表示人臉的大體形狀，而對應(yīng)于較小特征值 的特征向量則用于描述人臉的具體細(xì)節(jié)，或者從頻域來看，主分量表 示了人臉 的低頻部 分，而此 分量則描 述了人臉 的高頻部 分（ 源 碼 見EigenFace.m）。12105070average圖 3 特征臉，分別是第 1，2，10，50，70 分量，最后一張是平均臉。 圖片重建X要進(jìn)行圖片 X 的重建，首先對 X 投影到特征空間上，得到系數(shù) Y

13、= U T ( X m) ,然后選用一部分系數(shù)與特征向量進(jìn)行原始圖片的重建： X ' = mX + U (1 : t ) *Y (1 : t ) ，其中 1: t 表示取前 t 個(gè)元素。（見reconstruct.m）在圖 4 中，其中前兩張圖片來自訓(xùn)練樣本，第 3 張來自測試樣本，可以看到對于訓(xùn)練樣本，PCA 系數(shù)可以對圖片實(shí)現(xiàn)很好重建，而對于訓(xùn)練樣本以外的圖片重建效果很差。Original1550100150199圖 4 人臉圖像重建。第列張圖片是輸入原始圖，其它列圖片是重建結(jié)果，數(shù)字表示 t 的數(shù)目。 奇異值分解（SVD）1設(shè) A 是秩為 r 的 m*n（m>>n）維

14、矩陣，則存在兩個(gè)正交矩陣和一個(gè)對角陣：TA = a1 , a2 ,", ar = U VT其中 U = u0 , u1 ,", ur 1 ， V = v0 , v1 ,", vr 1 ， = diag (0 , 1 ,", r 1 ) ，且 UU= I ， VV T= I ，2 T m*m Tn*n T Ti 呈降序排列。其中 i為 AA和 A A 的非零特征值，ui 和 vi 分別是 AA 和 A Ai對應(yīng)于 2 的特征向量?？傻靡粋€(gè)推論：U = AV 1可以計(jì)算 AT A 的特征值 2 及相應(yīng)的正交歸一特征向量 v 后，可由推論知 AAT 的正交歸一

15、特i i征向量iu = 1iiAvi注意，協(xié)方差矩陣 CA= ( Ai AT ) / M 的特征值為： 2 / M 。 利用小矩陣計(jì)算大矩陣特征向量高階矩陣的特征向量可以轉(zhuǎn)化為求低階矩陣的特征向量：設(shè)：A 是秩為 r 的 m*n（m>>n）維矩陣， CX= AAT m*m ，是一個(gè)矩陣，現(xiàn)在要求 C 的X特征值及特 征向量，可 通過先求小 矩陣 AT A n*n 的特征 向量 v , v ,", v 和特征 值0 1r 10 , 1 ,", r 1 ，兩者之間有以下關(guān)系：i i iAT A v = vi i i左乘A AAT ( A v ) = ( A v )顯

16、然， C= AAT 的特征向量是 A v （注意沒有單位化）， , ,", 亦為其特征值。X i 0 1r 1結(jié)論： 與 1.2.6 的方法計(jì)算協(xié)方差矩陣的特征向量，特征值的結(jié)果是一致的，只是要注意 中的特征值要除以 M， 中的特征向量要單位化。 圖片歸一化圖片標(biāo)準(zhǔn)化通常是一個(gè)整體概念，要求把圖片歸一到均值為 0，方差為 1 下情況下。 這個(gè)概念類似于一般正態(tài)分布向標(biāo)準(zhǔn)正態(tài)分布的轉(zhuǎn)化：命題 4 若 X N (µ, 2 ) ，則 Z = X µ N (0,1)所以要對一組圖片中的一張 X i 進(jìn)行歸一化（標(biāo)準(zhǔn)化），只需要減去均值，除以方差就可以了。MT均值 mX =

17、 X iiM ，方差為 D = E ( X mX )( X mX ) 1.3 參考文獻(xiàn)1 鄧楠, 基于主成份分析的人臉識別, 西北大學(xué)碩士學(xué)位論文，2006.062 R.O. Duda, P.E. Hart, and D.G. Stork, Pattern Classification, seconded. John Wiley & Sons,2001.3 Sami Romdhani .Face Image Analysis using a Multiple Feature Fitting Strategy. PhD Thesis, University of Basel, Switz

18、erland, January 2005.4 Sami Romdhani .FACE RECOGNITION USING PRINCIPAL COMPONENTS ANALYSIS.1.4 附錄matlab 源碼 人臉識別% FaceRec.m% PCA 人臉識別修訂版，識別率 88% calc xmean,sigma and its eigen decomposition allsamples=;%所有訓(xùn)練圖像for i=1:40for j=1:5 a=imread(strcat('e:ORLs',num2str(i),'',num2str(j),'.

19、jpg');% imshow(a);b=a(1:112*92); % b 是行矢量 1×N，其中 N10304，提取順序是先列后行，即從上 到下，從左到右b=double(b);allsamples=allsamples; b;% allsamples 是一個(gè) M * N 矩陣，allsamples 中每一行數(shù) 據(jù)代表一張圖片，其中 M200end endsamplemean=mean(allsamples);% 平均圖片，1 × Nfor i=1:200 xmean(i,:)=allsamples(i,:)-samplemean;% xmean 是一個(gè) M 

20、15; N 矩陣，xmean每一行保存的數(shù)據(jù)是“每個(gè)圖片數(shù)據(jù)-平均圖片”end;% 獲取特征值及特征向量 sigma=xmean*xmean'% M * M 階矩陣 v d=eig(sigma);d1=diag(d);% 按特征值大小以降序排列dsort = flipud(d1);vsort = fliplr(v);%以下選擇 90%的能量dsum = sum(dsort); dsum_extract = 0; p = 0;while( dsum_extract/dsum < 0.9)p = p + 1;dsum_extract = sum(dsort(1:p);i=1;end%

21、 (訓(xùn)練階段)計(jì)算特征臉形成的坐標(biāo)系base = xmean' * vsort(:,1:p) * diag(dsort(1:p).(-1/2);% base 是 N×p 階矩陣，除以 dsort(i)(1/2)是對人臉圖像的標(biāo)準(zhǔn)化（使其方差為 1）% 詳見基于 PCA 的人臉識別算法研究p31% xmean' * vsort(:,i)是小矩陣的特征向量向大矩陣特征向量轉(zhuǎn)換的過程%while (i<=p && dsort(i)>0)% base(:,i) = dsort(i)(-1/2) * xmean' * vsort(:,i);

22、% base 是 N×p 階矩陣，除以 dsort(i)(1/2)是對人臉圖像的標(biāo)準(zhǔn)化（使其方差為 1）% 詳見基于 PCA 的人臉識別算法研究p31%i = i + 1;% xmean' * vsort(:,i)是小矩陣的特征向量向大矩陣特 征向量轉(zhuǎn)換的過程%end% 以下兩行 add by gongxun 將訓(xùn)練樣本對坐標(biāo)系上進(jìn)行投影,得到一個(gè) M*p 階矩陣 allcoor allcoor = allsamples * base; % allcoor 里面是每張訓(xùn)練人臉圖片在 M*p 子空間中的一個(gè)點(diǎn)， 即在子空間中的組合系數(shù)，accu = 0; % 下面的人臉識別過程

23、中就是利用這些組合系數(shù)來進(jìn)行識別% 測試過程for i=1:40for j=6:10 %讀入 40 x 5 副測試圖像 a=imread(strcat('e:ORLs',num2str(i),'',num2str(j),'.jpg'); b=a(1:10304);b=double(b);tcoor= b * base; %計(jì)算坐標(biāo)，是 1×p 階矩陣for k=1:200mdist(k)=norm(tcoor-allcoor(k,:);end;%三階近鄰dist,index2=sort(mdist);class1=floor( (ind

24、ex2(1)-1)/5 )+1; class2=floor(index2(2)-1)/5)+1; class3=floor(index2(3)-1)/5)+1;if class1=class2 && class2=class3 class=class1;elseif class1=class2 class=class1;elseif class2=class3class=class2;end;if class=i accu=accu+1;end;end;end;accuracy=accu/200 %輸出識別率 特征人臉% eigface.mfunction = eigface(

25、)% calc xmean,sigma and its eigen decomposition allsamples=;%所有訓(xùn)練圖像for i=1:40for j=1:5 a=imread(strcat('e:ORLs',num2str(i),'',num2str(j),'.jpg');% imshow(a);b=a(1:112*92); % b 是行矢量 1×N，其中 N10304，提取順序是先列后行，即從上 到下，從左到右b=double(b);allsamples=allsamples; b;% allsamples 是一個(gè) M

26、 * N 矩陣，allsamples 中每一行數(shù) 據(jù)代表一張圖片，其中 M200end endsamplemean=mean(allsamples);% 平均圖片，1 × Nfor i=1:200 xmean(i,:)=allsamples(i,:)-samplemean;% xmean 是一個(gè) M × N 矩陣，xmean每一行保存的數(shù)據(jù)是“每個(gè)圖片數(shù)據(jù)-平均圖片”end;% 獲取特征值及特征向量 sigma=xmean*xmean'% M * M 階矩陣 v d=eig(sigma);d1=diag(d);dsort = flipud(d1);vsort = f

27、liplr(v);%以下選擇 90%的能量dsum = sum(dsort); dsum_extract = 0; p = 0;while( dsum_extract/dsum < 0.9)p = p + 1;dsum_extract = sum(dsort(1:p);end p = 199;% (訓(xùn)練階段)計(jì)算特征臉形成的坐標(biāo)系%while (i<=p && dsort(i)>0)%base(:,i) = dsort(i)(-1/2) * xmean' * vsort(:,i);% base 是 N×p 階矩陣，除以dsort(i)(1/2

28、)是對人臉圖像的標(biāo)準(zhǔn)化，詳見基于 PCA 的人臉識別算法研究p31%i = i + 1;% xmean' * vsort(:,i)是小矩陣的特征向量向大矩 陣特征向量轉(zhuǎn)換的過程%endbase = xmean' * vsort(:,1:p) * diag(dsort(1:p).(-1/2);% 生成特征臉for (k=1:p),temp = reshape(base(:,k), 112,92); newpath = 'd:test' int2str(k) '.jpg' imwrite(mat2gray(temp), newpath);endavg = reshape(samplemean, 112,92);imwrite(mat2gray(avg), 'd:testaverage.jpg');% 將模型保存save('e:ORLmodel.mat', 'base', 'samplemean'); 人臉重建% Reconstru