基礎(chǔ)知識(shí)：https://blog.csdn.net/m0_37957160/article/details/107063223

https://mp.csdn.net/editor/html/107108718

PS：

先看一下PCA算法實(shí)現(xiàn)流程：（設(shè)有m條n維數(shù)據(jù)，意思就是有m個(gè)樣本，每個(gè)樣本有n維特征）

（1）將原始數(shù)據(jù)組成一個(gè)m行n列的矩陣X。

（2）將X的每一列（代表一個(gè)屬性字段），進(jìn)行零均值化處理，即減去這一列的均值。

（3）求出協(xié)方差矩陣。

（4）求出協(xié)方差矩陣的特征值及對(duì)應(yīng)的特征向量。

（5）將特征向量按對(duì)應(yīng)特征值大小，從上到下按行排列成矩陣，取前k行組成矩陣P。

（6）Y=PX即為降維到k維后的數(shù)據(jù)。（意思就是將原始數(shù)據(jù)旋轉(zhuǎn)到特征矩陣P所在的空間中，得到的數(shù)據(jù)Y即為降維后的結(jié)果）

----------------協(xié)方差矩陣--------

協(xié)方差矩陣通常這樣計(jì)算，先讓樣本矩陣中心化，即每一維度減去該維度的均值，使每一維度上的均值為0，然后直接用新得到的樣本矩陣乘上它的轉(zhuǎn)置，然后除以(N-1)即可。其實(shí)這種方法也是由前面的公式通道而來，只不過理解起來不是很直觀，但在抽象的公式推導(dǎo)時(shí)還是很常用的！

協(xié)方差總結(jié)：

1、協(xié)方差矩陣是對(duì)稱的

2、對(duì)角線是各個(gè)維度的方差

3、協(xié)方差矩陣表示各個(gè)維度之間的相關(guān)性

4、協(xié)方差矩陣中元素的正負(fù)值表示維度之間的正相關(guān)、負(fù)相關(guān)關(guān)系

---------------基于協(xié)方差矩陣的特征值特征向量求解------------------

