在科研、工程应用、生活中，我们所获取的数据往往包含着很多冗余信息，这些冗余信息往往对数据分析造成干扰，增加数据分析的复杂度。此时我们则需要对这些数据进行预处理，预处理的原则是：既能抓住其主要特征，又能剔除冗余信息，从而减少数据量。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维。
上代码：
-
void do_PCA(Mat src, Mat &pca, int k)
-
{
-
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)
-
{
-
Mat img = imread(
"lena.jpg", CV_LOAD_IMAGE_GRAYSCALE);
-
-
-
int k =
100;
//从图像原有的m行降为100行，列不变
-
Mat pca_m;
-
do_PCA(img, pca_m, k);
-
-
-
imshow(
"原图", img);
-
imshow(
"本文实现算法 PCA降维后", pca_m);
-
-
//使用Opencv的PCA函数接口
-
PCA pca(img, Mat(), CV_PCA_DATA_AS_COL, k);
-
Mat dst = pca.project(img);
//dst则为最终降维数据
-
-
imshow(
"Opencv PCA降维后", dst);
-
waitKey(
0);
-
}
运行上述代码，得到的结果如下。可以看到，本文实现的算法与Opencv函数的计算结果使一致的，说明我们的计算过程没错，鼓掌~再接再厉~
原图
本文实现的PCA算法的降维结果
Opencv现有的PCA算法的降维结果
微信公众号如下，欢迎关注，欢迎交流：
转载：https://blog.csdn.net/shandianfengfan/article/details/111351399