A tutorial on Principal Components Analysis | 主成分分析（PCA）教程

A tutorial on Principal Components Analysis

原著：Lindsay I Smith, A tutorial on Principal Components Analysis, February 26, 2002.

翻譯：houchaoqun.
時間：2017/01/18.
出處：http://blog.csdn.net/houchaoqun_xmu  |  http://blog.csdn.net/Houchaoqun_XMU/article/details/54601984

說明：本文在翻譯原文的過程當中，在保證不改變原文意思的前提下，對一些知識點進行了擴充並附上參考網址。本人剛開始嘗試翻譯文獻，經驗不足，有些理解可能不到位，有些翻譯可能不夠準確。若讀者們有發現有誤之處，還望海涵並加以指正，本人會及時加以指正。

• 本文中涉及到公式的編輯，使用的是「在線LaTeX公式編輯器」，將編輯完成後的公式（圖片的形式）複製到博文中便可。
• 原文可能缺乏一些圖，本人經過Matlab做圖並補充。
• 本人在原文的基礎上，額外提供了一些參考網址，有興趣的讀者們能夠詳細閱讀。
• 爲了更讀者更深入的理解PCA，本人還另外提供了案例，並用Matlab實現，僅供讀者們參考。
• 爲了更好的將理論運用到實際，本人打算在介紹完這篇PCA教程後，整理一篇基於PCA進行人臉識別的具體案例，因爲水平有限，須要改正和完善的地方但願讀者們提供寶貴的意見。

--------------------------------

主成分分析（PCA）教程

第1章：簡介（Introduction）

    本教程的目的是爲了讓讀者理解主成分分析（PCA）。PCA是一項頗有用的統計學技術，已經應用在人臉識別（Face Recognition）和圖像壓縮（Image Compression）等領域，而且做爲一項通用技術在高維數據中發現數據模型（降維技術）。

    開始介紹PCA以前，本教程首先介紹一些學習PCA相關的數學概念。包括，標準差，協方差，特徵向量及其對應的特徵值。這部分背景知識旨在使得PCA相關的章節更加簡潔明瞭，若是讀者對這部分的概念很熟悉，能夠選擇跳過這部份內容。

    本教程的一些相關案例旨在說明這些數學概念，讓讀者對這部份內容理解得更加深入。若是讀者想要更進一步的學習，「初等線性代數 第五版（Elementary Linear Algebra 5e）」這本數學工具書是一個很好的關於這些數學背景的資源。這本書的做者是Howard Anton，由John Wiley和Sons股份有限公司出版（國際標準圖書編號[ISBN] 0-471-85223-6）。

 

第2章：數學背景（Background Mathematics）

    本章將給出一些對於理解PCA原理（涉及相關的推導過程）所需的基礎數學技能。這部分知識都是相互獨立的，而且都有相關的例子。理解爲何可能使用到這個技術以及在數據方面，某個操做的結果告訴了咱們什麼，比記住準確的數學技術更重要。雖然並非這些技術都運用到PCA，可是最重要的技術（運用到PCA的技術）正是基於那些沒有明確須要的技術（沒有明確運用到PCA），因此對於初學者，仍是有必要好好理解一下各個知識點。

    本章包含兩大部分，第一部分介紹了關於統計學的知識，側重於介紹數據（分佈測量）是如何分佈的（標準差，方差，協方差以及協方差矩陣）。第二部分介紹了矩陣代數，側重於特徵向量和特徵值，矩陣的一些重要性質是理解PCA的基礎。

2.1 統計學（Statistics）

    關於統計學的整個主題是圍繞着你有一個大的數據集，而且你想要分析這個數據集中點跟點之間的關係。本章將側重介紹一些讀者可以親自操做（手算驗證）的測量指標（標準差，方差，協方差以及協方差矩陣），而且對於數據自己，這些測量指標告訴了咱們什麼（即，測量值反映數據的關係）。

2.1.1 標準差（Standard Deviation）

    爲了理解標準差，咱們須要一個數據集，統計學家一般使用來自整體（population）中的一個樣本（sample）。此處使用選舉投票（election polls）爲例，「整體」是這個國家的全部人，而統計學家測量（measure）的樣本是這個整體的一個子集。統計學的偉大之處在於僅僅經過測量（能夠經過電話調查或者相似的方法）整體中的一個樣本，就可以算出與整體最類似的測量結果（即，經過測量樣本進而解決整體）。在統計學這個章節，做者假設本文使用的數據集是來自更大整體的樣本（整體的數據集比樣本的數據集更大）。（本章）下文有一個關於樣本和整體更詳細信息的參考網址。

    例如，對於數據集X：

    咱們簡單的使用符號X表示整個數據集。若是想表示數據集中的某個元素，可使用下標表示指定的元素（值），如   表示數據集X的第三個元素，對應的值是4；  指的是這個序列中的第一個元素，值爲1；可是這裏並無像咱們在一些教科書上提到的 （由於此處的下標從1開始）。此外，符號n表示數據集X的元素個數。

    咱們能夠對一組數據集進行不少種計算。例如，咱們能夠計算這個樣本（數據集）的平均值（mean）。做者假設讀者們理解樣本X的均值的含義，此處僅給出計算公式，以下：

    注：符號   表示數據集X的均值。上述公式表示，「累加數據集X中的全部元素的值，而後除以元素總數n」。

   不幸的是，均值僅僅告訴咱們數據集的中間點。例如，如下兩組數據集具備相關的均值（mean = 10），可是這是兩組徹底不一樣的數據。

    所以，這兩組數據集有什麼不一樣呢？數據的分佈不一樣。數據集的標準差反映了數據是如何分佈的（個體間的離散程度）。

    咱們要怎樣計算標準差（Standard Deviation）呢？標準差（SD）的英文定義是，「The average distance from the mean of the data set to a point（數據集中每一個點到均值的平均距離）」。標準差的計算方式是，首先計算每一個點到樣本均值的距離的平方，而後進行累加，最後將累加的結果除以 n-1 再開根號，以下公式所示：

    其中，一般使用s表示一個樣本的標準差。讀者可能會問，「爲何是除以 (n-1) ，而不是 n？」。答案有些複雜，但總的來講，若是你的數據集是樣本數據集，也就是說若是你選擇的是一個真實世界中的子集（例如，對樣本數爲500我的的選舉投票案例進行測量），你就須要使用 (n-1) 。由於事實證實，對於樣本，使用(n-1)獲得的答案更接近於標準差的結果（若是你嘗試使用Matlab提供的cov函數，也是除以n-1）。若是你處理的數據集是整體，那你應該使用n。（ 其餘解釋參考[ 點擊進入連接 ]：如是整體，標準差公式根號內除以n，如是樣本,標準差公式根號內除以（n-1)；這樣能使咱們以較小的樣本集更好的逼近整體的標準差，即統計上所謂的「無偏估計」）。無論怎樣，若是你計算的不是樣本的標準差，並且整體的標準差，那麼公式根號內應該除以n，而不是(n-1)。對於這部份內容，有意願更進一步瞭解的讀者能夠參考網址： http://mathcentral.uregina.ca/RR/database/RR.09.95/weston2.html。這個網址用類似的方式介紹了標準差，還提供了一個實驗案例，解釋了標準差公式根號內的兩種分母（即，n和n-1）的區別，並對「樣本」和「整體」的區別進行了討論。

    以下表2.1所示的兩組數據，分別計算了它們的標準差。

   正如預期的，表2.1的第一部分的標準差（8.3266）更大是由於這組樣本數據相對於均值的分佈更離散。正如另外一個例子所示，樣本數據爲[ 10 10 10 10 ]，均值一樣爲10，可是它的標準差爲0。由於全部元素的值都同樣，不存在偏離均值的元素。

