在科研、工程應用、生活中，我們所獲取的數據往往包含着很多冗餘信息，這些冗餘信息往往對數據分析造成幹擾，增加數據分析的複雜度。此時我們則需要對這些數據進行預處理，預處理的原則是：既能抓住其主要特征，又能剔除冗餘信息，從而減少數據量。PCA降維就是這樣的一種數據預處理算法。

本文首先講解PCA降維的計算原理，再使用C 與Opencv來實現該算法，并與Opencv現有的PCA函數接口進行降維結果的對比。看到這裡，可能有人會問，Opencv都有現成的函數可以調用了，為什麼還要自己去寫呢？我想說的是，對于學習者來說，重複造輪子并不是壞事，它可以讓我們更加深刻地理解造輪子的過程，從而才有改進和創新的機會，相反，如果隻是使用别人造好的輪子，不深刻理解其構造原理，壓根就沒有改進創新的空間了。

1. 計算原理

假設有m行n列的數據，計為矩陣X0，其每一行數據看作一個一維行向量，那麼該數據本來有m個一維行向量，我們要使用PCA降維算法把其降為k個一維行向量(k < m)，計算過程如下：

(1) 求出每行數據的平均值。

(2) 去平均處理，把每行數據都減去本行數據的平均值。

(3) 去平均處理之後，同樣得到m行n列的數據，計為矩陣X1。

(4) 按以下公式計算X0的協方差矩陣，其中“*”表示矩陣乘法，得到的協方差矩陣Cov為m行m列矩陣：

(5) 計算協方差矩陣Cov的特征值與對應的特征向量。得到m個特征值，對應m個特征向量，其中每個特征向量的長度又為m，也即所有特征向量組成m行m列的矩陣。

(6) 将特征值按照從大到小排列，并根據排列後特征值的順序來按行從上到下排列特征向量(每個特征向量為一行)，使特征值與特征向量仍保持對應。

比如本來特征值依次為a1，a2，a3，a4，a5，對應的特征向量為：

v1

v2

v3

v4

v5

排序之後，a3>a1>a4>a5>a2，那麼特征向量的順序也作對應的調整：

v3

v1

v4

v5

v2

此時計特征向量組成的矩陣為V。

(7) 取矩陣V的前k行數據，得到k行m列的矩陣P，計算Y=P*X1，Y矩陣即為最終得到的k行n列的降維數據，從而實現把數據從m維降為n維。

上代碼：

voiddo_PCA(Matsrc,Mat&pca,intk) { if (src.rows <= 1 || src.cols <= 1) { return; } k = k > src.rows ? src.rows : k; Mat src_float; src.convertTo(src_float,CV_32F); //求每行的平均值 Mat col_mean = Mat::zeros(1, src.rows, CV_32FC1); float *col_mean_p = col_mean.ptr<float>(0); for (int i = 0; i < src_float.rows; i ) { float *p = src_float.ptr<float>(i); for (int j = 0; j < src_float.cols; j ) { col_mean_p[i] = p[j]; } } col_mean /= src_float.cols; //去平均值 for (int i = 0; i < src_float.rows; i ) { float *p = src_float.ptr<float>(i); for (int j = 0; j < src_float.cols; j ) { p[j] = p[j] - col_mean_p[i]; } } //計算協方差矩陣 Mat X = src_float*src_float.t(); X = X / src_float.cols; Mat eValuesMat, eVectorsMat; //調用opencv接口計算特征值與特征向量 eigen(X,eValuesMat,eVectorsMat);//這裡得到的特征值已經按照從大到小排序了，特征向量也與特征值相對應 Mat Vectors_k; eVectorsMat(Rect(0,0,eVectorsMat.cols,k)).copyTo(Vectors_k);//取特征向量的前k行,k*r pca=Vectors_k*src_float;//k*r*r*c=k*c pca=pca.clone();//确保矩陣連續，拷貝一份 printf("pca.rows=%d,pca.cols=%d\n",pca.rows,pca.cols); }

測試代碼：

void pca_test(void) { Matimg=imread("lena.jpg",CV_LOAD_IMAGE_GRAYSCALE); intk=100;//從圖像原有的m行降為100行，列不變 Matpca_m; do_PCA(img,pca_m,k); imshow("原圖", img); imshow("本文實現算法 PCA降維後", pca_m); //使用Opencv的PCA函數接口 PCApca(img,Mat(),CV_PCA_DATA_AS_COL,k); Matdst=pca.project(img);//dst則為最終降維數據 imshow("Opencv PCA降維後",dst); waitKey(0); }

運行上述代碼，得到的結果如下。可以看到，本文實現的算法與Opencv函數的計算結果是一緻的，說明我們的計算過程沒錯，鼓掌~再接再厲~

原圖

本文實現的PCA算法的降維結果

Opencv現有的PCA算法的降維結果