#include <vtkAutoInit.h>       
VTK_MODULE_INIT(vtkRenderingOpenGL);
VTK_MODULE_INIT(vtkInteractionStyle);
VTK_MODULE_INIT(vtkRenderingFreeType);
#include <string>
#include <iostream>
#include <pcl/io/io.h>
#include <pcl/io/pcd_io.h>
#include <pcl/point_cloud.h>
#include <pcl/point_types.h>		  
#include <vector>
#include <pcl/visualization/pcl_visualizer.h>
#include <Eigen/Core>
#include <pcl/common/transforms.h>using namespace std;
typedef pcl::PointXYZ PointType;
typedef pcl::Normal NormalType;int main(int argc, char **argv)
{pcl::PointCloud<PointType>::Ptr cloud(new pcl::PointCloud<PointType>());pcl::PointCloud<NormalType>::Ptr cloud_normal(new pcl::PointCloud<NormalType>());//導(dǎo)入點(diǎn)云std::string fileName(argv[1]);pcl::io::loadPCDFile(fileName, *cloud);//加載點(diǎn)云// PCA：計(jì)算主方向Eigen::Vector4f pcaCentroid;//質(zhì)心pcl::compute3DCentroid(*cloud, pcaCentroid);//估計(jì)質(zhì)心的坐標(biāo) // 齊次坐標(biāo)，（c0,c1,c2,1）cout << "質(zhì)心坐標(biāo)："<<pcaCentroid << endl;Eigen::Matrix3f covariance;pcl::computeCovarianceMatrixNormalized(*cloud, pcaCentroid, covariance);//計(jì)算歸一化協(xié)方差矩陣cout << "協(xié)方差矩陣：" << covariance << endl;// 計(jì)算主方向：特征向量和特征值Eigen::SelfAdjointEigenSolver<Eigen::Matrix3f> eigen_solver(covariance, Eigen::ComputeEigenvectors);Eigen::Matrix3f eigenVectorsPCA = eigen_solver.eigenvectors();//特征向量std::cout << "特征向量：" << eigenVectorsPCA << std::endl;//輸出特征向量Eigen::Vector3f eigenValuesPCA = eigen_solver.eigenvalues();//特征值std::cout << "特征值：" << eigenValuesPCA << std::endl;//輸出特征值//eigenVectorsPCA.col(2) = eigenVectorsPCA.col(0).cross(eigenVectorsPCA.col(1));//cross叉乘//校正主方向間垂直（特征向量方向： (e0, e1, e0 × e1) --- note: e0 × e1 = +/- e2）// 轉(zhuǎn)到參考坐標(biāo)系，將點(diǎn)云主方向與參考坐標(biāo)系的坐標(biāo)軸進(jìn)行對(duì)齊Eigen::Matrix4f transform(Eigen::Matrix4f::Identity());transform.block<3, 3>(0, 0) = eigenVectorsPCA.transpose();//轉(zhuǎn)置// R^(-1) = R^Ttransform.block<3, 1>(0, 3) = -1.0f * (transform.block<3, 3>(0, 0)) * (pcaCentroid.head<3>());// t^(-1) = -R^T * t//獲取向量的前n個(gè)元素：vector.head(n);//利用block()函數(shù)，可以從Matrix中取出一個(gè)小矩陣來進(jìn)行處理，使用的語法為matrix.block(i, j, p, q); //和 matrix.block<p, q>(i, j)。其中i，j是block左上角元素位于矩陣中的位置，p、q是block的大小，如果是小塊的話，pcl::PointCloud<PointType>::Ptr transformedCloud(new pcl::PointCloud<PointType>);//變換后的點(diǎn)云pcl::transformPointCloud(*cloud, *transformedCloud, transform);//轉(zhuǎn)換到原點(diǎn)時(shí)的主方向PointType o;o.x = 0.0;o.y = 0.0;o.z = 0.0;Eigen::Affine3f tra_aff(transform);Eigen::Vector3f pz = eigenVectorsPCA.col(0);//列Eigen::Vector3f py = eigenVectorsPCA.col(1);Eigen::Vector3f px = eigenVectorsPCA.col(2);pcl::transformVector(pz, pz, tra_aff);pcl::transformVector(py, py, tra_aff);pcl::transformVector(px, px, tra_aff);PointType pcaZ;pcaZ.x = 1000 * pz(0);pcaZ.y = 1000 * pz(1);pcaZ.z = 1000 * pz(2);PointType pcaY;pcaY.x = 1000 * py(0);pcaY.y = 1000 * py(1);pcaY.z = 1000 * py(2);PointType pcaX;pcaX.x = 1000 * px(0);pcaX.y = 1000 * px(1);pcaX.z = 1000 * px(2);//初始位置時(shí)的主方向PointType c;c.x = pcaCentroid(0);c.y = pcaCentroid(1);c.z = pcaCentroid(2);PointType pcZ;pcZ.x = 1000 * eigenVectorsPCA(0, 0) + c.x;pcZ.y = 1000 * eigenVectorsPCA(1, 0) + c.y;pcZ.z = 1000 * eigenVectorsPCA(2, 0) + c.z;PointType pcY;pcY.x = 1000 * eigenVectorsPCA(0, 1) + c.x;pcY.y = 1000 * eigenVectorsPCA(1, 1) + c.y;pcY.z = 1000 * eigenVectorsPCA(2, 1) + c.z;PointType pcX;pcX.x = 1000 * eigenVectorsPCA(0, 2) + c.x;pcX.y = 1000 * eigenVectorsPCA(1, 2) + c.y;pcX.z = 1000 * eigenVectorsPCA(2, 2) + c.z;//visualizationpcl::visualization::PCLVisualizer viewer;pcl::visualization::PointCloudColorHandlerCustom<PointType> color_handler(cloud, 255, 0, 0);//對(duì)輸入原始點(diǎn)云處理//這句話的意思是：對(duì)輸入為pcl::PointXYZ類型的點(diǎn)云，著色為紅色。其中，cloud表示真正處理的點(diǎn)云，color_handler表示處理結(jié)果.pcl::visualization::PointCloudColorHandlerCustom<PointType> tc_handler(transformedCloud, 0, 255, 0);viewer.addPointCloud(cloud, color_handler, "cloud");viewer.addPointCloud(transformedCloud, tc_handler, "transformCloud");//將點(diǎn)云transformedCloud,處理結(jié)果為tc_handler,添加到視圖中，其中,雙引號(hào)中的transformCloud,表示該點(diǎn)云的”標(biāo)簽“，//我們依然可以稱之為”名字“，之所以設(shè)置各個(gè)處理點(diǎn)云的名字，是為了在后續(xù)處理中易于區(qū)分; //v1表是添加到哪個(gè)視圖窗口（pcl中可設(shè)置多窗口模式）viewer.addArrow(pcaZ, o, 0.0, 0.0, 1.0, false, "arrow_Z");//給變化坐標(biāo)的綠色點(diǎn)云tc_handler添加箭頭viewer.addArrow(pcaY, o, 0.0, 1.0, 0.0, false, "arrow_Y");viewer.addArrow(pcaX, o, 1.0, 0.0, 0.0, false, "arrow_X");viewer.addArrow(pcZ, c, 0.0, 0.0, 1.0, false, "arrow_z");viewer.addArrow(pcY, c, 0.0, 1.0, 0.0, false, "arrow_y");viewer.addArrow(pcX, c, 1.0, 0.0, 0.0, false, "arrow_x");//viewer.setCameraPosition(100.0, 20.0, 35.0, 4.0, 5.0, 6.0);//設(shè)置坐標(biāo)原點(diǎn)viewer.addCoordinateSystem(100);//建立空間直角坐標(biāo)系//viewer.initCameraParameters();	//初始化相機(jī)參數(shù)，讓用戶在默認(rèn)的角度下觀察點(diǎn)云viewer.setBackgroundColor(1.0, 1.0, 1.0);while (!viewer.wasStopped()){viewer.spinOnce(10000);}return 0;
}

顯示結(jié)果如下：

?

總結(jié)

以上是生活随笔為你收集整理的PCA（3）：PCA实现C++代码的全部內(nèi)容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網(wǎng)站內(nèi)容還不錯(cuò)，歡迎將生活随笔推薦給好友。