2.1.2 方差（Variance）

    方差是另外一種測量數據集離散程度的方式。實際上，方差和標準差幾乎是相同的。公式以下：

    讀者們能夠發現方差僅僅是對標準差等式左右兩邊進行平方（方差的公式中沒有對其開根號）。對於樣本中的方差   是經常使用的符號。這兩個衡量值（方差和標準差）都是衡量數據的離散程度。標準差是最經常使用的衡量值，但方差一樣適用。本文在介紹了標準差以外，還介紹方差的緣由是爲下一個章節的協方差提供基礎（鋪墊，拋磚引玉）。

    練習：以下三組數據集所示，分別求出均值，標準差和方差。

• [ 12 23 34 44 59 70 98 ]
• [ 12 15 25 27 32 88 89 ]
• [ 15 35 78 82 90 96 97 ]

2.1.3 協方差（Covariance）

    上述咱們關注的兩種（標準差和方差）測量指標的對象僅僅是一維數據集（侷限性），好比房間裏全部人的身高或者教師批改試卷等。然而，許多數據集都不僅是一維，並且對於這些數據集，統計分析的目標一般是尋找這些維度之間是否存在任何關係。例如，咱們將班級上全部學生的身高以及學生考試的分數做爲數據集，而後使用統計分析方法去判斷學生的身高是否對他們的分數有任何影響。

    標準差和方差僅僅處理一維數據，所以你只能分別對數據集中的每一維計算出標準差（每一維對應一個標準差）。然而，找到一種與上述類似，而且可以度量各個維度偏離其均值程度的測量指標是很是有用的。

    【協方差】就是這樣的一種測量指標。協方差一般用來測量2維的數據集。若是你計算的協方差中的2個維度都是同一個變量，此時獲得的就是方差。（協方差正常處理的是二維數據，維數大於2的時候，就用協方差矩陣表示）所以，若是有一組三維數據集（x, y, z），那麼咱們須要分別計算x維度和y維度，x維度和z維度，y維度和z維度的協方差。計算x維度和x維度，y維度和y維度，z維度和z維度的協方差就至關於分別計算x，y，z維度的方差（這就是兩個維度都是同一個變量的狀況）。

    協方差的公式和方差的公式很類似，方差的公式（度量各個維度偏離其均值的程度）以下所示：

    其中，做者僅僅是拓展了平方項變成兩個部分（如，  表示成  *  ）。所以，協方差的公式以下：

    上述兩個公式（方差和標準差）除了分子中的第二個括號內將X用Y代替覺得，其餘都同樣。用英文定義就是，「For each data item, multiply the different between the x value and the mean of x, by the difference between the y value and the mean of y. Add all these up, and divide by (n -1)」（對於X和Y的每一個數據項，  與 （變量X的均值）的差乘以   與  （變量Y的均值）的差，而後將它們累加起來後除以 (n-1) ）。

    那麼，協方差是如何運做的呢？讓咱們使用一些實例數據。想象咱們走進世界中，收集到一些二維數據。咱們統計一羣學生用在學習  COSC241 這個課程的時間（小時）以及他們得到的分數。如今咱們有了二維的數據，第一維是H維，表示學生學習的時間（以小時爲單位），第二維是M維，表示學生得到的分數。以下圖所示，包含了咱們剛剛想象的數據以及H和M之間的協方差  的計算結果。

圖示：學生學習的時間H，得到的分數M

    那麼，這些數據告訴咱們什麼？協方差的符號（正好或者負號）比數值更加劇要。若是協方差的值是正數，代表這兩維（H和M）是正相關，也就是說，隨着學生學習的時間（H）的增長，學生得到的分數（M）也增長。若是值是負數，一個維度增長則另外一個維度減少（數據呈負相關）。對於這個例子，若是協方差最終獲得的結果是負值，也就是負相關，意味着隨着學生學習時間的增長，學生得到的分數則減小。最後還有一種狀況是，當協方差的值爲0，意味着這兩維（如，變量H和M）是相互獨立的（不相關）。

    實驗結果（上述數據）能夠很直觀的經過對原始數據做圖，以下圖所示，學生所得的分數M隨着學生學習的時間H的增長而增長。然而，數據的可視化（做圖）也僅僅侷限於2維和3維的數據。由於對於任意2維的數據均可以計算協方差的值，所以這項技術（協方差）常常用於發現高維數據之間的關係（這些高維數據是很難進行可視化的）。

    讀者們或許會問，「   等於   嗎？」。經過協方差公式告訴咱們，答案是「   等於   」，（展開協方差公式）這二者惟一的區別就是式子   用式子    替代。根據「乘法交換律」可知，兩個因數相乘，交換因數的位置，積不變，也就是說 和 獲得相同的結果。

