基本概念
主成分分析(Principal Component Analysis, PCA)是研究如何將多指標問題轉化為較少的綜合指標的一種重要的統計方法,它能將高維空間的問題轉化到低維空間去處理,使問題變得比較簡單、直觀,而且這些較少的綜合指標之間互不相關,又能提供原有指標的絕大部分信息。
PCA是一個無監督問題,不是基於標簽,而是基於方差。我們可以認為原始的數據在某一個維度上是非常密集的,PCA就是需要把這些密集的點給擴散大,相當於基於方差出找方差最大的方向,越大的方差方向數據點更會分分開,便於分類。
主成分分析除了降低多變量數據系統的維度以外,同時還簡化了變量系統的統計數字特征。主成分分析在對多變量數據系統進行最佳簡化的同時,還可以提供許多重要的系統信息,例如數據點的重心位置(或稱為平均水平),數據變異的最大方向,群點的散布范圍等。
主成分分析作為最重要的多元統計方法之一,其核心是K-L變換,在高維信號壓縮和特征提取等領域已獲得廣泛應用。主成分分析的基本思路可概述如下:借助一個正交變換矩陣,將分量相關的原隨機變量轉換成分量不相關的新變量。從代數角度,即將原變量的協方差矩陣轉換成對角矩陣。從幾何角度,將原變量系統變換成新的正交系統,使之指向樣本點散布最開的正交方向,進而對多維變量系統進行降維處理。按照特征提取的觀點,主成分分析相當於一種基於最小均方誤差的提取方法。
算法實現
用的是MATLAB自帶的一個簡單的數據集fisheriris,該數據集數據類別分為3類,setosa,versicolor,virginica。每類植物有50個樣本,共150個樣本代表150朵花瓣。每個樣本有4個屬性,分別為花萼長,花萼寬,花瓣長,花瓣寬。其中meas是150*4的矩陣代表着有150個樣本每個樣本有4個屬性描述,species代表着這150個樣本的分類。
三種顏色的三點分別對應於三種類型的花瓣,橫坐標和縱坐標分別為護板長度和花瓣寬度,實際我們可以看出對於這兩組數據來說,區分出來比較容易,群體與群體之間的差別比較明顯。這里我們想辦法將二維的特征降低到一維特征,使用matalb自帶的主成分分析函數pca進行主成分分析,使用的命令如下:
[coeff, score, latent, tsquared, explained, mu] = pca(meas(:,3:4));
我們可以看到圖中有紅色和綠色兩條短線,分別是PCA的兩組成分分量,結合圖中我們可以看出數據在紅色主成分分量分布的方差是比較大的,根據這組成分就能很好的將數據進行分類。所以對於花瓣長度和花瓣傳毒這樣一個二維數據而言,對其降維的話可以通過第一組成分分量代表數據集的分布。
使用pca函數得到的參數中第一個參數coeff保存的是主成分分析的主成分在原坐標系中坐標向量的內容,每一列對應一個主成分的分量,每一行對應於原坐標系坐標的分量。計算的coeff結果如下:
0.921777692631943 -0.387718822558475
0.387718822558475 0.921777692631943
第一列中0.921777692631943對應的是紅色短線的方向向量,0.387718822558475對應的是綠色短線的方向向量,該數據不包含中心未知的偏移量,這個信息包含在mu數據中。
Latent是原始數據協方差矩陣的特征值,從中可以看出各個主成分分量的占比,本數據集這個參數的結果如下:
3.66123804559049
0.0360460707406019
可以看出第一主成分對方差的貢獻是最大的。Explained是某個主成分對應方差占整體的百分比,其結果如下:
99.0250662484556
0.974933751544399
也就是說第一主成分占整體方差的99.025%,第二主成分只占0.975%。
然后使用花萼長,花萼寬,花瓣長,花瓣寬進行主成分分析,畫出鳶尾花數據集的PCA散點圖矩陣,其結果如圖所示。
從第一主成分來看可以很好的將數據分開,其數據直方圖混疊的程度比較低,但是成分2、成分3、成分4的散點圖以及矩陣困跌比較嚴重,很難將數據給分開。接下來看一下各個主成分對方差的貢獻率,可以看出第一主成分的貢獻率高達92.46%,第二主成分的貢獻率為5.3%,后面更少,可以判斷通過第主成分可以很好的將數據進行分類。
源代碼
%% 初始化工作空間
clc
clear
close all
%% 載入數據
load fisheriris
%% 二維數據示例
% 花瓣長度和花瓣寬度散點圖
figure;
speciesNum = grp2idx(species);
gscatter(meas(:,3),meas(:,4),speciesNum,['r','g','b'])
% 花瓣長度和寬度的PCA分析
[coeff, score, latent, tsquared, explained, mu] = pca(meas(:,3:4));
hold on
plot(mu(1)+[0 coeff(1,1)],mu(2)+[0 coeff(2,1)],'r')
plot(mu(1)+[0 coeff(1,2)],mu(2)+[0 coeff(2,2)],'g')
hold off
axis equal
axis([-1 7 -1 7])
xlabel('花瓣長度')
ylabel('花瓣寬度')
%% 四維數據示例
%% 主成分分析 Principal Component Analysis
% measRe = score*coeff'+repmat(mu,size(score,1),1);
[coeff, score, latent, tsquared, explained, mu] = pca(meas);
% PLotMatrix繪制散點圖矩陣
% 初始化繪圖窗口
hf = figure;
hf.Units = 'Pixels';
hf.Position = [50 50 800 800];
% 繪制散點圖矩陣
speciesNum = grp2idx(species);
[H,AX,BigAx] = gplotmatrix(score,[],speciesNum,['r','g','b']);
legend(AX(13+3),{'Setosa 山鳶尾','Versicolor 多色鳶尾','Virginica 弗吉尼亞鳶尾'},'Location','northwest','FontWeight','Bold','Fontsize',10)
title(BigAx,'鳶尾花數據PCA散點圖矩陣','FontWeight','Bold','Fontsize',16)
%橫坐標標題
xlabel(AX(4),{'Component 1';'成分1'},'FontWeight','Bold','Fontsize',12)
xlabel(AX(8+1),{'Component 2';'成分2'},'FontWeight','Bold','Fontsize',12)
xlabel(AX(12+2),{'Component 3','成分3'},'FontWeight','Bold','Fontsize',12)
xlabel(AX(16+3),{'Component 4','成分4'},'FontWeight','Bold','Fontsize',12)
% 縱坐標標題
ylabel(AX(1),{'Component 1';'成分1'},'FontWeight','Bold','Fontsize',12)
ylabel(AX(2),{'Component 2';'成分2'},'FontWeight','Bold','Fontsize',12)
ylabel(AX(3),{'Component 3','成分3'},'FontWeight','Bold','Fontsize',12)
ylabel(AX(4),{'Component 4','成分4'},'FontWeight','Bold','Fontsize',12)
% Biplot繪制向量圖
hf2 = figure;
vbls={'萼片長度','萼片寬度','花瓣長度','花瓣寬度'};
biplot(coeff(:,1:3),'Score',score(:,1:3),'Varlabels',vbls);