C++实战:从零实现PCA与LDA降维算法,掌握核心原理与工程细节
1. 项目概述从理论到代码亲手实现经典降维算法最近在整理一些机器学习的老项目发现很多朋友对PCA主成分分析和LDA线性判别分析这两个经典的降维算法“既熟悉又陌生”。熟悉是因为几乎每本教材都会提到陌生是当真正需要自己动手用C从零实现或者在一个实际项目中集成它们时又会遇到一堆细节问题。比如特征值分解到底怎么算类内散度矩阵和类间散度矩阵的公式背下来了但代码里维度对不上怎么办降维后的数据怎么可视化才有效这个项目实战就是来解决这些问题的。我将带你用纯C不依赖复杂的机器学习库一步步构建出PCA和LDA的核心算法并完成一个从数据加载、预处理、模型训练到结果可视化的完整流程。无论你是正在学习《模式识别》或《机器学习》课程的学生还是希望夯实算法底层实现能力的C开发者这个实战都能让你对降维的理解从“知道”跃升到“做到”。2. 核心算法原理与设计思路拆解2.1 为什么需要数据降维PCA与LDA的角色定位在我们开始敲代码之前必须搞清楚一个根本问题为什么要降维想象一下你面前有一张记录了100个特征的表格比如100种身体指标直接把它扔给一个分类器效果往往不好。这就是“维度灾难”数据在高维空间中变得极其稀疏容易过拟合计算成本也剧增。降维的核心目标就是在尽可能保留原始数据重要信息的前提下将数据投影到一个低维空间。PCA和LDA是解决这个问题的两把经典钥匙但开锁的思路截然不同。PCA是一种无监督的降维方法。它不考虑数据的标签类别只关注数据本身的方差。PCA认为方差大的方向承载了更多的信息。它的任务就是找到一组新的正交坐标轴主成分使得数据在这些轴上的投影方差最大。所以PCA常用于数据压缩、去噪和可视化比如把人脸图像的高维像素点压缩到几十个“特征脸”上。LDA则是一种有监督的降维方法。它充分利用了数据的类别标签。LDA的目标是找到一个投影方向使得类间散度不同类别数据中心之间的距离尽可能大同时类内散度同一类别数据点之间的离散程度尽可能小。简单说就是让投影后的数据“类内紧致类间分离”。因此LDA天生就是为了分类任务服务的在投影的同时就为后续的分类器创造了有利条件。选择PCA还是LDA取决于你的目标。如果你的数据没有标签或者你只想探索数据结构、减少存储开销选PCA。如果你的核心任务是分类并且拥有清晰的类别标签那么LDA通常是更好的起点。2.2 数学基石特征值分解与广义特征值分解两种算法的实现都绕不开线性代数中的核心运算。理解它们是写出正确代码的前提。对于PCA其核心是计算数据协方差矩阵的特征值和特征向量。给定一个中心化减去均值后的数据矩阵X(n_samples x n_features)其协方差矩阵C (X^T * X) / (n_samples-1)。我们需要求解特征方程C * w λ * w。这里特征值λ的大小代表了对应特征向量w方向上数据的方差特征向量w就是主成分的方向。我们按特征值从大到小排序选取前k个特征向量组成投影矩阵W降维操作就是X_reduced X * W。在C实现中我们将采用雅可比(Jacobi)迭代法来求解实对称矩阵的全部特征值和特征向量因为它稳定且易于实现虽然对于超大矩阵不是最快的但对于学习和小规模数据足够了。对于LDA问题转化为求解广义特征值问题S_b * w λ * S_w * w。其中S_b是类间散度矩阵S_w是类内散度矩阵。我们需要找到使广义瑞利商(w^T * S_b * w) / (w^T * S_w * w)最大化的向量w。这通常通过求解S_w^{-1} * S_b的特征向量来实现。但这里有个关键陷阱S_w经常是奇异的特别是当样本维度大于样本数时。在C实现中我们将采用正则化和PCA预降维两步走的策略来稳定求解。先对数据用PCA降到 (n_classes - 1) 维以内这是LDA能降到的最大维度确保S_w可逆再求解广义特征值问题。注意直接计算S_w^{-1} * S_b再求特征值在数值计算上并不稳定。更稳健的做法是同时对角化S_b和S_w或者求解(S_w εI)^{-1} * S_b其中ε是一个很小的正则化参数。我们的代码将采用后一种方法。3. 项目结构与核心类设计3.1 环境准备与第三方库选择本项目使用纯C17标准旨在展示算法本质因此原则上避免使用Eigen、OpenCV等重型线性代数库。但是为了数据可视化这个最终环节我们将轻度依赖一个图形库来绘制结果。这里我选择matplotlib-cpp它是一个C调用Python matplotlib的桥接库配置简单能轻松绘制高质量的散点图、折线图。对于矩阵运算的核心部分我们将自己实现。你需要准备一个支持C17的编译器如GCC 7 Clang 5 MSVC 2017。Python环境用于matplotlib-cpp后端。安装Python的matplotlib、numpy库。项目目录结构设计如下Cpp_DimensionalityReduction/ ├── include/ │ ├── pca.h │ ├── lda.h │ ├── matrix.h // 自定义矩阵类封装基础运算 │ └── utils.h // 数据读取、预处理等工具函数 ├── src/ │ ├── pca.cpp │ ├── lda.cpp │ ├── matrix.cpp │ └── utils.cpp ├── lib/ // 存放matplotlib-cpp等第三方头文件 ├── datasets/ // 存放示例数据如Iris鸢尾花数据集 ├── examples/ // 示例代码 │ ├── example_pca_iris.cpp │ └── example_lda_iris.cpp └── CMakeLists.txt3.2 核心类接口设计我们的设计遵循清晰、易用的原则。首先在matrix.h中定义一个简单的Matrix类包含基础的行、列、数据存储以及矩阵乘法、转置、切片等操作。这不是工业级实现但足以支撑我们的算法演示。PCA 类设计 (pca.h):class PCA { public: PCA() default; explicit PCA(int n_components); // 拟合模型计算主成分 void fit(const Matrix X); // 将数据转换到主成分空间 Matrix transform(const Matrix X) const; // 拟合并转换的一步操作 Matrix fit_transform(const Matrix X); // 获取主成分特征向量 Matrix get_components() const { return components_; } // 获取解释方差比 std::vectordouble get_explained_variance_ratio() const { return explained_variance_ratio_; } private: int n_components_; Matrix mean_; // 训练数据的均值用于中心化 Matrix components_; // 主成分特征向量每行是一个成分 std::vectordouble explained_variance_ratio_; };LDA 类设计 (lda.h):class LDA { public: LDA() default; explicit LDA(int n_components); // 拟合模型需要数据和对应的标签 void fit(const Matrix X, const std::vectorint y); // 转换数据 Matrix transform(const Matrix X) const; Matrix fit_transform(const Matrix X, const std::vectorint y); // 获取投影矩阵 Matrix get_components() const { return components_; } private: int n_components_; Matrix mean_; // 全局均值 Matrix components_; // LDA投影方向 std::vectorMatrix class_means_; // 每个类别的均值向量 // 核心求解函数 Matrix solve_eigen(const Matrix Sb, const Matrix Sw); };实操心得将fit和transform分离是机器学习库的常见设计模式。它允许你用一个数据集训练fit模型然后用同一个模型去转换训练集或新的测试集transform保证了处理的一致性。fit_transform为便捷操作提供了捷径。4. 核心算法C实现详解4.1 矩阵运算基础实现 (matrix.cpp)由于我们“造轮子”一个可靠的Matrix类是基石。重点实现内存管理使用std::vectordouble连续存储行优先。基本运算加减、乘、转置。矩阵乘法是性能关键我们实现一个简单的三重循环版本对于学习目的足够清晰。切片与访问方便地获取行、列或子矩阵。特殊矩阵生成单位矩阵、零矩阵。这里贴一个矩阵乘法的核心实现注意循环顺序对缓存友好ikj顺序Matrix Matrix::operator*(const Matrix other) const { if (cols_ ! other.rows_) { throw std::invalid_argument(Matrix dimensions mismatch for multiplication.); } Matrix result(rows_, other.cols_, 0.0); for (int i 0; i rows_; i) { for (int k 0; k cols_; k) { double aik data_[i * cols_ k]; if (std::fabs(aik) 1e-10) continue; // 稀疏性优化 for (int j 0; j other.cols_; j) { result.data_[i * other.cols_ j] aik * other.data_[k * other.cols_ j]; } } } return result; }4.2 PCA的核心实现协方差矩阵与雅可比迭代在pca.cpp的fit函数中我们按以下步骤进行步骤1数据中心化。计算每个特征维度的均值然后从原始数据中减去。这是PCA的必要前提确保主成分分析围绕原点进行。void PCA::fit(const Matrix X) { int n_samples X.rows(); int n_features X.cols(); // 计算均值 mean_ compute_mean(X); // 实现一个按列求均值的函数 Matrix X_centered X - mean_; // 需要重载矩阵减法运算符 // ... 后续步骤 }步骤2计算协方差矩阵。公式为C (X_centered^T * X_centered) / (n_samples - 1)。注意这里除的是n_samples-1以获得样本协方差的无偏估计。步骤3特征值分解雅可比迭代法。雅可比法的思想是通过一系列平面旋转正交变换逐步将对称矩阵非对角元素化为0。实现要点迭代直到非对角元素绝对值最大值小于某个阈值如1e-10。每次选择绝对值最大的非对角元素a[p][q]。计算旋转角度theta。对矩阵A和特征向量矩阵V初始为单位阵进行旋转变换。void jacobi_eigen(Matrix A, Matrix eigenvectors, std::vectordouble eigenvalues, double epsilon 1e-10) { int n A.rows(); eigenvectors Matrix::identity(n); int max_iter 1000; for (int iter 0; iter max_iter; iter) { // 寻找最大非对角元 double max_off 0.0; int p 0, q 0; for (int i 0; i n; i) { for (int j i1; j n; j) { if (std::fabs(A(i, j)) max_off) { max_off std::fabs(A(i, j)); p i; q j; } } } if (max_off epsilon) break; // 已对角化 // 计算旋转角度 double app A(p, p); double aqq A(q, q); double apq A(p, q); double theta 0.5 * std::atan2(2 * apq, aqq - app); double c std::cos(theta); double s std::sin(theta); // 对A进行旋转变换 // ... 更新A的第p行、第q行和第p列、第q列 // 对特征向量矩阵V进行同样旋转 for (int k 0; k n; k) { double vkp eigenvectors(k, p); double vkq eigenvectors(k, q); eigenvectors(k, p) c * vkp - s * vkq; eigenvectors(k, q) s * vkp c * vkq; } } // 提取对角线元素作为特征值 for (int i 0; i n; i) { eigenvalues[i] A(i, i); } }步骤4特征值排序与主成分选取。将特征值及对应的特征向量按降序排序。根据设定的n_components_选取前k个特征向量组成投影矩阵components_每行是一个主成分方向。同时计算解释方差比每个主成分的方差特征值占总方差特征值之和的比例这有助于我们决定保留多少个成分。// 对特征值和特征向量进行排序降序 // 使用索引排序避免直接交换特征向量列 std::vectorint idx argsort(eigenvalues, false); // false表示降序 Matrix sorted_components(n_features, n_components_); for (int i 0; i n_components_; i) { int col_idx idx[i]; // 将特征向量矩阵的第col_idx列赋值给sorted_components的第i列 double lambda eigenvalues[col_idx]; total_variance lambda; for (int j 0; j n_features; j) { sorted_components(j, i) eigenvectors(j, col_idx); } explained_variance_ratio_[i] lambda / total_variance_all; // 总方差需提前计算 } components_ sorted_components;步骤5数据转换。transform函数非常简单将中心化后的数据矩阵与投影矩阵components_相乘即可。Matrix PCA::transform(const Matrix X) const { Matrix X_centered X - mean_; return X_centered * components_; // 注意我们的components_是 (n_features x k) }4.3 LDA的核心实现散度矩阵与广义特征值求解LDA的fit函数接收数据X和标签y。步骤1计算全局均值与类别均值。遍历数据统计每个类别的样本索引分别计算各类别的均值向量。步骤2计算类内散度矩阵S_w。S_w是每个类别散度矩阵S_i的和。对于第i类S_i Σ (x - μ_i) * (x - μ_i)^T其中求和针对该类所有样本xμ_i是该类均值。Matrix Sw(n_features, n_features, 0.0); for (int i 0; i n_classes; i) { Matrix Si(n_features, n_features, 0.0); for (int idx : class_indices[i]) { Matrix x X.get_row(idx); // 获取样本行向量需转为列向量 Matrix x_centered x - class_means_[i]; // Si x_centered^T * x_centered; (注意维度x_centered是1xn_features) // 更高效的做法是逐元素累加外积 for (int p 0; p n_features; p) { for (int q 0; q n_features; q) { Si(p, q) x_centered(0, p) * x_centered(0, q); } } } Sw Sw Si; }步骤3计算类间散度矩阵S_b。S_b Σ n_i * (μ_i - μ) * (μ_i - μ)^T其中n_i是第i类的样本数μ是全局均值向量。Matrix Sb(n_features, n_features, 0.0); for (int i 0; i n_classes; i) { Matrix mean_diff class_means_[i] - global_mean; // Sb n_i * (mean_diff^T * mean_diff) for (int p 0; p n_features; p) { for (int q 0; q n_features; q) { Sb(p, q) class_counts[i] * mean_diff(0, p) * mean_diff(0, q); } } }步骤4解决小样本问题S_w奇异。这是LDA实现中最容易出错的地方。当特征维数高、样本少时S_w是奇异的不可逆。我们采用“PCA LDA”的两阶段法先用PCA将数据降至n_features_reduced min(n_features, n_samples - n_classes)维。实际上LDA最多能降到n_classes - 1维。在PCA降维后的空间里计算S_w和S_b此时S_w大概率是满秩的。步骤5求解广义特征值问题。我们求解(S_w εI)^{-1} * S_b的特征向量。其中ε是一个小的正则化参数如1e-4I是单位矩阵目的是稳定求逆。Matrix solve_eigen(const Matrix Sb, const Matrix Sw) { int n Sw.rows(); double epsilon 1e-4; Matrix Sw_reg Sw; for (int i 0; i n; i) { Sw_reg(i, i) epsilon; // 正则化 } // 计算 Sw_reg 的逆这里实现一个简单的高斯消元求逆实际可用更稳定的方法如LU分解 Matrix Sw_inv inverse(Sw_reg); Matrix A Sw_inv * Sb; // 对A进行特征值分解同样使用雅可比法 Matrix eigenvectors; std::vectordouble eigenvalues; jacobi_eigen(A, eigenvectors, eigenvalues); // 选取前k个最大特征值对应的特征向量 // ... 排序和选取逻辑与PCA类似 return components; // 投影矩阵 }步骤6数据转换。与PCA类似但注意LDA的转换通常不包含中心化步骤因为S_b和S_w的计算已经隐含了均值处理或者我们使用全局均值进行中心化。转换公式X_lda (X - global_mean) * W_lda。注意事项LDA投影矩阵W_lda的列向量投影方向一般不是正交的这与PCA不同。此外LDA对数据分布做了很强的假设各类数据服从同协方差的正态分布实际数据若不满足效果可能打折扣。5. 实战演示以鸢尾花(Iris)数据集为例5.1 数据准备与预处理我们使用经典的Iris数据集它包含3类鸢尾花Setosa, Versicolor, Virginica每类50个样本每个样本4个特征花萼长宽、花瓣长宽。这是一个非常适合演示LDA的小数据集。首先实现一个简单的CSV数据读取函数utils.cppstruct Dataset { Matrix data; std::vectorint labels; std::vectorstd::string feature_names; std::vectorstd::string label_names; }; Dataset load_iris(const std::string filename) { std::ifstream file(filename); // ... 解析CSV文件将前4列作为数据最后一列作为标签映射为0,1,2 // 简单起见这里省略具体的文件解析代码可以使用字符串流分割 return dataset; }数据加载后一个良好的习惯是进行标准化Standardization即对每个特征减去其均值并除以标准差使其均值为0方差为1。这能消除不同特征量纲的影响对于基于距离度量的算法包括PCA的方差、LDA的散度尤为重要。void standardize(Matrix X) { int n X.cols(); for (int j 0; j n; j) { double mean compute_column_mean(X, j); double std compute_column_std(X, j, mean); if (std 1e-10) std 1.0; // 防止除零 for (int i 0; i X.rows(); i) { X(i, j) (X(i, j) - mean) / std; } } }5.2 PCA降维与可视化我们先将150x4的数据降到2维以便在平面上可视化。#include pca.h #include utils.h #include matplotlibcpp.h namespace plt matplotlibcpp; int main() { Dataset iris load_iris(datasets/iris.csv); Matrix X iris.data; std::vectorint y iris.labels; // 数据标准化 standardize(X); // PCA降维 PCA pca(2); Matrix X_pca pca.fit_transform(X); // 可视化 std::vectordouble x_coords, y_coords; std::vectorint colors; // 用颜色区分类别 for (int i 0; i X_pca.rows(); i) { x_coords.push_back(X_pca(i, 0)); y_coords.push_back(X_pca(i, 1)); colors.push_back(y[i]); } plt::scatter(x_coords, y_coords, 10, colors); plt::title(Iris Dataset PCA Projection (2D)); plt::xlabel(Principal Component 1); plt::ylabel(Principal Component 2); plt::show(); // 打印解释方差比 auto ratios pca.get_explained_variance_ratio(); std::cout Explained variance ratio (PC1, PC2): ratios[0] , ratios[1] std::endl; std::cout Total variance explained by 2 components: ratios[0] ratios[1] std::endl; return 0; }运行后你将看到一幅散点图。通常PCA的前两个主成分能解释原始数据90%以上的方差。观察散点图你会发现不同类别的点有混合因为PCA不考虑标签信息。5.3 LDA降维与可视化同样降到2维。注意LDA最多能降到n_classes - 1 2维这里正合适。#include lda.h // ... 同样的数据加载和标准化 int main() { // ... 加载数据标准化 LDA lda(2); Matrix X_lda lda.fit_transform(X, y); // 这里需要传入标签y // 可视化代码与PCA类似 // ... plt::title(Iris Dataset LDA Projection (2D)); plt::xlabel(Linear Discriminant 1); plt::ylabel(Linear Discriminant 2); plt::show(); }运行LDA可视化你会看到截然不同的结果三个类别的数据点在二维空间中被清晰地分开了。这正是LDA“类内紧致类间分离”目标的直观体现。对于Iris这种线性可分性较好的数据LDA的投影效果通常优于PCA。5.4 结果对比与算法评估如何定量评估降维效果对于无监督的PCA我们可以看累计解释方差比。对于有监督的LDA一个直接的评估方法是在降维后的低维空间上用一个简单的分类器如K近邻进行分类计算准确率。我们可以实现一个简单的KNN分类器进行评估double knn_accuracy(const Matrix X_train, const std::vectorint y_train, const Matrix X_test, const std::vectorint y_test, int k3) { int correct 0; for (int i 0; i X_test.rows(); i) { // 计算X_test[i]与所有X_train的欧氏距离 std::vectorstd::pairdouble, int distances; // (距离, 训练样本索引) for (int j 0; j X_train.rows(); j) { double dist 0.0; for (int d 0; d X_train.cols(); d) { double diff X_test(i, d) - X_train(j, d); dist diff * diff; } distances.emplace_back(std::sqrt(dist), j); } // 按距离排序取前k个 std::sort(distances.begin(), distances.end()); // 统计k个最近邻的类别标签 std::unordered_mapint, int vote; for (int n 0; n k n distances.size(); n) { int label y_train[distances[n].second]; vote[label]; } // 找出票数最多的类别 int predicted std::max_element(vote.begin(), vote.end(), [](auto a, auto b){ return a.second b.second; })-first; if (predicted y_test[i]) correct; } return static_castdouble(correct) / X_test.rows(); }然后将原始数据、PCA降维后数据、LDA降维后数据分别划分训练集和测试集如80%-20%用KNN评估分类准确率。你大概率会发现原始4维数据准确率不错PCA-2D数据准确率可能略有下降丢失了部分判别信息而LDA-2D数据准确率可能接近甚至超过原始数据因为它直接优化了分类目标。6. 常见问题、调试技巧与性能优化6.1 编译与链接问题问题1使用matplotlib-cpp时链接错误。确保你的CMakeLists.txt正确找到了Python和matplotlib。一个简单的CMake配置如下cmake_minimum_required(VERSION 3.10) project(CppDimensionalityReduction) set(CMAKE_CXX_STANDARD 17) # 查找Python find_package(PythonLibs REQUIRED) include_directories(${PYTHON_INCLUDE_DIRS}) # 添加matplotlib-cpp头文件路径 include_directories(${PROJECT_SOURCE_DIR}/lib) add_executable(pca_demo examples/example_pca_iris.cpp src/pca.cpp src/matrix.cpp src/utils.cpp) target_link_libraries(pca_demo ${PYTHON_LIBRARIES}) # 对于Linux/macOS可能还需要链接pthread、m等 target_link_libraries(pca_demo pthread m)问题2矩阵运算维度不匹配导致崩溃。这是实现中最常见的错误。务必在每一个矩阵乘法、加减法操作前用断言或异常检查维度。在我们的Matrix类运算符重载中第一步就是检查维度兼容性。6.2 算法实现中的数值稳定性问题特征值分解不收敛或结果错误。雅可比法不收敛增加最大迭代次数如5000并确保你的收敛阈值epsilon设置合理如1e-10。检查旋转角度计算atan2(2*apq, aqq-app)是否正确。特征值顺序混乱确保排序逻辑正确。argsort函数应返回按特征值降序排列的索引。LDA中S_w求逆失败即使正则化了如果数据极端病态求逆仍可能不稳定。考虑使用更稳健的求解广义特征值的方法如计算S_w的Cholesky分解S_w L * L^T然后求解L^{-1} * S_b * (L^{-1})^T的标准特征值问题。或者可以转而求解S_b * w λ * (S_w εI) * w的广义特征值问题使用专门的数值库如ARPACK的子空间迭代法。对于学习项目我们的正则化PCA预处理通常足够。6.3 性能瓶颈分析与优化我们实现的雅可比法是O(n³)的且只适用于对称矩阵。对于特征维度很高n1000的数据它会非常慢。在实际生产环境中你绝不会自己写这个。但了解原理后你可以轻松地替换为高效的库调用使用Eigen库推荐用于真实项目// PCA拟合部分可以简化为 #include Eigen/Dense Eigen::MatrixXd X_eigen ...; // 将数据转为Eigen矩阵 Eigen::MatrixXd centered X_eigen.rowwise() - X_eigen.colwise().mean(); Eigen::MatrixXd cov (centered.adjoint() * centered) / (X_eigen.rows() - 1); Eigen::SelfAdjointEigenSolverEigen::MatrixXd eigensolver(cov); if (eigensolver.info() ! Eigen::Success) { /* 处理错误 */ } Eigen::MatrixXd components eigensolver.eigenvectors().rightCols(n_components); // 注意Eigen特征值升序排列Eigen使用高度优化的算法速度快几个数量级且经过充分测试。对于LDAEigen也提供了广义特征值求解器#include Eigen/Dense // 假设Sb, Sw是Eigen矩阵 Eigen::GeneralizedSelfAdjointEigenSolverEigen::MatrixXd ges; ges.compute(Sb, Sw epsilon * Eigen::MatrixXd::Identity(n, n)); Eigen::MatrixXd components ges.eigenvectors().rightCols(n_components); // 取对应最大特征值的向量实操心得自己动手实现算法是学习的绝佳途径它能让你深刻理解每个步骤的细节和陷阱。但在实际项目中除非有极特殊的定制需求否则强烈建议使用Eigen、Intel MKL或Armadillo这些成熟的线性代数库。它们的代码经过千锤百炼在速度、精度和稳定性上远超我们手写的版本。6.4 扩展与进阶思考核PCA (Kernel PCA)对于非线性数据标准的线性PCA无能为力。核PCA通过核函数将数据隐式映射到高维特征空间再在那个空间进行线性PCA。实现关键在于计算核矩阵的特征值分解。增量PCA (Incremental PCA)当数据太大无法一次性装入内存时可以使用增量PCA它允许你分批次处理数据逐步更新主成分。二次判别分析 (QDA)LDA假设所有类别共享同一个协方差矩阵。如果放松这个假设每个类别有自己的协方差矩阵就得到了QDA它是一个非线性分类器决策边界是二次的。流形学习算法如t-SNE、UMAP这些是现代高维数据可视化的利器它们的目标是保持数据点之间的局部邻近关系而非全局方差。你可以尝试用C实现t-SNE虽然挑战不小作为下一个进阶项目。这个从零开始的C PCA与LDA实现项目就像一次深入算法腹地的探险。你不仅得到了两个可用的降维工具更重要的是你亲手触摸了特征值分解、广义特征值问题、矩阵计算这些机器学习基石背后的每一个齿轮。下次当你调用sklearn.decomposition.PCA时你会对屏幕上跳出的结果多一份了然于胸的底气。编程实现算法的过程就是将这些数学公式内化为工程直觉的过程这种收获是单纯阅读论文或使用库函数无法比拟的。