2.1.4 協方差矩陣（The Covariance Matrix）

    回顧協方差的相關知識，它一般測量的是二維數據。假如咱們有一組大於2維的數據，那麼將有多個協方差能夠計算（有多種組合形式的協方差）。例如，對於一個三維的數據（x, y, z），咱們可以計算 ,   和  。實際上，對於一組n維數據，咱們能夠計算   種不一樣組合的協方差。一種有用的方式是，「計算不一樣維度之間全部可能狀況的協方差值，而後將它們放入一個矩陣」，構成協方差矩陣。（本教程假設讀者們對矩陣很熟悉而且知道如何去定義矩陣）。所以，一組n維數組的協方差矩陣的定義以下：

其中， 表示一個n行n列的矩陣，  表示第x維。這一串複雜的公式（ugly looking formula）的意思是，若是你有一組n維的數據集，那麼獲得的是一個n行n列的協方差矩陣，矩陣中的每個元素值是經過計算兩個不一樣維度的協方差獲得的（矩陣的位置跟協方差的組合有對應關係）。例如，第2行，第3列位置的元素值是經過計算第2維和第3維之間的協方差獲得的。

    舉一個例子，咱們將虛構一組的3維數據集的協方差矩陣（使用x, y, z表示這三維）。獲得一個3行3列的協方差矩陣，對應的元素值以下所示：

    一些關鍵的注意點：一個注意點是，沿着協方差矩陣的主對角線方向，你會發現元素值是某一維和它自己的協方差，也就是這個維的方差。另外一個注意點是，由 可知，協方差矩陣是一個對稱陣。（協方差還有另外一種定義形式，參考 百度百科）

   練習：

    以下所示的一組2維數據集，計算出x和y維之間的協方差，並根據計算結果說明一下數據。

    計算以下所示的一組3維數據集的協方差矩陣。

2.2 矩陣代數（Matrix Algebra）

    本章節提供矩陣代數的背景知識有助於下文對PCA的理解。做者將側重的介紹對於給定一個矩陣（方陣）的特徵向量和特徵值。一樣的，做者假設讀者們已經瞭解了一些矩陣的基礎知識。

2.2.1 特徵向量（Eigenvectors）

   如讀者們所知，假如兩個矩陣的行數和列數知足矩陣相乘的條件（例如，C=AB，當矩陣A的列數等於矩陣B的行數時，A與B能夠相乘），就能夠將他們相乘。特徵向量就是其中的一個特例，考慮以下2個「矩陣乘以向量」的式子。

式1：沒有特徵向量

式2：有一個特徵向量

    在第一個乘法式子（式1）裏，矩陣乘以向量的結果並非形如一個整數乘以原始向量（沒有特徵向量）。而在第二個乘法式子（式2）裏相乘獲得的結果是4乘以原始向量 。爲何是這樣子呢？緣由以下：原始向量 是一個2維空間的向量，（在x-y座標系下）向量 表示一個從原點（0, 0）指向（3, 2）的一個箭頭。式2左邊的另外一個乘法因子是一個方陣（行數等於列數的矩陣）能夠表示爲一個變換矩陣。若是將這個變換矩陣左乘以（左乘和右乘有區別）一個向量，所得的結果是由原始狀態（向量）經過變換矩陣獲得的另外一個向量。

    特徵向量來自於變換（矩陣）的性質（ 線性變換的特徵向量是一個非簡併的向量，其方向在該變換下不變。該向量在此變換下縮放的比例稱爲其特徵值）。（A * V）想象有這麼一個變換矩陣A做爲左乘因子，將向量V映射到直線 y=x 上。而後你能夠注意到是否有一個向量坐落於直線方程 y=x 上，這就是它自身的映射。此時，向量V（能夠任意擴大它的倍數，由於向量的長度並不影響結果）就是變換矩陣A的一個特徵向量。

    特徵向量還有其餘什麼性質呢？讀者們首先應該知道特徵向量必須是由一個方陣求解獲得。但也不是全部的方陣都有特徵向量。對於給定一個 n * n 的矩陣，而且該矩陣存在特徵向量，那麼這個矩陣有n個特徵向量。相似的，對於給定一個 3 * 3 的矩陣（假設存在特徵向量），那麼這個矩陣有3個特徵向量。

   特徵向量的另外一個性質是，即便變換矩陣A乘以原始向量V =  以前，先將V擴大2倍，再將矩陣與擴大後的向量 V' 相乘，獲得的結果仍然是變換矩陣的特徵向量，且一樣擴大了2倍，以下圖所示。

向量與矩陣相乘前，乘以倍數2

變換矩陣乘以擴大了2倍的原始向量

    這是由於，若是你對原始向量擴大n倍，意味着將向量變得更長但不改變向量的方向。最後，由變換矩陣獲得的全部特徵向量都是正交的，無論數據集有多少維，兩兩都相互垂直。這部分知識很是重要，由於這意味着你能夠經過這些正交的特徵向量來數據表達，而不用根據 x-y 座標系。咱們將在下文的PCA章節進行討論。（參考： 基與正交基,特徵值與特徵向量）

    另外一個讀者們應該知道的重點是，「數學家們在求解特徵向量時，他們更願意去找到這些模正好爲1的特徵向量」。這是由於，如讀者們所知，向量的大小並不影響這個向量是否爲特徵向量，可是方向就有影響。所以，爲了標準化特徵向量，當咱們求得一個特徵向量時，一般將它的模縮放爲1，使得全部特徵向量都有相同的模（模爲1）。以下，經過一個例子加以驗證。（此部分更詳細的內容可參考： PCA數學原理）

向量   是一個特徵向量，這個向量的模爲： ，而後咱們將原始向量   除以模 ，使得特徵向量的模變爲1，以下：

    那麼，咱們如何着手去找到這些神祕的特徵向量呢？不幸的是，可以簡單的進行求解的條件是你處理的對象是一個很小的矩陣，就像一個規模不大於 3*3 的方陣。對於求解規模較大方陣的特徵向量，一般的方法是經過一些複雜的迭代方法進行求解（這些方法超出了本教程的範圍，做者沒作進一步的介紹）。若是讀者們須要經過程序求解方陣的特徵向量，能夠找一個數學庫（例如，Matlab，你只須要調用相關的函數便可，具體的實現過程這個庫已經幫你完成了）。此外，有一個叫作「newmat」的數學工具包僅供讀者參考： http://webnz.com/robert/。

    想要進一步瞭解「關於如何求解通常狀況下的特徵向量，正交性」的讀者們，能夠參考由Howard Anton編著，John Wiley & Sons股份有限公司出版的教科書「基礎線性代數 第五版」（ISBN 0-471-85223-6）。

