PCA（主成分分析）原理、步骤及代码实现

1、PCA处理步骤：

假设一组数据：m条n维数据

例如点云：
x1, y1, z1            x1, x2, ...xn
x2, y2, z2      ->    y1, y2, ...yn
...                   z1, z2, ...zn
xn, yn zn

1. 将原始数据组成n行m列矩阵X（每一行为同一字段）
   2.将X的每一行进行零均值化（即每一行数据均减去改行均值）
   3.求出均值化后矩阵的协方差矩阵
   4.求出协方差矩阵对应的特征值、特征向量
   5.将特征向量依据特征值的大小从上到下，按行排列成矩阵，取前K行组成矩阵P
   6.Y = PX，Y即为降维到K维后的数据。

2、数学原理可以参考此处

3、C++借助Eigen库代码实现

//基于Eigen库实现PCA，此代码经本人调试完成，未出现错误
//示例数据：5个2维数据   
/* 
1 1
1 3
2 3
4 4
2 4
*/

#include <iostream>
#include <algorithm>
#include<cstdlib>
#include <fstream>
#include <Eigen/Dense>

using namespace std;
using namespace Eigen;

MatrixXd featurnormal(MatrixXd &X)
{
	// Eigen:基本功能
	// VectorXd    等价于Matrix<double,Dynamic,1>;       //列向量，一列
	// RowVectorXd 等价于Matrix<double,1,Dynamic>;       //行向量，一行
	// MatrixXd    等价于Matrix<double,Dynamic,Dynamic>;

	//计算每一维的均值
	MatrixXd X1 = X.transpose();  //先来个普通转置【行列变换】

	VectorXd meanval = X1.rowwise().mean();//colwise按列求均值。//rowwise按行求均值。 //X1.rowwise()不可用
	cout << "meanval" << endl << meanval << endl;
	
	//样本均值化为0	
	X1.colwise() -= meanval;   //按列减，每一列与均值做差
	
	return X1;
}

void ComputeCov(MatrixXd &X1, MatrixXd &C)
{
	//计算协方差矩阵设我们有m个n维数据记录，将其按列排成n乘m的矩阵X，设C =X*XT/m
    //此处除以m，具体原因你深思[可参考此文](https://blog.csdn.net/hustqb/article/details/78394058)
		
	C = X1*X1.adjoint();//相当于原始数据XT*X    adjiont()求伴随矩阵 相当于求矩阵的转置
	/*cout << "C:" << endl << C<< endl;*/
	C = C.array() / X1.cols();//C.array()矩阵的数组形式  ///逐元素除以m

	cout << "C:" << endl << C << endl;
}

void ComputEig(MatrixXd &C, MatrixXd &vec, MatrixXd &val)
{
	//计算特征向量和特征值 SelfAdjointEigenSolver自动将计算得到的特征向量和特征值排序
	SelfAdjointEigenSolver<MatrixXd> eig(C);
	//所求特征向量依据对应的特征值从小到大排列。
	vec = eig.eigenvectors();
	val = eig.eigenvalues();
	
	//故应该依靠翻转将特征向量依据特征值从大到小排列。【如下】
	vec.colwise().reverse();
	cout << "***************************************" << endl;
	cout << "vec.colwise().reverse():" << endl << vec.colwise().reverse() << endl;
	cout << "***************************************" << endl;
	val.colwise().reverse();
	cout << "val.colwise().reverse():" << endl << val.colwise().reverse() << endl;
	cout << "***************************************" << endl;
}

// 计算维度
int ComputDim(MatrixXd &val)
{
	int dim;
	double sum = 0;
	for (int i = val.rows() - 1; i >= 0; --i)
	{
		sum += val(i, 0);
		dim = i;
		if (sum / val.sum() >= 0.8)//达到所需要的要求时
			break;
	}
	//cout << val.rows() - dim << endl;
	return val.rows() - dim;
}

int main()
{
	//文件操作   文件路径使用： \\ 或者 /
	ifstream fin("F:\\PCLprincipal\\pca070808\\Project1\\m5n2.txt");
	if (fin.is_open())
	{
		cout << "载入成功" << endl;
	}
	ofstream fout("X1X2outpot1.txt");

	//定义所需要的量
	const int m = 5, n = 2, z = 2;
	MatrixXd X(m, n), C(z, z);  //C为协方差阵
	MatrixXd vec, val;

	//读取数据
	double in[z];
	for (int i = 0; i < m; ++i)
	{
		for (int j = 0; j < n; ++j)
		{
			fin >> in[j];  //单值读入
		}
		for (int j = 1; j <= n; ++j)
			X(i, j - 1) = in[j - 1];
	}

	cout << endl;
	cout << "X为原始的数据：" << endl << X << endl;
	cout << endl;
	
	// 1 均值0标准化
	MatrixXd X1 = featurnormal(X);
	cout << "X1均值0标准化的数据：" << endl << X1.transpose() << endl;   ///.transpose()方便为了查看 N*3
	cout << endl;

	// 2 计算协方差
	ComputeCov(X1, C);

	// 3 计算特征值和特征向量
	ComputEig(C, vec, val);

	// 4 计算损失率，确定降低的维数
	int dim = ComputDim(val);
	//int dim = 1;
	// 5计算结果

	//但当涉猎
	//MatrixXd res001 = vec.transpose().rightCols(dim).transpose();
	//MatrixXd res001 = vec.bottomRows(dim);
	//MatrixXd res001 = vec.transpose().leftCols(dim).transpose();
	/*cout <<"X.colwise().reverse()"<<endl<< X.colwise().reverse() << endl;*/ //垂直翻转，优秀的很，

	MatrixXd res001 = vec.colwise().reverse().topRows(dim);  
	//按列翻转取前dim行 特征向量，实现依据特征值从大到小排列获取特征向量。
	
	cout << "res001:" << endl << res001 << endl;  //2*3
	MatrixXd res = res001*X1; //2*3  *  N*3  Y=PX 此Y即为降维之后的结果。

	//方便查看
	cout << "res为用PCA算法之后的数据：" << endl << res.transpose() << endl;
	
	// 6输出结果
	fout << "the result is" << res.rows() << "x" << res.cols() << "after pca algorithm" << endl;
	fout << res.transpose();
	fout.close();
	
	system("pause");
	return 0;
}

• 矩阵：依据矩阵有实数矩阵与复数矩阵数据类型分为：
• 转置：共轭转置、普通转置
  1.普通转置：实数矩阵的共轭转置矩阵就是转置矩阵【行与列对换】
  2.共轭转置：复数矩阵的共轭转置矩阵就是将形如a+bi的数变成a-bi,实数的共轭是它本身.
• Eigen库详解1
• Eigen库详解2
• Eigen库详解3
• Eigen库详解4