PCA算法与Opencv实现

转载： http://blog.sina.com.cn/s/blog_4b9b714a0100hc8o.html

对于PCA，一直都是有个概念，没有实际使用过，今天终于实际使用了一把，发现PCA还是挺神奇的。

在OPENCV中使用PCA非常简单，只要几条语句就可以了。

1、初始化数据

//每一行表示一个样本

CvMat* pData = cvCreateMat( 总的样本数, 每个样本的维数, CV_32FC1 );

CvMat* pMean = cvCreateMat(1, 样本的维数, CV_32FC1);

//pEigVals中的每个数表示一个特征值

CvMat* pEigVals = cvCreateMat(1, min(总的样本数,样本的维数), CV_32FC1);

//每一行表示一个特征向量

CvMat* pEigVecs = cvCreateMat( min(总的样本数,样本的维数), 样本的维数, CV_32FC1);

2、PCA处理,计算出平均向量pMean,特征值pEigVals和特征向量pEigVecs

cvCalcPCA( pData, pMean, pEigVals, pEigVecs, CV_PCA_DATA_AS_ROW );

3、选出前P个特征向量(主成份),然后投影,结果保存在pResult中，pResult中包含了P个系数

CvMat* pResult = cvCreateMat( 总的样本数, PCA变换后的样本维数(即主成份的数目), CV_32FC1 );

cvProjectPCA( pData, pMean, pEigVecs, pResult );

4、 重构,结果保存在pRecon中

CvMat* pRecon = cvCreateMat( 总的样本数, 每个样本的维数, CV_32FC1 );

cvBackProjectPCA( pResult, pMean, pEigVecs, pRecon );

5、重构误差的计算

计算pRecon和pData的"差"就可以了.

使用时如果是想用PCA判断“是非”问题，则可以先用正样本计算主成分，判断时，对需要判断得数据进行投影，然后重构，计算重构出的数据与原数据的差异，如果差异在给定范围内，可以认为“是”。

如果相用PCA进行分类，例如对数字进行分类，则先用所有数据（0-9的所有样本）计算主成分，然后对每一类数据进行投影，计算投影的系数，可简单得求平均。即对每一类求出平均系数。分类时，将需要分类得数据进行投影，得到系数，与先前计算出得每一类得平均系数进行比较，可判为最接近得一类。当然这只是最简单得使用方法。

PCA（principal component analysis，主成分分析），又称为k-l变换，我想是大家用的最多的降维手段，对于PCA的理解，我想大神们都各有各的绝招，可以应用的场合也非常多。下面就介绍一下OpenCV中PCA这个类，因为常用，所以这个类相对OpenCV而言显得比较独立，放在了core这部分中。

PCA类的成员函数包括构造函数、运算符重载（）、project、backProject这几个函数，还包括成员变量eigenvectors、eigenvalues、mean。使用也很方便。比如我要计算一组向量的PCA，我们只需要定义个PCA实例，获得主成分，调用project测试新样本，也可以再调用backProject重建原始向量，是project的一个逆运算。

在实际使用我我发现了一个bug，构造函数和运算符重载（）函数的第三个参数只能使用0或者1（分别对应CV_PCA_DATA_AS_ROW和CV_PCA_DATA_AS_COL），但调用对应的宏却编译无法通过，不知道大家是否也是一样的情况

下面结合一个实例分析：

2001年全国各地区消费情况指数

地区 食品X1 衣着X2 家庭设备X3 医疗保健X4 交通通信X5 娱乐教育X6 居住X7

北京 101.5 100.4 97.0 98.7 100.8 114.2 104.2
天津 100.8 93.5 95.9 100.7 106.7 104.3 106.4
河北 100.8 97.4 98.2 98.2 99.5 103.6 102.4
山西 99.4 96.0 98.2 97.8 99.1 98.3 104.3
内蒙古 101.8 97.7 99.0 98.1 98.4 102.0 103.7
辽宁 101.8 96.8 96.4 92.7 99.6 101.3 103.4
吉林 101.3 98.2 99.4 103.7 98.7 101.4 105.3
黑龙江 101.9 100.0 98.4 96.9 102.7 100.3 102.3
上海 100.3 98.9 97.2 97.4 98.1 102.1 102.3
江苏 99.3 97.7 97.6 101.1 96.8 110.1 100.4
浙江 98.7 98.4 97.0 99.6 95.6 107.2 99.8
安徽 99.7 97.7 98.0 99.3 97.3 104.1 102.7
福建 97.6 96.5 97.6 102.5 97.2 100.6 99.9
江西 98.0 98.4 97.1 100.5 101.4 103.0 99.9
山东 101.1 98.6 98.7 102.4 96.9 108.2 101.7
河南 100.4 98.6 98.0 100.7 99.4 102.4 103.3
湖北 99.3 96.9 94.0 98.1 99.7 109.7 99.2
湖南 98.6 97.4 96.4 99.8 97.4 102.1 100.0
广东 98.2 98.2 99.4 99.3 99.7 101.5 99.9
广西 98.5 96.3 97.0 97.7 98.7 112.6 100.4
海南 98.4 99.2 98.1 100.2 98.0 98.2 97.8
重庆 99.2 97.4 95.7 98.9 102.4 114.8 102.6
四川 101.3 97.9 99.2 98.8 105.4 111.9 99.9
贵州 98.5 97.8 94.6 102.4 107.0 115.0 99.5
云南 98.3 96.3 98.5 106.2 92.5 98.6 101.6
西藏 99.3 101.1 99.4 100.1 103.6 98.7 101.3
陕西 99.2 97.3 96.2 99.7 98.2 112.6 100.5
甘肃 100.0 99.9 98.2 98.3 103.6 123.2 102.8
青海 102.2 99.4 96.2 98.6 102.4 115.3 101.2
宁夏 100.1 98.7 97.4 99.8 100.6 112.4 102.5
新疆 104.3 98.7 100.2 116.1 105.2 101.6 102.6 

我们用PCA来分析这些向量的对消费水平的贡献，代码实现如下：

#include <opencv2/core/core.hpp>
#include <iostream>
#include <fstream>

using namespace cv;
using namespace std;

void main(int argc, char** argv)
{
 const int SAMPLE_NUM = 31;
 const int VECTOR_NUM = 7;
 fstream fp;
 string temp;
 char* filename = argv[1];
 Mat pcaSet(SAMPLE_NUM, VECTOR_NUM, CV_32F);
 // 读取数据，生成pcaSet
 fp.open(filename, ios::in);
 for (int i=0; i<SAMPLE_NUM; i++)
 {
  getline(fp, temp, ' ');
  for (int j=0; j<VECTOR_NUM; j++)
  {
   fp >> pcaSet.at<float>(i, j);
  }
 }
 fp.close();

//  cout << pcaSet << endl;
 PCA pca(pcaSet, Mat(), 0/*CV_PCA_DATA_AS_ROW*/);
 cout << pca.eigenvalues << endl;
 cout << endl;
 cout << pca.eigenvectors << endl;
}

我这里只把前3维的结果显示如下：

eigenvalues：[43.18，14.60，9.21]、

eigenvectors：[0.013，0.034，-.0.100， -0.130， 0.207，0.963，-0.020；

    0.157，0.038，0.121，0.893，0.395，0.046，0.060；

    0.214，0.018，-0.001，-0.404，0.813，-0.228，0.277]