2.2.2 特徵值

    特徵值跟特徵向量緊密相關，實際上，讓咱們回顧一下「2.2.1 特徵向量」這節的2個例子，咱們能夠注意到這2個例子中，「相同的方陣A，分別乘以原始向量   以及擴大2倍後的向量  ，最後獲得的特徵向量前面的係數是同樣的」。在這2個例子中，獲得的係數都是4，4就是這個這個特徵向量對應的特徵值。在方陣A乘以向量V以前，不管咱們先對向量V乘上一個多大的係數（非零：線性代數中規定特徵向量不能夠爲零向量）獲得的V'，最後獲得的都會是4乘以 V' 。（ 特徵值的定義：設A爲n階矩陣，若存在常數λ及n維非零向量x，使得Ax=λx，則稱λ是矩陣A的特徵值，x是A屬於特徵值λ的特徵向量）

    讀者們能夠發現特徵向量和特徵值是成對出現的。當你使用一個理想的程序庫去計算特徵向量時，同時你也會獲得相應的特徵值。

   練習：

    對於下列矩陣：

    分別對下列的5個向量，判斷是否爲上述方陣的特徵向量，若是是，求出該特徵向量向量對應的特徵值。

，  ，  ，  ， 

    【補充】能夠經過Matlab提供的「eig()」函數求解特徵向量及其對應的特徵值，以下所示：

close all;
clc;
clear;

A = [3 0 1; -4 1 2; -6 0 -2];
[d, v] = eig(A)

A * d - v * d

    如上實驗結果所示，方陣A求得2個特徵向量 [0; 1; 0] 和[0; -1; 0] 對應的特徵值都是1。而上述例子中所給的5個向量中，只有[0;1; 0]是方陣A的特徵向量。值得注意的是：對於上述例子的求解方法，應該根據特徵值和特徵向量的定義進行求解，即根據是否知足「A*Vector = egienvalue*Vector」的形式，若是知足，則向量Vector是方陣A的特徵向量，反之不是。

第3章：主成分分析（PCA）

    最後，咱們來到了主成分分析（PCA）這節。主成分分析是個什麼東西呢？它是一種「識別數據中的模式，經過強調數據的異同（similarities and differences）的方式來表達數據」的方法。因爲很難發現高維數據中的模式，對於圖形表示法是不可用的，而PCA正是一種強大的數據分析工具。

    PCA的另外一個主要優點是，一旦你發現數據中的這些模式，你就能夠經過減小維數來壓縮數據（保留主成分的數據），並且不會損失太多的信息。這項技術（優點）運用在圖像壓縮上，後面章節會具體介紹。

    本章將爲讀者們介紹「在一組數據集上實現主成分分析」所需的相關步驟。本文並不會詳細介紹PCA這項技術爲何能起做用（下文中，做者沒有具體給出計算過程，是直接給出結果），可是做者會嘗試着提供「每一個要點（步驟）會發生什麼（有什麼做用）」的解釋，所以當讀者們嘗試着使用到PCA這項技術時，可以作出明確的決定。

3.1 理論方法（PCA的六大步驟）

步驟一：獲取數據集

    爲了簡單起見，本文將使用做者本身構造的數據集。這是一組2維的數據集，選擇這樣一組數據集是由於可以對其做圖，更直觀的顯示PCA分析過程當中的每個步驟的效果。以下左圖所示是本章的原始數據集Data和減去均值後的數據集DataAdjust，右圖是對DataAdjust進行做圖的效果圖。

         

步驟二：減去均值

    爲了讓PCA正常工做，必須將原始數據集中的每一維的全部元素值減去該維的均值（數據的標準化）。（此處的數據集是2維數據，即x-y）也就是說，全部的 值都減去 （ 就是x維上全部元素值的均值），全部的   值都減去  。這個過程使得咱們獲得一組均值爲0的數據集。

步驟三：計算協方差矩陣

    這個步驟正如2.1.4節討論過的計算協方差矩陣。由於這是一組2維數據集，因此將會獲得一個 2 X 2 的協方差矩陣。在這裏也是同樣（無心外），所以做者直接給出結果，以下：

    根據上述結果，因爲協方差矩陣   非對角線上的元素值是正數，咱們能夠推出變量 x 和 y 是正相關（同增同減）。

步驟四：計算協方差矩陣的特徵向量及其對應的特徵值

    由於協方差矩陣是一個方陣，因此咱們能夠計算協方差矩陣的特徵向量及其對應的特徵值，這一步驟是很是重要的。特徵向量及其對應的特徵值告訴咱們關於數據有用的信息，待會兒將展現給讀者們看。在此期間，計算所得的特徵向量和特徵值以下所示：

    ，     

注意：原文中給的特徵向量的符號跟Matlab求得的不同（這裏給出Matlab求得的解），基於Matlab求解特徵值和特徵向量的代碼以下所示：

cov = [0.616555556, 0.61544444; 0.616555556, 0.716555556];  % 協方差矩陣
[eigenvectors, eigenvalues] = eig(cov);  % 求解特徵向量及其對應的特徵值

值得注意的是，此處求解獲得的這些特徵向量都是單位特徵向量，即 ，這兩個特徵向量的模都爲1。這一點對於PCA是很是重要的，幸運的是，大多數數學工具包求解特徵向量時，獲得的結果都是單位特徵向量。
    那麼，這又是什麼意思呢？以下圖所示爲標準化後的數據集（各維變量減去對應的均值），讀者們能夠發現這些數據中存在一個明顯的模式（quite a strong pattern）。正如協方差矩陣所預期的同樣，x和y兩組變量的值呈正相關（同增同減）。在圖的上方，做者也畫出了對應的兩個特徵向量所表示的直線方程，它們以虛線的形式出現於座標系的對角線。正如前面提到的同樣，這一組特徵向量相互垂直（正交），分別是圖中斜率大於0以及斜率小於0的這兩條虛線。但更重要的是，特徵向量爲咱們提供了數據模式（此處的模式就是直線方程）的相關信息。首先讓咱們先關注這條通過數據集（標準化後）中心的特徵向量對應的直線（斜率小於0那條直線），看起來像畫了一條最擬合的直線嗎（很明顯不是）？這個特徵向量告訴咱們這兩組數據集（應該指的是這兩維的數據，x-y）是如何關聯到這條直線。（這就是咱們要找的直線，數據模式）而第二個特徵向量爲了咱們提供了另外一條直線（另外一種數據模式），斜率大於0的那條線，座標上全部的點沿着這條主線分佈，但都以必定程度的距離偏離在這條主線的兩側。

    因此，經過求解協方差矩陣的特徵向量，咱們已經可以提取出描述數據集的直線。如今，咱們剩下的步驟包括，對數據進行轉化，使其可以表達在特徵向量對應的直線上。

步驟五：選擇主成分，造成特徵向量

    接下來就是「數據壓縮和降維」的到來。若是讀者們仔細觀察前幾章講到的特徵向量和特徵值，注意到求解獲得的特徵值是不同的。實際上，最大的特徵值對應的特徵向量就是數據集的主成分。在本文的例子裏，最大特徵值對應的特徵向量就是x-y座標系中，斜率大於0的那條直線（此處論文好像有點問題，讀者們可查看原文進一步驗證）。這就是表達數據集維度之間的關係最有意義，最具表明性的模式。

    一般來講，一旦從協方差矩陣中求解獲得特徵向量，下一步就是從大到小對特徵值進行排序（此處要對應好特徵向量和特徵值），也就是對數據成分主次的排序（特徵值越大，數據在該維就越有意義，信息量越大）。如今，若是讀者們有需求的話，就能夠忽略那些成分較小的維度。雖然你丟失了一些數據信息，可是忽略的這些維度的特徵值很小，丟失的信息不多。若是你忽略一些成分（維度），最終獲得的數據集的維數將比原始數據集的維數少（這就是數據降維）。（To be precise）確切地講就是，若是你處理的原始數據集是n維，你將計算獲得n個特徵向量和特徵值（有對應關係），而後你選擇特徵值前p（p < n）大的特徵向量，忽略掉其餘維度的數據，最後獲得的數據集就是一組p維的數據集。

    至此，讀者們如今須要作的事情就是構造一個特徵向量（fearue vector），其實就是由向量構成的一個矩陣。  就是取前k大的特徵值對應的eigenvectors，而後根據以下形式，將這k個向量以列的形式構成一個矩陣。

    回到本文例子的數據集，經過求解，咱們有2個特徵向量   ，  ，意味着咱們有2種選擇。以下所示，咱們能夠同時選擇這2個特徵向量做爲 ：

  

    咱們也能夠選擇忽略特徵值較小的那維數據，僅僅保留特徵值較大的那維數據，以下所示：

    做者將在下一節展現以上兩種結果。

（注：以上的 FeatureVector ，下文運用的時候，須要對其先進行轉置）

步驟六：產生新的數據集（降維後）

    這是PCA的最後一個步驟，也是最簡單的。一旦咱們選好了主成分（eigenvectors），並構造了特徵向量（ ）後。咱們首先對  進行轉置（獲得1行2列或者2行2列的 ），而後左乘以原始數據（已標準化且轉置後，即2行n列），以下所示：

其中，  表示提取前k維主成分後構成的   並進行轉置後的矩陣（  是按列進行存儲，通過轉置後的  變成按行進行存儲）；  表示原始數據減去均值後再進行轉置獲得的結果，  矩陣的每一列表示一組具體的數據項，如 ，每一行表示獨立的一個維（如，x維或者y維）。讀者可能會由於這個步驟中出現的忽然「轉置」感到困惑，做者爲此感到很抱歉，但若是讀者們先對 和 進行轉置，這裏的給出的公式將更加簡單。而不是在 和 的右上方表上一個轉置符號T（其實，翻譯這裏的時候，我不是很清楚這句話的用意）；  表示通過PCA算法後最終獲得的數據集，其中每一列表示數據項 ，每一行表示具體的維（x維或者y維）。

    行文至此，（降維後的   ）咱們將得到什麼呢？僅僅根據咱們選擇的前k維特徵向量，咱們就能夠獲得原始數據。案例中的數據集有2維，包括x維和y維。讀者們能夠選擇本身喜歡的兩個維度來表示這些數據。若是這些維度是相互垂直（perpendicular）的，那麼以這種方式表達數據會更有效。這就是爲何特徵向量須要兩兩互相垂直會這麼重要。咱們已經改變了數據的表示形式（原來是基於x-y座標系），如今是表示在2個特徵向量上。當新的數據集進了降維的狀況下，咱們已是忽略掉一些特徵向量（ ），只保留了咱們選擇的前k個特徵向量（ ）。

    爲了在數據上展現這些，做者已經對每一種可能的特徵向量（ ）完成了最後的變換（  左乘 ）。同時做者也對每一種狀況的結果進行了轉置，將數據轉換回原始漂亮的表格格式（就是將更漂亮，更直觀的表示格式）。做者還對最後獲得的新數據進行做圖，體現這些數據是怎樣關聯到成分（提取的前k個主成分）。

• 第一種狀況（保留了全部的特徵向量）：以下左圖表示通過 變換後的 ，右圖表示在新的座標系（）下繪製 這些數據座標點，能夠理解爲原始數據和座標系通過旋轉必定的角度獲得的。不難理解，這種狀況下，沒有損失任何數據。

• 第二種狀況（只選擇最大特徵值對應的特徵向量）：以下圖所示爲降維後的數據集（從2維降到1維，和預期的同樣）。讀者們若是將這組數據與第一種狀況的數據做比較，你會發現這組數據就是第一種狀況那組數據的第一列，即對應 x 維的數據。

    因此，若是你想對這組數據做圖，以下圖所示，獲得的就是一維空間，並且圖上 y=0 直線附近的這些座標點剛好是 x 維空間上的位置（此時獲得的座標系，至關於y=0，主成分都在x軸上的數據）。咱們已經拋掉另外一維空間信息（y維），也就是另外一個特徵向量（值較小的那個特徵向量被拋棄了）。

    那麼，到這裏咱們完成了什麼呢？咱們基本上已經對數據進行了變換，使得數據表達成它們之間的模式，本例中這些模式是「描述這些數據之間的關係的最擬合的直線」。這是很是有用的，由於咱們如今已經將數據「根據這些直線（也就是特徵向量eigenvector）對原始數據貢獻信息的程度」進行組合分類。最初，咱們有 x 和 y 座標軸（完整的數據信息），但每一個數據點的 x 值和 y 值實際上並不能確切告訴咱們這個數據跟其餘數據之間的關係（y值提供的信息量不多）。如今（降維後），數據點的座標值（只有 x 值）告訴咱們這個數據點適合於哪條趨勢線（trend lines）。本例中，兩個eigenvector都使用（轉換）的狀況下，咱們僅僅簡單的改變了數據，用求解獲得的兩個特徵向量（eigenvectors）代替本來的 x-y 座標系而已。可是，將特徵值較小（y維，貢獻小）的那一維去掉，只保留與特徵值較大的那一維（x維）相關的數據的這種狀況，就達到很好的降維效果，並且在丟失少許數據信息的狀況下，用一維就可以很好的表示原始數據。

3.1.1 恢復舊數據

    若是讀者們使用 PCA 進行數據壓縮，必定都會想要「恢復原始數據」，這個問題想必受到不少讀者的關注，畢竟數據通過PCA變換後，還是丟失了一些信息（下文有提供相關的案例）。這部份內容摘抄於 http://www.vision.auc.dk/ sig/Teaching/Flerdim/Current/hotelling/hotelling.html。（這是原文提供的網址，可是好像進不去了- -）

    那麼，咱們該如何恢復原始數據呢？在這以前，讀者們應該知道只有在「將全部特徵向量（eigenvectors）共同構成最終的變換矩陣」這種狀況下，咱們才能精確的恢復原始數據。若是最終的變換矩陣已是通過降維（丟掉一些特徵向量，好比上個例子中的 y 軸被丟棄），那麼原始數據已是丟失掉一些信息了，儘管是少許的信息。

    回顧一下以前的「最終變換矩陣」，以下：

    爲了恢復原始數據，咱們對上述公式作一下變形，獲得以下公式：

其中，  是   的逆（inverse）。當咱們將全部的eigenvectors做爲feature vector時（由於eigenvector和feature vector都翻譯成「特徵向量」，因此此處，我直接使用英文表示，加以區分），你會發現，原來矩陣 feature vector 的逆剛好等於feature vector的轉置（當且僅當知足「feature vector矩陣的元素是數據集對應的單位特徵向量，eigenvectors」這個條件的時候，成立。證實過程能夠參考， 這裏）。這使得恢復數據變得更加容易，由於上述的表達式獲得了改善，以下所示：

可是，爲了恢復到最原始數據，咱們還須要對每一個維度的元素都加上相應的均值（由於，咱們一開始的時候就對原始數據作了標準化，也就是各維度的元素都減去了均值）。因此，以下所示，是更完整的表達式：

這個公式對「並無把全部的 eigenvectors 做爲feature vector」這種狀況也一樣適用。因此即便當你丟棄一些eigenvector時，上述等式仍然能夠獲得正確的轉換。

    做者在本文中並不打算使用完整的feature vector（使用全部的eigenvector）演示數據恢復，由於這樣獲得的結果剛好就是咱們一開始處理的數據集（最原始的數據）。然而，做者將作以下實驗：「僅保留特徵值較大的維度空間的狀況，展現數據是如何丟失的」，以下圖所示（降維後的數據），將圖中的數據與最原始的數據做比較，你會發現只有沿着主成分eigenvector方向上的變量（x維）被保留了下來，沿着另外一個eigenvector方向上的變量（y維）已經不見了。

    練習：

• 求解協方差矩陣獲得的特徵向量可以帶給咱們什麼（有什麼做用）？
• 在PCA流程中，咱們可以在「哪一個步驟」決定壓縮數據？會有什麼影響？
• 找一個關於「PCA運用到面部識別」的案例並用圖示法表示主要的eigenvectors（前k大的特徵值對應的eigenvectors），能夠經過「Eigenfaces」關鍵字進行搜索。

第4章：計算機視覺領域的應用（Application to Computer Vision）

    本章將對「PCA在計算機視覺中的應用」進行概述。首先介紹圖像一般是怎樣表示的，而後介紹PCA可以幫助咱們對這些圖像作什麼樣的處理。在本章中，關於「面部識別（facial recognition）」的相關信息是來自於「Face Recognition: Eigenface, Elastic Matching, and Neural Nets」, Jun Zhang et al. Proceedings of the IEEE, Vol. 85, No. 9, September 1997。 表徵信息來自於「Digital Image Processing」, Rafael C. Gonzalez and Paul Wintz, Addison-Wesley Publishing Company, 1987（這一樣是關於「通常狀況下， K-L變換更進一步的知識」一篇的很好的參考文獻）。圖像壓縮的相關資料來自於「 http://www.vision.auc.dk/ sig/Teaching/Flerdim/Current/hotelling/hotelling.html」（又進不去0- -，原文是2002年的，可能有些網址已經發生了改變，後續有找到的話，再更新分享給你們），這個參考網址也爲咱們提供一些關於「使用不一樣數量的eigenvector的圖像重建」的案例。

4.1 （圖像）表示方法（Representation）

    在計算機視覺領域，當咱們使用矩陣方法時，必定要考慮圖像的表示形式。一個 N*N 的方陣能夠被表示成一個  維的向量，以下所示：

其中，圖像中的像素以行爲單位，一行接着一行被放置，造成一個一維的圖像。例如，前N個元素（  - ）將被做爲圖像的第一行，接下來的N個元素（  - ）做爲第二行，以此類推。向量中的值反應圖像的強度信息，或許就是一個簡單的灰度值（對於灰度圖像，就是對應灰度值）。

補充，以下矩陣所示，表示一個圖像進行數字化的矩陣表示形式：

4.2 基於PCA尋找模式（模型）（PCA to find patterns）

    假設咱們有20張圖像，每張圖像由N個像素的高和N個像素的寬組成（N*N的矩陣）。對於每張圖像，咱們可使用上一節的方法將其表示爲一個圖像向量。而後，咱們能夠將全部的圖像（如今一張圖像對應一個向量）放到一個大矩陣裏面，形如：

    以此做爲咱們使用PCA算法的第一步（如今處理的對象是一個由20張圖像構成的大矩陣）。一旦使用PCA，咱們要作的就是求解協方差矩陣獲得特徵向量（eigenvectors）。這（PCA）爲何有用呢？假設咱們想實現面部識別，原始的數據集是人臉圖像。問題是，給定一張新的人臉圖像，識別出這張新圖像對應原始數據圖像中的哪一類人臉圖像，也就是根據人臉圖像信息進行分類（注意，這張新圖像並非來自於咱們一開始所給的那20張圖像）？這個問題在計算機視覺的解法是，基於PCA分析獲得的新座標系下，測量新的人臉圖像與已知的20張人臉圖像之間的區別，而不是在原來的座標系。

    事實證實通過PCA算法獲得的新的座標系更有利於識別人臉，這是由於PCA算法告訴咱們原始圖像（數據）之間的差別（differences and similarities）。PCA算法肯定了數據中的統計模型。

    由於全部的向量都是  維，因此咱們將獲得 個特徵向量（eigenvector）。實際上，咱們也能夠丟棄一些意義不大的eigenvectors（只保留特徵值前k大對應的eigenvectors），識別的效果一樣不錯。

4.3 基於PCA的圖像壓縮（PCA for image compression）

    使用PCA算法進行圖像壓縮又稱爲Hotelling 變換或者K-L變換。假如咱們有20張圖像，每張 個像素。咱們能夠構造 個向量，每一個向量20維，每一維對應這20張圖像中同一個像素的強度值，下文我（博主）將補充說明。這一點與上一個例子的大不一樣，上一個例子是構成的向量vector中的每一個元素都是對應不一樣的像素，而如今這個例子構成的每一個向量（ 個）的元素對應20張圖像相同的一個像素值。

---------------

補充說明： 個20維的向量

• 第1個向量形如：
• 第2個向量形如：
• 第個向量形如：

其中，每一個向量中有20個像素值，好比第1個向量  的20個元素值分別對應從第一張圖像到第20張圖像的第一個像素值 ，以此類推。

---------------

    如今開始對這組數據使用PCA算法。咱們將獲得20個特徵向量（eigenvectors），由於每一個vector都是20維。爲了壓縮數據，咱們可使用特徵值前15大對應的15個eigenvector做爲變換矩陣。最後獲得的數據集只剩下15維，減小了1/4的空間。然而，但須要重構原始數據集時，獲得的圖像將丟失一些信息。PCA是一項有損信息的壓縮技術，由於解壓後的圖像並不跟原始數據徹底同樣，一般更糟（差別更大）。

5. 博主補充：PCA實例

• 原始數據集：，使用PCA算法將這組二維數據降維成一維數據。
• 由於這個矩陣的每行已是零均值，這裏咱們直接求協方差矩陣，過程以下：

• 而後求其特徵值和特徵向量，具體求解方法再也不詳述，能夠參考相關資料。求解後特徵值以下所示：

• 其對應的特徵向量分別是：

•  其中對應的特徵向量分別是一個通解， 和  可取任意實數。那麼標準化後的特徵向量（是一個單位特徵向量）爲 ：

• 所以咱們的矩陣P是：

• 能夠驗證協方差矩陣C的對角化：

• 最後咱們用P的第一行乘以數據矩陣，就獲得了降維後的表示：

• 降維投影結果以下圖：

注：本例（來自於文章「 PCA數學原理」）在求解協方差矩陣的過程當中，協方差的分母是用n表示。可是咱們在本教程中提到，處理樣本而非整體的時候，咱們更傾向於選擇 (n-1) ，即 。在此，咱們保留以上例子，但爲了讓例子更貼近本教程，博主使用Matlab對該例子進行求解，協方差的分母使用的是 (n-1)，代碼和實驗結果以下所示：

• 數據集仍是上個例子的數據集：
• 求解協方差矩陣，代碼以下：

close all;
clc;
clear;

%% PCA案例
X = [-1,-1,0,2,0];
Y = [-2,0,0,1,1];
mean_X = mean(X);
mean_Y = mean(Y);
n = 5;
sum_xy = 0;
sum_yx = 0;
sum_xx = 0;
sum_yy = 0;

for i = 1:5
    sum_xy = sum_xy + (X(i) - mean_X) * (Y(i) - mean_Y);
    sum_yx = sum_yx + (Y(i) - mean_Y) * (X(i) - mean_X);
    sum_xx = sum_xx + (X(i) - mean_X) * (X(i) - mean_X);
    sum_yy = sum_yy + (Y(i) - mean_Y) * (Y(i) - mean_Y);
end

cov_xx = sum_xx / (n-1)
cov_xy = sum_xy / (n-1)
cov_yy = sum_yy / (n-1)
cov_yx = sum_yx / (n-1)

    -- 也可使用Matlab提供的cov函數直接獲得協方差矩陣：

 

       -- 協方差矩陣：

• 根據協方差矩陣求解特徵向量和特徵值：

cov_test = cov(X,Y);    % 協方差矩陣
[eigenvectors eigenvalues] = eig(cov_test);  %特徵向量和特徵值

       -- 特徵向量（經過Matlab的eig求得的特徵向量已是標準化的）：

       -- 特徵值：0.5 和 2.5

• 選擇最大特徵值（此處爲2.5）對應的特徵向量做爲feature vector：
• 根據公式 ，求解降維後的新數據集：（從2維降到1維）

X = [-1,-1,0,2,0];
Y = [-2,0,0,1,1];

% 特徵向量
featurevector = [0.7071;0.7071];     % 直接提取最大特徵值對應的特徵向量
RowFeatureVector = featurevector';   % 轉置
RowAdjustData = [X;Y];
FinalData = RowFeatureVector * RowAdjustData

附錄A

實現代碼：

    本文代碼是基於Scilab（一款開源的軟件）實現的。與Matlab相似， SCILAB也是一種科學工程計算軟件，其數據類型豐富，能夠很方便地實現各類矩陣運算與圖形顯示，能應用於科學計算、數學建模、信號處理、決策優化、線性/非線性控制等各個方面。做者使用如下代碼生成了本文中的全部案例。代碼中除了第一部分的宏（macro）之外，其他部分都是做者本身完成的。

    注：這部分代碼，本人還沒有進一步驗證，後續驗證了再更新。

%  This macro taken from
%  http://www.cs.montana.edu/˜harkin/courses/cs530/scilab/macros/cov.sci
%  No alterations made
%  Return the covariance matrix of the data in x, where each column of x
%  is one dimension of an n-dimensional data set. That is, x has x columns
%  and m rows, and each row is one sample.
% 
%  For example, if x is three dimensional and there are 4 samples.
%  x = [1 2 3;4 5 6;7 8 9;10 11 12]
%  c = cov (x)
function [c]=cov (x)
    % Get the size of the array
    sizex=size(x);
    % Get the mean of each column
    meanx = mean (x, 'r');
    % For each pair of variables, x1, x2, calculate
    % sum ((x1 - meanx1)(x2-meanx2))/(m-1)
    for var = 1:sizex(2),
        x1 = x(:,var);
        mx1 = meanx (var);
        for ct = var:sizex (2),
            x2 = x(:,ct);
            mx2 = meanx (ct);
            v = ((x1 - mx1)' * (x2 - mx2))/(sizex(1) - 1);
            cv(var,ct) = v;
            cv(ct,var) = v;
            % do the lower part of c also.
        end
    end
    c=cv;
    
% This a simple wrapper function to get just the eigenvectors
% since the system call returns 3 matrices
function [x]=justeigs (x)
% This just returns the eigenvectors of the matrix
    [a, eig, b] = bdiag(x);
    x= eig;
%  this function makes the transformation to the eigenspace for PCA
%  parameters:
%  adjusteddata = mean-adjusted data set
%  eigenvectors = SORTED eigenvectors (by eigenvalue)
%  dimensions = how many eigenvectors you wish to keep
% 
%  The first two parameters can come from the result of calling
%  PCAprepare on your data.
%  The last is up to you.
function [finaldata] = PCAtransform(adjusteddata,eigenvectors,dimensions)
    finaleigs = eigenvectors(:,1:dimensions);
    prefinaldata = finaleigs'*adjusteddata';
    finaldata = prefinaldata';
%  This function does the preparation for PCA analysis
%  It adjusts the data to subtract the mean, finds the covariance matrix,
%  and finds normal eigenvectors of that covariance matrix.
%  It returns 4 matrices
%  meanadjust = the mean-adjust data set
%  covmat = the covariance matrix of the data
%  eigvalues = the eigenvalues of the covariance matrix, IN SORTED ORDER
%  normaleigs = the normalised eigenvectors of the covariance matrix,
%  IN SORTED ORDER WITH RESPECT TO
%  THEIR EIGENVALUES, for selection for the feature vector.
% 
%  NOTE: This function cannot handle data sets that have any eigenvalues
%  equal to zero. It’s got something to do with the way that scilab treats
%  the empty matrix and zeros.

function [meanadjusted,covmat,sorteigvalues,sortnormaleigs] = PCAprepare (data)
% Calculates the mean adjusted matrix, only for 2 dimensional data
means = mean(data,'r');
meanadjusted = meanadjust(data);
covmat = cov(meanadjusted);
eigvalues = spec(covmat);
normaleigs = justeigs(covmat);
sorteigvalues = sorteigvectors(eigvalues', eigvalues');
sortnormaleigs = sorteigvectors(eigvalues', normaleigs);


%  This removes a specified column from a matrix
%  A = the matrix
%  n = the column number you wish to remove
function [columnremoved] = removecolumn(A,n)
inputsize = size(A);
numcols = inputsize(2);
temp = A(:,1:(n-1));
for var = 1:(numcols - n)
temp(:,(n+var)-1) = A(:,(n+var));
end,
columnremoved = temp;


%  This finds the column number that has the
%  highest value in it’s first row.
function [column] = highestvalcolumn(A)
inputsize = size(A);
numcols = inputsize(2);
maxval = A(1,1);
maxcol = 1;
for var = 2:numcols
if A(1,var) > maxval
maxval = A(1,var);
maxcol = var;
end,
end,
column = maxcol

%  This sorts a matrix of vectors, based on the values of
%  another matrix
% 
%  values = the list of eigenvalues (1 per column)
%  vectors = The list of eigenvectors (1 per column)
% 
%  NOTE: The values should correspond to the vectors
%  so that the value in column x corresponds to the vector
%  in column x.
function [sortedvecs] = sorteigvectors(values,vectors)
    inputsize = size(values);
    numcols = inputsize(2);
    highcol = highestvalcolumn(values);
    sorted = vectors(:,highcol);
    remainvec = removecolumn(vectors,highcol);
    remainval = removecolumn(values,highcol);
    for var = 2:numcols
        highcol = highestvalcolumn(remainval);
        sorted(:,var) = remainvec(:,highcol);
        remainvec = removecolumn(remainvec,highcol);
        remainval = removecolumn(remainval,highcol);
    end
    sortedvecs = sorted;


%  This takes a set of data, and subtracts
%  the column mean from each column.
function [meanadjusted] = meanadjust(Data)
    inputsize = size(Data);
    numcols = inputsize(2);
    means = mean(Data,'r');
    tmpmeanadjusted = Data(:,1) - means(:,1);
    for var = 2:numcols
        tmpmeanadjusted(:,var) = Data(:,var) - means(:,var);
    end
    meanadjusted = tmpmeanadjusted

附錄B：本文相關的Matlab做圖代碼

1. 【協方差部分】H和M的關係圖

%% ctrl + R 註釋
%% ctrl + T 取消註釋

close all;
clc;
clear;
fontsize = 11;
figure;
% x = 0:0.00025:1;
% y = x.*sin(10*3.14.*x)+2;
% plot(x,y)

x = [9 15 25 14 10 18 0 16 5 19 16 20];
y = [39 56 93 61 50 75 32 85 42 70 66 80];
plot(x, y, 'o');
xlabel({'H / 小時'}, 'FontSize',fontsize);
ylabel({'M / 分數'}, 'FontSize',fontsize);

h = legend(['H(學習時間)和M(得到分數)的關係', sprintf('\n')], 'location', 'best');
grid on;

set(h, 'Box', 'off');

參考文獻

1. PCA數學原理：http://www.360doc.com/content/13/1124/02/9482_331688889.shtml
2. 杜洪，夏欣，琚生根，王能。基於PCA圖像壓縮算法研究與實現：http://www.cnki.com.cn/Article/CJFDTotal-SCDX201405009.htm
3. 本文對應的原文：【PDF】http://pan.baidu.com/s/1nvNPJ4H   |   【文章出處】https://arxiv.org/abs/1404.1100
4. 博主的CSDN主頁：http://blog.csdn.net/houchaoqun_xmu
5. 【Latex】製做演示文檔或者課程報告-常見用法（一）
6. 【Latex】製做演示文檔或者課程報告 - 製做演示文檔（二）
7. 【Latex】製做演示文檔或者課程報告 - 製做課程報告（三）
8. TSP_旅行商問題 - 蠻力法DFS（一）
9. TSP_旅行商問題 - 模擬退火算法（三）
10. TSP_旅行商問題 - 遺傳算法（四）