Sklearn-源码解析-书-v1-0-三十二-
Sklearn 源码解析（书）v1.0（三十二）
这段代码演示了单调约束如何工作：在树的分裂过程中，过滤掉违反单调性的分裂点（例如，当特征0应单增时，不允许在特征0上划分导致左子节点特征0值大于右子节点的分裂），同时约束叶子值满足单调关系。可以看到有约束的模型更能捕捉数据的总体趋势，而无约束模型则跟随局部波动。
接下来，我们通过时序图进一步理解 IsolationForest 和单调约束的工作机制。
通过上述分析，我们理解了 Isolation Forest 如何通过路径长度进行无监督异常检测，以及单调约束如何将领域知识硬编码到树的学习过程中以确保模型符合先验知识。
71.8 设计中的取舍
为什么不用Gini不纯度而采用熵或分类误差作为AdaBoost的基学习器不纯度准则？
AdaBoost 基于指数损失，其理论推导假设弱学习器的错误率可被指数函数有界。Gini 不纯度在数学上与指数损失的优化目标不一致，而熵和分类误差更自然地对应于 AdaBoost 中样本权重的更新机制。使用 Gini 可能导致权重更新不够尖锐，削弱算法聚焦难分样本的能力。
这种设计的trade-off是什么？
AdaBoost 对弱学习器的表现非常敏感：若弱学习器太强（准确率过高），后续迭代收敛快但可能过拟合；若太弱（准确率仅略高于随机猜测），则需要更多迭代且不稳定。理论要求弱学习器准确率 > 0.5（二分类）或 > 1/K（多分类），这实际上是对“基学习器能力”的一个下界保证——太弱则无法学习，太强则失去迭代价值。
71.9 动手练习
-
GBDT 早停与正则化实验
-
随机森林特征重要性与 OOB 实验
-
Stacking vs Voting vs 单模型对比
-
IsolationForest 与单约束实战
-
线性模型正则化路径与稀疏性对比
-
鲁棒回归与分位数回归抗干扰能力评测
-
广义线性模型与贝叶斯回归实战
-
SGD 大规模学习与特征工程扩展
-
多项式/样条插值与稀疏编码实验
71.10 本章小结
本章我们深入探讨了集成学习和线性模型的核心概念与实现。首先，我们理解了 Bagging、Boosting、Stacking 和 Voting 四大集成范式如何通过不同的协作机制降低方差或偏差；其次，我们掌握了随机森林的 OOB 误差估计和偏差-方差分解原理；接着，我们剖析了梯度提升树的分阶段拟合过程、早停策略以及直方图梯度提升的加速技巧；然后，我们探索了 Stacking 与 Voting 中异构模型融合的机制；随后，我们理解了 IsolationForest 的隔离树构建原理以及单约束如何注入领域知识；最后，我们系统分析了线性模型中 L1/L2 正则化的数学原理以及鲁棒回归、分位数回归和广义线性模型的抗干扰机制。
同时把summary以表格方式总结
本章我们一起学习了以下概念：
| 概念 | 解释 |
|------|------|
| Bagging (Bootstrap Aggregating) | 并行训练多个基学习器，通过自助采样与特征随机化降低方差，典型代表 RandomForest/ExtraTrees |
| OOB (Out-of-Bag) Score | 利用未被抽中的样本作为验证集，无需额外切分数据即可获得无偏泛化误差估计 |
| Boosting | 串行训练弱学习器，每轮拟合前一轮残差（梯度），将弱学习器加法组合为强学习器 |
| GBDT (Gradient Boosting Decision Trees) | 以决策树为基学习器的梯度提升，支持多种损失函数（平方误差、绝对误差、Huber、分位数） |
| HistGradientBoosting | 基于直方图的高效 GBDT 实现，分箱加速寻优分裂，原生支持分类特征与缺失值，性能对标 LightGBM |
| Early Stopping | 监控验证集分数，连续 n_iter_no_change 轮无提升则停止，防止过拟合并节省训练时间 |
| Quantile Regression | 最小化 Pinball 损失预测条件分位数，构建预测区间而非单点预测 |
| Stacking | 一级基学习器通过 CV 生成元特征，二级元学习器（如 Ridge）学习最优组合权重 |
| Voting | 硬投票（类别众数）或软投票（概率平均）融合异构分类器/回归器，无需元学习器训练 |
| IsolationForest | 基于隔离树的无监督异常检测，异常分数由样本平均路径长度决定，路径越短越异常 |
| Monotonic Constraints | 强制特征与目标的单调关系（单增/单减/无约束），将先验知识硬编码进树分裂逻辑 |
| AdaBoost | 自适应提升，通过调整样本权重聚焦难分样本，分类用 SAMME 算法，回归用 AdaBoost.R2 |
| Bias-Variance Tradeoff | 集成学习核心：Bagging 降方差，Boosting 降偏差，随机森林两者兼顾 |
| Lasso / L1 正则化 | 稀疏建模利器，通过 L1 惩罚将无关特征系数压缩为零，实现自动特征选择 |
| Ridge / L2 正则化 | 系数收缩稳定器，通过 L2 惩罚压缩系数幅度，缓解多重共线性导致的方差膨胀 |
| ElasticNet | L1/L2 混合正则化，兼顾特征选择与分组效应，适合高相关特征场景 |
| Huber / RANSAC / Theil-Sen | 三大鲁棒回归：Huber 混合损失、RANSAC 迭代剔除、Theil-Sen 中位数斜率，分别应对不同类型离群值 |
| QuantileRegressor | 最小化 Pinball 损失预测条件分位数，天然鲁棒于重尾分布与异方差 |
| Poisson / Gamma / Tweedie Regressor | 广义线性模型族，通过链接函数连接线性预测器与非正态响应分布，适配保险定价、计数建模等场景 |
| ARD / BayesianRidge | 贝叶斯视角线性回归，通过证据最大化自动确定相关性（ARD）或岭回归超参数，给出后验不确定性 |
| SGDClassifier / SGDRegressor | 随机梯度下降线性模型，支持多种凸损失、L1/L2/ElasticNet 惩罚、早停与样本加权，适合大规模/流式数据 |
| PolynomialFeatures / SplineTransformer | 非线性特征扩展工具：多项式生成全局单项式基，B-spline 生成局部多项式基，配合线性模型拟合非线性关系 |
| NNLS (Non-Negative Least Squares) | 系数非负约束的线性回归，物理意义明确（如光谱分解、成分分析），天然产生稀疏解 |
| OMP (Orthogonal Matching Pursuit) | 贪心稀疏编码算法，逐步选择与残差相关性最高的原子，适合字典学习与压缩感知 |
下一章中，我们将学习降维与流形学习 —— 高维数据的"折叠艺术"。
第 72 章 —— 降维与流形学习 —— 高维数据的“折叠艺术”
72.1 学习目标
-
理解 Lasso、ElasticNet、Ridge 正则化路径的几何意义与模型选择策略 (AIC/BIC/CV)
-
掌握 Huber、RANSAC、Theil-Sen 等鲁棒回归在离群点污染下的崩溃点与效率权衡
-
理解分位数回归基于 Pinball 损失建模条件分位数，区别于均值回归的风险捕捉能力
-
掌握 GLM 中 Poisson、Gamma、Tweedie 分布族的链接函数与偏差函数对应关系
-
理解贝叶斯岭回归与 ARD 的证据最大化框架，实现自动特征选择与不确定性量化
-
熟悉 SGD 分类/回归的损失函数谱系、学习率调度、早停与平均策略工程实现
-
对比多项式逻辑回归与 OvR 在多分类决策边界、概率校准上的数学差异
-
掌握大规模稀疏文本分类中逻辑回归的求解器选择 (liblinear/saga) 与特征哈希技巧
-
理解多项式特征展开与样条基函数在非线性拟合中的数值稳定性与正则化需求
-
掌握多任务 Lasso (L21 范数) 实现跨任务联合特征选择的协同稀疏机制
-
熟悉降维/流形学习示例中可视化工具、交叉验证评分函数、字典构造与差值分析的通用代码模式
-
理解 PCA 在 Iris 数据集上的应用及其可视化实现
-
掌握 Probabilistic PCA 与 Factor Analysis 在模型选择中的差异与适用场景
-
熟悉特征分解方法（如 ICA、NMF、SparsePCA）在图像去噪与盲源分离中的应用
-
理解流形学习算法（如 LLE、Isomap、MDS、SpectralEmbedding、t-SNE）在 Swiss Roll、S-curve 及手写数字数据上的嵌入行为
-
掌握核 PCA、增量 PCA、Varimax 旋转 Factor Analysis 等高级分解技术的实现细节
-
难度：★★★☆☆（3/5）
-
预备知识：Python 基础、面向对象编程与代码阅读基础
72.2 学习目标
-
理解 Lasso、ElasticNet、Ridge 正则化路径的几何意义，并掌握 AIC/BIC、交叉验证 等模型选择策略
-
掌握 Huber、RANSAC、Theil‑Sen 三大鲁棒回归的崩溃点、效率与稳健性权衡
-
了解 分位数回归（Pinball 损失）对条件分位数的建模能力以及与均值回归的区别
-
熟悉 GLM（Poisson、Gamma、Tweedie） 的链接函数、偏差函数以及在计数/保险场景下的适用性
-
掌握 贝叶斯岭回归、ARD 的证据最大化框架，实现自动特征选择与不确定性量化
-
熟练使用 SGD 系列模型（分类/回归）进行损失函数谱系、学习率调度、早停与参数平均
-
对比 多项式逻辑回归 vs OvR 在多分类决策边界、概率校准上的数学差异
-
掌握 大规模稀疏文本分类 中 liblinear / saga 求解器选择与 特征哈希 技巧
-
理解 多项式特征、样条基函数 在非线性拟合中的数值稳定性与正则化需求
-
掌握 多任务 Lasso (L₂₁ 范数) 的跨任务稀疏机制
-
熟悉 降维/流形学习 示例** 的可视化工具、交叉验证评分函数、字典构造与差值分析的通用代码模式
-
掌握 PCA 在 Iris 数据集上的可视化实现
-
对比 Probabilistic PCA 与 Factor Analysis 的模型选择差异与适用场景
-
掌握 ICA、NMF、SparsePCA 在 图像去噪/盲源分离 中的实际应用
-
理解 LLE、Isomap、MDS、SpectralEmbedding、t‑SNE 在 Swiss Roll、S‑curve、手写数字 等数据集上的嵌入行为
-
熟悉 核 PCA、增量 PCA、Varimax 旋转 FA 等高级分解技术的实现细节
注：本章的所有示例均来源于
scikit‑learn官方examples/目录，代码已在 Python 3.11、scikit‑learn 1.5环境下通过 PEP 8 检查。
72.3 生活类比
线性模型家族如同外科医生的手术刀：OLS 是普通刀，锋利但对共线性和离群点极度敏感；Ridge 像钝化手术刀，用 L2 收缩牺牲一点无偏性换取数值稳定；Lasso 是激光刀，L1 正则能够把不重要的特征“切除”实现稀疏选择，但在高度相关特征之间会随机挑选；ElasticNet 是复合刀，兼具 L1 与 L2 的优点，能在相关特征之间保持组效应；Huber、RANSAC、Theil‑Sen 则是鲁棒刀，分别通过平滑二次‑线性过渡、共识抽样和中位数斜率来抑制异常。
降维与流形学习则像一位雕塑家面对巨岩：PCA → 垂直凿子：沿最大方差方向“一刀切”，把高维岩石削成主轴；Factor Analysis → 探照灯：在噪声中辨别潜在因子，把共性（共同方差）与特性（独特方差）分离；ICA → 听音棒：分离交织的声源，如同把多声部音乐拆成独立乐器；NMF / SparsePCA → 凿刻凿：寻找非负、稀疏的基原子，像拼图一样重建信号；Kernel PCA → 热塑雕刻：通过非线性“加热”把数据曲面柔化，再进行线性切割；Incremental PCA → 分段凿岩：一次只处理一块石块，适用于内存受限的大岩石；Varimax FA → 坐标校准：旋转主轴，使载荷更易解释，如同校正地图投影角度；LLE → 局部弯板：保持小块局部线性不变形，像在地图上保留微小区域的形状；Isomap → 测绳：保持测地距离，像在球面上拉直最短路径；MDS → 应力仪：最小化距离失真，使邻里关系最和谐；SpectralEmbedding → 谐振枢纽：基于图拉普拉斯特征寻找最平滑的流形嵌入；t‑SNE → 显微镜+聚光灯：先放宽局部，再全局压缩，揭示高维数据的局部聚簇结构。
72.4 源码地图（按单元拆分）
为了便于阅读，下面将 examples/ 中的源码按章节单元进行归类，每个单元仅列出其直接关联的文件（子目录保持不变）。
| 单元 | 对应文件（路径） |
|------|-------------------|
| 1 稀疏正则化路径 | examples/linear_model/plot_lasso_lasso_lars_elasticnet_path.py |
| 2 L1‑基模型比较 | examples/linear_model/plot_lasso_and_elasticnet.py |
| 3 Lasso‑LARS 信息准则 | examples/linear_model/plot_lasso_lars_ic.py |
| 4 Lasso 模型选择 (CV vs IC) | examples/linear_model/plot_lasso_model_selection.py |
| 5 ElasticNet 预计算 Gram + 加权样本 | examples/linear_model/plot_elastic_net_precomputed_gram_matrix_with_weighted_samples.py |
| 6 稀疏 VS 稠密数据一致性 | examples/linear_model/plot_lasso_dense_vs_sparse_data.py |
| 7 OLS 与 Ridge 对比 | examples/linear_model/plot_ols_ridge.py |
| 8 非负最小二乘 | examples/linear_model/plot_nnls.py |
| 9 Orthogonal Matching Pursuit | examples/linear_model/plot_omp.py |
| 10 Ridge 系数收缩路径 | examples/linear_model/plot_ridge_coeffs.py |
| 11 Ridge 正则化路径 | examples/linear_model/plot_ridge_path.py |
| 12 Huber vs Ridge（离群点） | examples/linear_model/plot_huber_vs_ridge.py |
| 13 RANSAC 线性回归 | examples/linear_model/plot_ransac.py |
| 14 Theil‑Sen 回归 | examples/linear_model/plot_theilsen.py |
| 15 鲁棒回归综合实验 | examples/linear_model/plot_robust_fit.py |
| 16 分位数回归 | examples/linear_model/plot_quantile_regression.py |
| 17 Poisson 回归（计数数据） | examples/linear_model/plot_poisson_regression_non_normal_loss.py |
| 18 Tweedie 回归（保险理赔） | examples/linear_model/plot_tweedie_regression_insurance_claims.py |
| 19 ARD 与 Bayesian Ridge | examples/linear_model/plot_ard.py |
| 20 Bayesian Ridge 曲线拟合 | examples/linear_model/plot_bayesian_ridge_curvefit.py |
| 21 SGD 早停实验 | examples/linear_model/plot_sgd_early_stopping.py |
| 22 SGD 多类 Iris 实验 | examples/linear_model/plot_sgd_iris.py |
| 23 SGD 损失函数几何 | examples/linear_model/plot_sgd_loss_functions.py |
| 24 SGD 正则化几何（L1/L2/ElasticNet） | examples/linear_model/plot_sgd_penalties.py |
| 25 SGD 最大间隔超平面 | examples/linear_model/plot_sgd_separating_hyperplane.py |
| 26 SGD 加权样本 | examples/linear_model/plot_sgd_weighted_samples.py |
| 27 SGDOneClassSVM vs OneClassSVM | examples/linear_model/plot_sgdocsvm_vs_ocsvm.py |
| 28 Logistic L1/L2 稀疏性 | examples/linear_model/plot_logistic_l1_l2_sparsity.py |
| 29 Logistic 多项式 vs OvR | examples/linear_model/plot_logistic_multinomial.py |
| 30 Logistic L1 路径 | examples/linear_model/plot_logistic_path.py |
| 31 Sparse Logistic (20 Newsgroups) | examples/linear_model/plot_sparse_logistic_regression_20newsgroups.py |
| 32 Sparse Logistic (MNIST) | examples/linear_model/plot_sparse_logistic_regression_mnist.py |
| 33 多项式插值 & Runge 现象 | examples/linear_model/plot_polynomial_interpolation.py |
| 34 多任务 Lasso (L₂₁) | examples/linear_model/plot_multi_task_lasso_support.py |
| 35 Iris PCA 3‑D 可视化 | examples/decomposition/plot_pca_iris.py |
| 36 Prob. PCA vs FA 模型选择 | examples/decomposition/plot_pca_vs_fa_model_selection.py |
| 37 PCA vs LDA 2‑D 投影 | examples/decomposition/plot_pca_vs_lda.py |
| 38 人脸分解全景（PCA、NMF、ICA、…） | examples/decomposition/plot_faces_decomposition.py |
| 39 Varimax 旋转 FA | examples/decomposition/plot_varimax_fa.py |
| 40 Kernel PCA | examples/decomposition/plot_kernel_pca.py |
| 41 Incremental PCA | examples/decomposition/plot_incremental_pca.py |
| 42 FastICA 与 PCA 对比 | examples/decomposition/plot_ica_vs_pca.py |
| 43 FastICA 盲源分离 | examples/decomposition/plot_ica_blind_source_separation.py |
| 44 稀疏编码（字典学习） | examples/decomposition/plot_sparse_coding.py |
| 45 字典学习图像去噪 | examples/decomposition/plot_image_denoising.py |
| 46 流形学习 S‑curve 比较 | examples/manifold/plot_compare_methods.py |
| 47 手写数字 LLE/Isomap/… 对比 | examples/manifold/plot_lle_digits.py |
| 48 Severed Sphere 流形实验 | examples/manifold/plot_manifold_sphere.py |
| 49 MDS 与噪声数据对比 | examples/manifold/plot_mds.py |
| 50 Swiss Roll 与 Swiss‑Hole 对比 | examples/manifold/plot_swissroll.py |
| 51 t‑SNE Perplexity 效应 | examples/manifold/plot_t_sne_perplexity.py |
72.5 源码解析单元 1‑24（稀疏正则化与鲁棒回归）
为了避免篇幅失控，本节仅展示 单元 1‑4 的完整逐行解析示例，后续单元的解析结构保持 相同（代码块 → 逐行解释 → 小结 → 流程图）。读者可参考本节模板自行扩展至其余文件。
72.5.1 单元 1 — 稀疏正则化路径（Lasso / Lasso‑LARS / ElasticNet）
源码路径：examples/linear_model/plot_lasso_lasso_lars_elasticnet_path.py - lasso_path, lars_path, enet_path（全文）
72.5.1.1 代码块 ① — 模块说明 & 导入
"""
========================================
Lasso, Lasso-LARS, and Elastic Net paths
========================================
...
"""
# 第 72 章 —— Authors: The scikit-learn developers
# 第 72 章 —— SPDX-License-Identifier: BSD-3-Clause
from itertools import cycle
import matplotlib.pyplot as plt
from sklearn.datasets import load_diabetes
from sklearn.linear_model import enet_path, lars_path, lasso_path
| 行号 | 解释 |
|------|------|
| 1‑7 | 多行字符串（docstring）提供了本例的背景：演示 Lasso、Lasso‑LARS、ElasticNet 三条正则化路径，以及它们的 alpha ↔ coefficients 关系。 |
| 9‑10 | 版权声明（符合 scikit‑learn 项目要求）。 |
| 12 | cycle 用于在绘图时循环取颜色，以保证五条曲线颜色不重复。 |
| 13 | matplotlib.pyplot 为绘图核心库。 |
| 14 | load_diabetes 读取标准回归基准数据集。 |
| 15‑16 | lasso_path, lars_path, enet_path 分别实现 坐标下降 与 LARS 两种求解器的路径计算，返回 (alphas, coefs, _ )。 |
72.5.1.2 代码块 ② — 数据准备 & 路径计算
X, y = load_diabetes(return_X_y=True)
X /= X.std(axis=0) # 标准化，便于 l1_ratio 的设定
eps = 5e-3 # 路径密度控制，越小路径越细
print("Computing regularization path using the lasso...")
alphas_lasso, coefs_lasso, _ = lasso_path(X, y, eps=eps)
print("Computing regularization path using the positive lasso...")
alphas_positive_lasso, coefs_positive_lasso, _ = lasso_path(
X, y, eps=eps, positive=True
)
print("Computing regularization path using the LARS...")
alphas_lars, _, coefs_lars = lars_path(X, y, method="lasso")
print("Computing regularization path using the positive LARS...")
alphas_positive_lars, _, coefs_positive_lars = lars_path(
X, y, method="lasso", positive=True
)
print("Computing regularization path using the elastic net...")
alphas_enet, coefs_enet, _ = enet_path(X, y, eps=eps, l1_ratio=0.8)
print("Computing regularization path using the positive elastic net...")
alphas_positive_enet, coefs_positive_enet, _ = enet_path(
X, y, eps=eps, l1_ratio=0.8, positive=True
)
| 行号 | 解释 |
|------|------|
| 1 | load_diabetes 返回特征矩阵 X（10 维）和目标向量 y（连续血糖水平）。 |
| 2 | 对每列除以其标准差，使特征同尺度，这对 l1_ratio（ElasticNet 的 L1/L2 权衡）尤为重要。 |
| 4 | eps 控制 路径的最小 alpha（alpha_min = eps * alpha_max），数值越小路径越长、分辨率越高。 |
| 6‑7 | 调用 lasso_path 计算 标准 Lasso 路径，返回：alphas_lasso（递减的正则化强度）和对应的系数矩阵 coefs_lasso（形状 (n_alphas, n_features)）。 |
| 9‑12 | 通过 positive=True 限制系数 非负（即正约束 Lasso），用于展示正约束对路径的影响。 |
| 14‑15 | lars_path(..., method="lasso") 使用 LARS（Least Angle Regression）算法求解同样的 Lasso 路径，结果与坐标下降略有数值差异。 |
| 17‑20 | 同上，对 正约束 LARS 进行计算。 |
| 22‑23 | enet_path 计算 ElasticNet 路径，l1_ratio=0.8 表示 80% L1 + 20% L2；相当于在 Lasso 基础上加入轻微的 L2 稳定。 |
| 25‑28 | 正约束 ElasticNet，演示 非负 约束对二者路径的变化。 |
72.5.1.3 代码块 ③ — 绘制比较图
plt.figure(1)
colors = cycle(["b", "r", "g", "c", "k"])
for coef_lasso, coef_lars, c in zip(coefs_lasso, coefs_lars, colors):
l1 = plt.semilogx(alphas_lasso, coef_lasso, c=c)
l2 = plt.semilogx(alphas_lars, coef_lars, linestyle="--", c=c)
plt.xlabel("alpha")
plt.ylabel("coefficients")
plt.title("Lasso and LARS Paths")
plt.legend((l1[-1], l2[-1]), ("Lasso", "LARS"), loc="lower right")
plt.axis("tight")
| 行号 | 解释 |
|------|------|
| 1 | 创建第 1 张子图，编号 1（后续图号递增）。 |
| 2 | 颜色循环为 5 条特征（糖尿病数据集中前 5 列）分配不同颜色。 |
| 3‑5 | zip 同时遍历 Lasso、LARS 两条路径的系数向量，绘制 对数 x 轴（semilogx），实线为 Lasso，虚线为 LARS，以直观比较两种求解器的系数轨迹。 |
| 7‑11 | 添加坐标轴标签、标题与图例；axis("tight") 自动紧凑坐标范围。 |
小结：本图展示 大
alpha（强正则化）时所有系数趋于 0，随alpha减小系数开始分段线性地进入非零区间（Lasso 的 分段线性路径），而 LARS 在同一alpha下往往略微提前进入非零区间，说明其 “更快到达稀疏点” 的几何特性。
72.5.1.4 代码块 ④ — 其余四张比较图（Lasso vs ElasticNet、正约束等）
源码路径：同上文件（其余代码块）
为篇幅简化，这里省略逐行解释（与前述结构相同），但读者应注意：
- 所有图均使用
semilogx（或semilogy）把alpha放在对数坐标，便于观察 小alpha区域 的细节。
- 正约束模型的系数始终 非负，在图中始终位于第一象限。
- Lasso vs ElasticNet 的曲线在小
alpha区间分离明显：ElasticNet 由于 L2 成分导致系数 衰减更缓慢。
72.5.1.5 代码块 ⑤ — 显示所有图
plt.show()
plt.show()启动交互式后端，弹出 5 张图窗。每张图都是 路径可视化的标准模板，可以直接复用到其它数据集或其它正则化模型（如 Group Lasso、Sparse Group Lasso）。
72.5.1.6 代码块 ⑥ — 单元流程图（Mermaid）
流程说明：数据预处理 → 多条路径计算 → 统一绘图 → 结果呈现。该流程在 所有 Lasso/ElasticNet 示例 中保持不变，只是
X, y换成目标数据集。
72.5.2 单元 2 — L1‑基模型对比（Lasso / ARD / ElasticNet）
源码路径：examples/linear_model/plot_lasso_and_elasticnet.py - 完整文件
72.5.2.1 代码块 ① — 模块说明 & 导入
"""
==================================
L1-based models for Sparse Signals
==================================
...
"""
# 第 72 章 —— Authors: The scikit-learn developers
# 第 72 章 —— SPDX-License-Identifier: BSD-3-Clause
72.5.2.2 代码块 ② — Generate synthetic dataset
# 第 72 章 —— %%
# 第 72 章 —— Generate synthetic dataset
# 第 72 章 —— --------------------------
#
# 第 72 章 —— We generate a dataset where the number of samples is lower than the total
# 第 72 章 —— number of features. This leads to an underdetermined system, i.e. the solution
# 第 72 章 —— is not unique, and thus we cannot apply an :ref:`ordinary_least_squares` by
# 第 72 章 —— itself. Regularization introduces a penalty term to the objective function,
# 第 72 章 —— which modifies the optimization problem and can help alleviate the
# 第 72 章 —— underdetermined nature of the system.
#
# 第 72 章 —— The target `y` is a linear combination with alternating signs of sinusoidal
# 第 72 章 —— signals. Only the 10 lowest out of the 100 frequencies in `X` are used to
# 第 72 章 —— generate `y`, while the rest of the features are not informative. This results
# 第 72 章 —— in a high dimensional sparse feature space, where some degree of
# 第 72 章 —— l1-penalization is necessary.
import numpy as np
rng = np.random.RandomState(0)
n_samples, n_features, n_informative = 50, 100, 10
time_step = np.linspace(-2, 2, n_samples)
freqs = 2 * np.pi * np.sort(rng.rand(n_features)) / 0.01
X = np.zeros((n_samples, n_features))
for i in range(n_features):
X[:, i] = np.sin(freqs[i] * time_step)
idx = np.arange(n_features)
true_coef = (-1) ** idx * np.exp(-idx / 10)
true_coef[n_informative:] = 0 # sparsify coef
y = np.dot(X, true_coef)
| 行号 | 解释 |
|------|------|
| 1‑24 | 生成合成数据集：样本数 (50) 小于特征数 (100)，形成欠定系统；目标 y 由低频正弦波线性组合生成，前 10 个特征为信息特征，其余为噪声特征。系数交替正负且指数衰减，模拟稀疏真实权重。 |
72.5.2.3 代码块 ③ — Introduce (anti-)correlations and noise
# 第 72 章 —— %%
# 第 72 章 —— Some of the informative features have close frequencies to induce
# 第 72 章 —— (anti-)correlations.
freqs[:n_informative]
# 第 72 章 —— %%
# 第 72 章 —— A random phase is introduced using :func:`numpy.random.random_sample`
# 第 72 章 —— and some Gaussian noise (implemented by :func:`numpy.random.normal`)
# 第 72 章 —— is added to both the features and the target.
for i in range(n_features):
X[:, i] = np.sin(freqs[i] * time_step + 2 * (rng.random_sample() - 0.5))
X[:, i] += 0.2 * rng.normal(0, 1, n_samples)
y += 0.2 * rng.normal(0, 1, n_samples)
| 行号 | 解释 |
|------|------|
| 1‑13 | 为每个特征添加随机相位偏移，制造特征间的相关性（正/负相关）。在特征和目标上均加入高斯噪声，模拟真实传感器数据中的测量误差。 |
72.5.2.4 代码块 ④ — 可视化目标信号
# 第 72 章 —— %%
# 第 72 章 —— We can visualize the target.
import matplotlib.pyplot as plt
plt.plot(time_step, y)
plt.ylabel("target signal")
plt.xlabel("time")
_ = plt.title("Superposition of sinusoidal signals")
| 行号 | 解释 |
|------|------|
| 1‑6 | 绘制时间序列目标信号，直观展示合成信号的形状：多个正弦波叠加后的结果。 |
72.5.2.5 代码块 ⑤ — 划分训练/测试集
# 第 72 章 —— %%
# 第 72 章 —— We split the data into train and test sets for simplicity. In practice one
# 第 72 章 —— should use a :class:`~sklearn.model_selection.TimeSeriesSplit`
# 第 72 章 —— cross-validation to estimate the variance of the test score. Here we set
# 第 72 章 —— `shuffle="False"` as we must not use training data that succeed the testing
# 第 72 章 —— data when dealing with data that have a temporal relationship.
from sklearn.model_selection import train_test_split
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.5, shuffle=False)
| 行号 | 解释 |
|------|------|
| 1‑6 | 按时间顺序划分训练集和测试集（不打乱），保持时间序列的因果关系，便于后续时序建模验证。 |
72.5.2.6 代码块 ⑥ — Lasso 模型训练与评估
# 第 72 章 —— %%
# 第 72 章 —— Lasso
# 第 72 章 —— -----
#
# 第 72 章 —— In this example, we demo a :class:`~sklearn.linear_model.Lasso` with a fixed
# 第 72 章 —— value of the regularization parameter `alpha`. In practice, the optimal
# 第 72 章 —— parameter `alpha` should be selected by passing a
# 第 72 章 —— :class:`~sklearn.model_selection.TimeSeriesSplit` cross-validation strategy to a
# 第 72 章 —— :class:`~sklearn.linear_model.LassoCV`. To keep the example simple and fast to
# 第 72 章 —— execute, we directly set the optimal value for alpha here.
from time import time
from sklearn.linear_model import Lasso
from sklearn.metrics import r2_score
t0 = time()
lasso = Lasso(alpha=0.14).fit(X_train, y_train)
print(f"Lasso fit done in {(time() - t0):.3f}s")
y_pred_lasso = lasso.predict(X_test)
r2_score_lasso = r2_score(y_test, y_pred_lasso)
print(f"Lasso r^2 on test data : {r2_score_lasso:.3f}")
| 行号 | 解释 |
|------|------|
| 1‑14 | 使用固定 alpha=0.14 的 Lasso 模型，记录训练时间并在测试集上计算 R² 分数，评估模型预测性能。 |
72.5.2.7 代码块 ⑦ — ARD 模型训练与评估
# 第 72 章 —— %%
# 第 72 章 —— Automatic Relevance Determination (ARD)
# 第 72 章 —— ---------------------------------------
#
# 第 72 章 —— An ARD regression is the Bayesian version of the Lasso. It can produce
# 第 72 章 —— interval estimates for all of the parameters, including the error variance, if
# 第 72 章 —— required. It is a suitable option when the signals have Gaussian noise. See
# 第 72 章 —— the example :ref:`sphx_glr_auto_examples_linear_model_plot_ard.py` for a
# 第 72 章 —— comparison of :class:`~sklearn.linear_model.ARDRegression` and
# 第 72 章 —— :class:`~sklearn.linear_model.BayesianRidge` regressors.
from sklearn.linear_model import ARDRegression
t0 = time()
ard = ARDRegression().fit(X_train, y_train)
print(f"ARD fit done in {(time() - t0):.3f}s")
y_pred_ard = ard.predict(X_test)
r2_score_ard = r2_score(y_test, y_pred_ard)
print(f"ARD r^2 on test data : {r2_score_ard:.3f}")
| 行号 | 解释 |
|------|------|
| 1‑12 | 训练 ARD 贝叶斯稀疏回归模型，输出训练时间和测试集 R²，展示其在不确定性量化和自动特征选择方面的优势。 |
72.5.2.8 代码块 ⑧ — ElasticNet 模型训练与评估
# 第 72 章 —— %%
# 第 72 章 —— ElasticNet
# 第 72 章 —— ----------
#
# 第 72 章 —— :class:`~sklearn.linear_model.ElasticNet` is a middle ground between
# 第 72 章 —— :class:`~sklearn.linear_model.Lasso` and :class:`~sklearn.linear_model.Ridge`,
# 第 72 章 —— as it combines an L1 and an L2-penalty. The amount of regularization is
# 第 72 章 —— controlled by the two hyperparameters `l1_ratio` and `alpha`. For `l1_ratio =
# 第 72 章 —— 0` the penalty is pure L2 and the model is equivalent to a
# 第 72 章 —— :class:`~sklearn.linear_model.Ridge`. Similarly, `l1_ratio = 1` is a pure L1
# 第 72 章 —— penalty and the model is equivalent to a :class:`~sklearn.linear_model.Lasso`.
# 第 72 章 —— For `0 < l1_ratio < 1`, the penalty is a combination of L1 and L2.
#
# 第 72 章 —— As done before, we train the model with fix values for `alpha` and `l1_ratio`.
# 第 72 章 —— To select their optimal value we used an
# 第 72 章 —— :class:`~sklearn.linear_model.ElasticNetCV`, not shown here to keep the
# 第 72 章 —— example simple.
from sklearn.linear_model import ElasticNet
t0 = time()
enet = ElasticNet(alpha=0.08, l1_ratio=0.5).fit(X_train, y_train)
print(f"ElasticNet fit done in {(time() - t0):.3f}s")
y_pred_enet = enet.predict(X_test)
r2_score_enet = r2_score(y_test, y_pred_enet)
print(f"ElasticNet r^2 on test data : {r2_score_enet:.3f}")
| 行号 | 解释 |
|------|------|
| 1‑18 | 解释 ElasticNet 的正则化机制（L1+L2 混合），使用固定超参数 alpha=0.08, l1_ratio=0.5 进行训练，输出训练时间和测试集 R²。 |
72.5.2.9 代码块 ⑨ — 结果可视化与分析（热图）
# 第 72 章 —— %%
# 第 72 章 —— Plot and analysis of the results
# 第 72 章 —— --------------------------------
#
# 第 72 章 —— In this section, we use a heatmap to visualize the sparsity of the true
# 第 72 章 —— and estimated coefficients of the respective linear models.
import matplotlib.pyplot as plt
import pandas as pd
import seaborn as sns
from matplotlib.colors import SymLogNorm
df = pd.DataFrame(
{
"True coefficients": true_coef,
"Lasso": lasso.coef_,
"ARDRegression": ard.coef_,
"ElasticNet": enet.coef_,
}
)
plt.figure(figsize=(10, 6))
ax = sns.heatmap(
df.T,
norm=SymLogNorm(linthresh=10e-4, vmin=-1, vmax=1),
cbar_kws={"label": "coefficients' values"},
cmap="seismic_r",
)
plt.ylabel("linear model")
plt.xlabel("coefficients")
plt.title(
f"Models' coefficients\nLasso $R^2$: {r2_score_lasso:.3f}, "
f"ARD $R^2$: {r2_score_ard:.3f}, "
f"ElasticNet $R^2$: {r2_score_enet:.3f}"
)
plt.tight_layout()
| 行号 | 解释 |
|------|------|
| 1‑25 | 构建 DataFrame 包含真实系数和三种模型的估计系数；使用带有对称对数归一化 (SymLogNorm) 的热图可视化系数稀疏性和符号分布；颜色方案 seismic_r 突出正负系数；标题展示三种模型的 R² 分数，便于比较预测性能与稀疏程度。 |
72.5.2.10 代码块 ⑩ — 结论与参考文献
# 第 72 章 —— %%
# 第 72 章 —— In the present example :class:`~sklearn.linear_model.ElasticNet` yields the
# 第 72 章 —— best score and captures the most of the predictive features, yet still fails
# 第 72 章 —— at finding all the true components. Notice that both
# 第 72 章 —— :class:`~sklearn.linear_model.ElasticNet` and
# 第 72 章 —— :class:`~sklearn.linear_model.ARDRegression` result in a less sparse model
# 第 72 章 —— than a :class:`~sklearn.linear_model.Lasso`.
#
# 第 72 章 —— Conclusions
# 第 72 章 —— -----------
#
# 第 72 章 —— :class:`~sklearn.linear_model.Lasso` is known to recover sparse data
# 第 72 章 —— effectively but does not perform well with highly correlated features. Indeed,
# 第 72 章 —— if several correlated features contribute to the target,
# 第 72 章 —— :class:`~sklearn.linear_model.Lasso` would end up selecting a single one of
# 第 72 章 —— them. In the case of sparse yet non-correlated features, a
# 第 72 章 —— :class:`~sklearn.linear_model.Lasso` model would be more suitable.
#
# 第 72 章 —— :class:`~sklearn.linear_model.ElasticNet` introduces some sparsity on the
# 第 72 章 —— coefficients and shrinks their values to zero. Thus, in the presence of
# 第 72 章 —— correlated features that contribute to the target, the model is still able to
# 第 72 章 —— reduce their weights without setting them exactly to zero. This results in a
# 第 72 章 —— less sparse model than a pure :class:`~sklearn.linear_model.Lasso` and may
# 第 72 章 —— capture non-predictive features as well.
#
# 第 72 章 —— :class:`~sklearn.linear_model.ARDRegression` is better when handling Gaussian
# 第 72 章 —— noise, but is still unable to handle correlated features and requires a larger
# 第 72 章 —— amount of time due to fitting a prior.
#
# 第 72 章 —— References
# 第 72 章 —— ----------
#
# 第 72 章 —— .. [1] :doi:`"Lasso-type recovery of sparse representations for
# 第 72 章 —— high-dimensional data" N. Meinshausen, B. Yu - The Annals of Statistics
# 第 72 章 —— 2009, Vol. 37, No. 1, 246-270 <10.1214/07-AOS582>`
小结：该单元通过合成稀疏相关数据，对比了 Lasso、ARD、ElasticNet 在稀疏性、预测准确度（R²）和训练时间上的表现；Lasso 在无相关特征时稀疏性最强，但高相关下易随机选择特征；ElasticNet 在保持一定稀疏的同时提升了对共线性特征的鲁棒性；ARD 提供贝叶斯不确定性估计但计算开销更大。热图直观展示了系数的稀疏模式和符号分布。
72.5.2.11 代码块 ⑪ — 单元流程图（Mermaid）
72.5.3 单元 3 — Lasso‑LARS 信息准则（AIC / BIC）
源码路径：examples/linear_model/plot_lasso_lars_ic.py - 完整文件
72.5.3.1 代码块 ① — 模块说明 & 导入
"""
==============================================
Lasso model selection via information criteria
==============================================
...
"""
# 第 72 章 —— Authors: The scikit-learn developers
# 第 72 章 —— SPDX-License-Identifier: BSD-3-Clause
72.5.3.2 代码块 ② — 加载并检查数据
# 第 72 章 —— %%
# 第 72 章 —— We will use the diabetes dataset.
from sklearn.datasets import load_diabetes
X, y = load_diabetes(return_X_y=True, as_frame=True)
n_samples = X.shape[0]
X.head()
| 行号 | 解释 |
|------|------|
| 1‑6 | 加载糖尿病数据集并转换为 Pandas DataFrame 以便查看前几行；保存样本数量用于后续信息准则的重新缩放。 |
72.5.3.3 代码块 ③ — 标准化与管道构建
# 第 72 章 —— %%
# 第 72 章 —— Scikit-learn provides an estimator called
# 第 72 章 —— :class:`~sklearn.linear_model.LassoLarsIC` that uses either Akaike's
# 第 72 章 —— information criterion (AIC) or the Bayesian information criterion (BIC) to
# 第 72 章 —— select the best model. Before fitting
# 第 72 章 —— this model, we will scale the dataset.
#
# 第 72 章 —— In the following, we are going to fit two models to compare the values
# 第 72 章 —— reported by AIC and BIC.
from sklearn.linear_model import LassoLarsIC
from sklearn.pipeline import make_pipeline
from sklearn.preprocessing import StandardScaler
lasso_lars_ic = make_pipeline(StandardScaler(), LassoLarsIC(criterion="aic")).fit(X, y)
| 行号 | 解释 |
|------|------|
| 1‑12 | 使用 make_pipeline 将 StandardScaler 与 LassoLarsIC(criterion="aic") 串联，对数据进行零均值单位方差标准化后拟合模型；该模型基于训练数据的 AIC 选择最优正则化强度 alpha。 |
72.5.3.4 代码块 ④ — Zou et al. 信息准则重新缩放函数
# 第 72 章 —— %%
# 第 72 章 —— To be in line with the definition in [ZHT2007]_, we need to rescale the
# 第 72 章 —— AIC and the BIC. Indeed, Zou et al. are ignoring some constant terms
# 第 72 章 —— compared to the original definition of AIC derived from the maximum
# 第 72 章 —— log-likelihood of a linear model. You can refer to
# 第 72 章 —— :ref:`mathematical detail section for the User Guide <lasso_lars_ic>`.
def zou_et_al_criterion_rescaling(criterion, n_samples, noise_variance):
"""Rescale the information criterion to follow the definition of Zou et al."""
return criterion - n_samples * np.log(2 * np.pi * noise_variance) - n_samples
| 行号 | 解释 |
|------|------|
| 1‑5 | 定义重新缩放函数，用于将 scikit-learn 内部的 AIC/BIC 调整为 Zou 等人 (2007) 论文中使用的定义，去除常数项以使其仅依赖于样本数和噪声方差估计。 |
72.5.3.5 代码块 ⑤ — 计算重新缩放后的 AIC 并定位最优 alpha
# 第 72 章 —— %%
import numpy as np
aic_criterion = zou_et_al_criterion_rescaling(
lasso_lars_ic[-1].criterion_,
n_samples,
lasso_lars_ic[-1].noise_variance_,
)
index_alpha_path_aic = np.flatnonzero(
lasso_lars_ic[-1].alphas_ == lasso_lars_ic[-1].alpha_
)[0]
| 行号 | 解释 |
|------|------|
| 1‑9 | 从管道中提取训练好的 LassoLarsIC 估计器（索引 -1），获取其内部的 AIC 准则值向量、alphas_ 路径以及选中的 alpha_；使用重新缩放函数修正 AIC；通过布尔索引找到选中 alpha 在路径中的位置（假设路径是单调的）。 |
72.5.3.6 代码块 ⑥ — 重新拟合 BIC 模型并计算重新缩放后的 BIC
# 第 72 章 —— %%
lasso_lars_ic.set_params(lassolarsic__criterion="bic").fit(X, y)
bic_criterion = zou_et_al_criterion_rescaling(
lasso_lars_ic[-1].criterion_,
n_samples,
lasso_lars_ic[-1].noise_variance_,
)
index_alpha_path_bic = np.flatnonzero(
lasso_lars_ic[-1].alphas_ == lasso_lars_ic[-1].alpha_
)[0]
| 行号 | 解释 |
|------|------|
| 1‑12 | 修改管道参数将准则切换为 BIC，重新拟合模型；同样提取 BIC 准则、重新缩放以及定位最优 BIC 对应的 alpha 在路径中的索引。 |
72.5.3.7 代码块 ⑦ — 验证 AIC 与 BIC 最优 alpha 是否一致
# 第 72 章 —— %%
# 第 72 章 —— Now that we collected the AIC and BIC, we can as well check that the minima
# 第 72 章 —— of both criteria happen at the same alpha. Then, we can simplify the
# 第 72 章 —— following plot.
index_alpha_path_aic == index_alpha_path_bic
| 行号 | 解释 |
|------|------|
| 1‑3 | 检查 AIC 和 BIC 选中的 alpha 是否对应路径中的同一索引；若为真，则两种准则在本数据集上选取了相同的正则化强度，可简化后续绘图（仅需绘制一条竖线）。 |
72.5.3.8 代码块 ⑧ — 绘制 AIC/BIC 曲线并标记选中点
# 第 72 章 —— %%
# 第 72 章 —— Finally, we can plot the AIC and BIC criterion and the subsequent selected
# 第 72 章 —— regularization parameter.
import matplotlib.pyplot as plt
plt.plot(aic_criterion, color="tab:blue", marker="o", label="AIC criterion")
plt.plot(bic_criterion, color="tab:orange", marker="o", label="BIC criterion")
plt.vlines(
index_alpha_path_bic,
aic_criterion.min(),
aic_criterion.max(),
color="black",
linestyle="--",
label="Selected alpha",
)
plt.legend()
plt.ylabel("Information criterion")
plt.xlabel("Lasso model sequence")
_ = plt.title("Lasso model selection via AIC and BIC")
| 行号 | 解释 |
|------|------|
| 1‑12 | 使用 matplotlib 绘制重新缩放后的 AIC（蓝色圆点）和 BIC（橙色圆点）随模型序列（即路径中 alpha 的递增顺序）的变化；添加黑色虚线垂直线标记选中的 alpha（即两个准则的最小点对应的位置）；添加图例、坐标轴标签和标题。 |
小结：本例演示了如何使用
LassoLarsIC基于 AIC 或 BIC 自动选择正则化强度alpha。通过重新缩放使其与理论文献保持一致，并展示了在糖尿病数据集上两种准则往往指向相同的最优点。该方法的优势在于仅需一次模型拟合（线性时间），适用于样本量较大且满足线性正态假设的场景；其局限在于对噪声方差的估计敏感，且在高维稀疏或强共线性情况下可能不如交叉验证鲁棒。
72.5.3.9 代码块 ⑨ — 单元流程图（Mermaid）
72.5.4 单元 4 — Lasso 模型选择（AIC/BIC vs 交叉验证）
源码路径：examples/linear_model/plot_lasso_model_selection.py - 完整文件
72.5.4.1 代码块 ① — 模块说明 & 导入
"""
=================================================
Lasso model selection: AIC-BIC / cross-validation
=================================================
...
"""
# 第 72 章 —— Authors: The scikit-learn developers
# 第 72 章 —— SPDX-License-Identifier: BSD-3-Clause
72.5.4.2 代码块 ② — 加载并查看数据
# 第 72 章 —— %%
# 第 72 章 —— Dataset
# 第 72 章 —— -------
# 第 72 章 —— In this example, we will use the diabetes dataset.
from sklearn.datasets import load_diabetes
X, y = load_diabetes(return_X_y=True, as_frame=True)
X.head()
| 行号 | 解释 |
|------|------|
| 1‑5 | 加载糖尿病数据集为 Pandas DataFrame，展示前几行以确认特征名称和数据类型。 |
72.5.4.3 代码块 ③ — 添加随机特征以增大特征数
# 第 72 章 —— %%
# 第 72 章 —— In addition, we add some random features to the original data to
# 第 72 章 —— better illustrate the feature selection performed by the Lasso model.
import numpy as np
import pandas as pd
rng = np.random.RandomState(42)
n_random_features = 14
X_random = pd.DataFrame(
rng.randn(X.shape[0], n_random_features),
columns=[f"random_{i:02d}" for i in range(n_random_features)],
)
X = pd.concat([X, X_random], axis=1)
# 第 72 章 —— Show only a subset of the columns
X[X.columns[::3]].head()
| 行号 | 解释 |
|------|------|
| 1‑9 | 生成 14 个服从标准正态分布的随机特征并与原始数据横向拼接，以增加特征维度并考察 Lasso 在高维噪声特征中的特征选择能力；最后展示每隔三列的子集以验证拼接成功。 |
72.5.4.4 代码块 ④ — 基于信息准则的模型选择（AIC）
# 第 72 章 —— %%
# 第 72 章 —— Selecting Lasso via an information criterion
# 第 72 章 —— --------------------------------------------
# 第 72 章 —— :class:`~sklearn.linear_model.LassoLarsIC` provides a Lasso estimator that
# 第 72 章 —— uses the Akaike information criterion (AIC) or the Bayes information
# 第 72 章 —— criterion (BIC) to select the optimal value of the regularization
# 第 72 章 —— parameter alpha.
#
# 第 72 章 —— Before fitting the model, we will standardize the data with a
# 第 72 章 —— :class:`~sklearn.preprocessing.StandardScaler`. In addition, we will
# 第 72 章 —— measure the time to fit and tune the hyperparameter alpha in order to
# 第 72 章 —— compare with the cross-validation strategy.
#
# 第 72 章 —— We will first fit a Lasso model with the AIC criterion.
import time
from sklearn.linear_model import LassoLarsIC
from sklearn.pipeline import make_pipeline
from sklearn.preprocessing import StandardScaler
start_time = time.time()
lasso_lars_ic = make_pipeline(StandardScaler(), LassoLarsIC(criterion="aic")).fit(X, y)
fit_time = time.time() - start_time
| 行号 | 解释 |
|------|------|
| 1‑16 | 记录开始时间；构建包含标准化和 LassoLarsIC(criterion="aic") 的管道；拟合模型并计算拟合耗时，为后续与交叉验证的时间开销比较提供基础。 |
72.5.4.5 代码块 ⑤ — 提取 AIC 曲线并识别最优 alpha
# 第 72 章 —— %%
# 第 72 章 —— We store the AIC metric for each value of alpha used during `fit`.
results = pd.DataFrame(
{
"alphas": lasso_lars_ic[-1].alphas_,
"AIC criterion": lasso_lars_ic[-1].criterion_,
}
).set_index("alphas")
alpha_aic = lasso_lars_ic[-1].alpha_
| 行号 | 解释 |
|------|------|
| 1‑6 | 将管道中模型的 alphas_ 路径和对应的 AIC 准则值存入 DataFrame，以 alpha 为索引；提取所选中的 alpha（即 AIC 最小点对应的值）保存至变量 alpha_aic。 |
72.5.4.6 代码块 ⑥ — 基于信息准则的模型选择（BIC）
# 第 72 章 —— %%
# 第 72 章 —— Now, we perform the same analysis using the BIC criterion.
lasso_lars_ic.set_params(lassolarsic__criterion="bic").fit(X, y)
results["BIC criterion"] = lasso_lars_ic[-1].criterion_
alpha_bic = lasso_lars_ic[-1].alpha_
| 行号 | 解释 |
|------|------|
| 1‑4 | 修改管道参数将准则切换为 BIC，重新拟合模型；将 BIC 准则值加入之前的 DataFrame；提取所选中的 alpha（BIC 最小点）保存至变量 alpha_bic。 |
72.5.4.7 代码块 ⑦ — 突出显示 DataFrame 中的最小值（AIC/BIC）
# 第 72 章 —— %%
# 第 72 章 —— We can check which value of `alpha` leads to the minimum AIC and BIC.
def highlight_min(x):
x_min = x.min()
return ["font-weight: bold" if v == x_min else "" for v in x]
results.style.apply(highlight_min)
| 行号 | 解释 |
|------|------|
| 1‑4 | 定义一个用于 Pandas Styler 的函数，将 DataFrame 中每列的最小值加粗显示；应用于包含 AIC 和 BIC 的结果表格，以直观标出每种准则下的最优 alpha。 |
72.5.4.8 代码块 ⑧ — 绘制 AIC/BIC 曲线并标记选中点
# 第 72 章 —— %%
# 第 72 章 —— Finally, we can plot the AIC and BIC values for the different alpha values.
# 第 72 章 —— The vertical lines in the plot correspond to the alpha chosen for each
# 第 72 章 —— criterion. The selected alpha corresponds to the minimum of the AIC or BIC
# 第 72 章 —— criterion.
ax = results.plot()
ax.vlines(
alpha_aic,
results["AIC criterion"].min(),
results["AIC criterion"].max(),
label="alpha: AIC estimate",
linestyles="--",
color="tab:blue",
)
ax.vlines(
alpha_bic,
results["BIC criterion"].min(),
results["BIC criterion"].max(),
label="alpha: BIC estimate",
linestyle="--",
color="tab:orange",
)
ax.set_xlabel(r"$\alpha$")
ax.set_ylabel("criterion")
ax.set_xscale("log")
ax.legend()
_ = ax.set_title(
f"Information-criterion for model selection (training time {fit_time:.2f}s)"
)
| 行号 | 解释 |
|------|------|
| 1‑16 | 使用 pandas 内置绘图绘制 AIC 和 BIC 随 alpha（对数尺度）的变化曲线；添加蓝色和橙色的虚线垂直线分别标记 AIC 和 BIC 选中的 alpha；设置坐标轴标签、图例以及标题，标题中包含基于信息准则的模型选择训练耗时。 |
72.5.4.9 代码块 ⑨ — 信息准则方法的讨论（优点与局限）
# 第 72 章 —— %%
# 第 72 章 —— Model selection with an information-criterion is very fast. It relies on
# 第 72 章 —— computing the criterion on the in-sample set provided to `fit`. Both criteria
# 第 72 章 —— estimate the model generalization error based on the training set error and
# 第 72 章 —— penalize this overly optimistic error. However, this penalty relies on a
# 第 72 章 —— proper estimation of the degrees of freedom and the noise variance. Both are
# 第 72 章 —— derived for large samples (asymptotic results) and assume the model is
# 第 72 章 —— correct, i.e. that the data are actually generated by this model.
#
# 第 72 章 —— These models also tend to break when the problem is badly conditioned (more
# 第 72 章 —— features than samples). It is then required to provide an estimate of the
# 第 72 章 —— noise variance.
| 行号 | 解释 |
|------|------|
| 1‑10 | 阐释信息准则（AIC/BIC）模型选择的优势：仅需一次拟合，计算快速；其理论基础是基于训练数据的似然，通过惩罚项校正过度乐观的训练误差以估计泛化性能。同时指出其局限：对自由度和噪声方差的估计依赖渐近假设（大样本），在特征多于样本（p > n）或模型误设时可能失效，此时需要外部噪声方差估计。 |
72.5.4.10 代码块 ⑩ — 基于交叉验证的模型选择（坐标下降 LassoCV）
# 第 72 章 —— %%
# 第 72 章 —— Selecting Lasso via cross-validation
# 第 72 章 —— ------------------------------------
# 第 72 章 —— The Lasso estimator can be implemented with different solvers: coordinate
# 第 72 章 —— descent and least angle regression. They differ with regards to their
# 第 72 章 —— execution speed and sources of numerical errors.
#
# 第 72 章 —— In scikit-learn, two different estimators are available with integrated
# 第 72 章 —— cross-validation: :class:`~sklearn.linear_model.LassoCV` and
# 第 72 章 —— :class:`~sklearn.linear_model.lassoLarsCV` that respectively solve the
# 第 72 章 —— problem with coordinate descent and least angle regression.
#
# 第 72 章 —— In the remainder of this section, we will present both approaches. For both
# 第 72 章 —— algorithms, we will use a 20-fold cross-validation strategy.
#
# 第 72 章 —— Lasso via coordinate descent
# 第 72 章 —— ............................
# 第 72 章 —— Let's start by making the hyperparameter tuning using
# 第 72 章 —— :class:`~sklearn.linear_model.LassoCV`.
from sklearn.linear_model import LassoCV
start_time = time.time()
model = make_pipeline(StandardScaler(), LassoCV(cv=20)).fit(X, y)
fit_time = time.time() - start_time
| 行号 | 解释 |
|------|------|
| 1‑16 | 描述交叉验证在 Lasso 中的两种求解器实现（坐标下降与 LARS）；记录开始时间；构建包含标准化和 LassoCV(cv=20) 的管道；拟合模型并计算训练耗时，为后续与信息准则方法比较提供基础。 |
72.5.4.11 代码块 ⑪ — 绘制坐标下降 LassoCV 的均方误差路径
# 第 72 章 —— %%
import matplotlib.pyplot as plt
ymin, ymax = 2300, 3800
lasso = model[-1]
plt.semilogx(lasso.alphas_, lasso.mse_path_, linestyle=":")
plt.plot(
lasso.alphas_,
lasso.mse_path_.mean(axis=-1),
color="black",
label="Average across the folds",
linewidth=2,
)
plt.axvline(lasso.alpha_, linestyle="--", color="black", label="alpha: CV estimate")
plt.ylim(ymin, ymax)
plt.xlabel(r"$\alpha$")
plt.ylabel("Mean square error")
plt.legend()
_ = plt.title(
f"Mean square error on each fold: coordinate descent (train time: {fit_time:.2f}s)"
)
| 行号 | 解释 |
|------|------|
| 1‑12 | 提取管道中训练好的 LassoCV 模型；使用 semilogx 绘制每折的 MSE 路径（虚线）；绘制所有折的平均 MSE（粗黑线）；添加垂直虚线标记通过交叉验证选择的最优 alpha；限制 y 轴范围以聚焦在相关区域；添加坐标轴标签、图例以及标题，标题中包含基于坐标下降的交叉验证训练耗时。 |
72.5.4.12 代码块 ⑫ — 基于交叉验证的模型选择（LARS LassoLarsCV）
# 第 72 章 —— %%
# 第 72 章 —— Lasso via least angle regression
# 第 72 章 —— ................................
# 第 72 章 —— Let's start by making the hyperparameter tuning using
# 第 72 章 —— :class:`~sklearn.linear_model.LassoLarsCV`.
from sklearn.linear_model import LassoLarsCV
start_time = time.time()
model = make_pipeline(StandardScaler(), LassoLarsCV(cv=20)).fit(X, y)
fit_time = time.time() - start_time
| 行号 | 解释 |
|------|------|
| 1‑8 | 记录开始时间；构建包含标准化和 LassoLarsCV(cv=20) 的管道；拟合模型并计算训练耗时。 |
72.5.4.13 代码块 ⑬ — 绘制 LARS LassoLarsCV 的均方误差路径
# 第 72 章 —— %%
lasso = model[-1]
plt.semilogx(lasso.cv_alphas_, lasso.mse_path_, ":")
plt.semilogx(
lasso.cv_alphas_,
lasso.mse_path_.mean(axis=-1),
color="black",
label="Average across the folds",
linewidth=2,
)
plt.axvline(lasso.alpha_, linestyle="--", color="black", label="alpha CV")
plt.ylim(ymin, ymax)
plt.xlabel(r"$\alpha$")
plt.ylabel("Mean square error")
plt.legend()
_ = plt.title(f"Mean square error on each fold: Lars (train time: {fit_time:.2f}s)")
| 行号 | 解释 |
|------|------|
| 1‑12 | 提取管道中训练好的 LassoLarsCV 模型；使用 semilogx 绘制每折的 MSE 路径（点线）；绘制所有折的平均 MSE（粗黑线）；添加垂直虚线标记通过交叉验证选择的最优 alpha；设置坐标轴范围、标签、图例以及标题，标题中包含基于 LARS 的交叉验证训练耗时。 |
72.5.4.14 代码块 ⑬ — 交叉验证方法的总结
# 第 72 章 —— %%
# 第 72 章 —— Summary of cross-validation approach
# 第 72 章 —— ....................................
# 第 72 章 —— Both algorithms give roughly the same results.
#
# 第 72 章 —— Lars computes a solution path only for each kink in the path. As a result, it
# 第 72 章 —— is very efficient when there are only of few kinks, which is the case if
# 第 72 章 —— there are few features or samples. Also, it is able to compute the full path
# 第 72 章 —— without setting any hyperparameter. On the opposite, coordinate descent
# 第 72 章 —— computes the path points on a pre-specified grid (here we use the default).
# 第 72 章 —— Thus it is more efficient if the number of grid points is smaller than the
# 第 72 章 —— number of kinks in the path. Such a strategy can be interesting if the number
# 第 72 章 —— of features is really large and there are enough samples to be selected in
# 第 72 章 —— each of the cross-validation fold. In terms of numerical errors, for heavily
# 第 72 章 —— correlated variables, Lars will accumulate more errors, while the coordinate
# 第 72 章 —— descent algorithm will only sample the path on a grid.
#
# 第 72 章 —— Note how the optimal value of alpha varies for each fold. This illustrates
# 第 72 章 —— why nested-cross validation is a good strategy when trying to evaluate the
# 第 72 章 —— performance of a method for which a parameter is chosen by cross-validation:
# 第 72 章 —— this choice of parameter may not be optimal for a final evaluation on
# 第 72 章 —— unseen test set only.
| 行号 | 解释 |
|------|------|
| 1‑14 | 比较 LARS 和坐标下降两种交叉验证实现的效率特点：LARS 在路径折点少时更快且可无超参数全路径计算；坐标下降在预指定网格点少时更优，适用于高维特征且每折有足够样本进行选择；讨论了数值误差来源：LARS 在高相关下误差累积，而坐标下降仅在网格点采样。强调了每折 alpha 的变化凸显了嵌套交叉验证的必要性，以避免在最终评估时使用偏频的超参数选择。 |
72.5.4.15 代码块 ⑭ — 结论：信息准则 vs 交叉验证
# 第 72 章 —— %%
# 第 72 章 —— Conclusion
# 第 72 章 —— ----------
# 第 72 章 —— In this tutorial, we presented two approaches for selecting the best
# 第 72 章 —— hyperparameter `alpha`: one strategy finds the optimal value of `alpha`
# 第 72 章 —— by only using the training set and some information criterion, and another
# 第 72 章 —— strategy is based on cross-validation.
#
# 第 72 章 —— In this example, both approaches are working similarly. The in-sample
# 第 72 章 —— hyperparameter selection even shows its efficacy in terms of computational
# 第 72 章 —— performance. However, it can only be used when the number of samples is large
# 第 72 章 —— enough compared to the number of features.
#
# 第 72 章 —— That's why hyperparameter optimization via cross-validation is a safe
# 第 72 章 —— strategy: it works in different settings.
| 行号 | 解释 |
|------|------|
| 1‑10 | 总结两种模型选择策略：信息准则（AIC/BIC）基于训练数据一次性计算，速度快但依赖模型假设和噪声方差估计；交叉验证（CV）通过多次拟合估计泛化性能，计算开销更大但更鲁棒，尤其在样本量不足、特征高维或模型假设不满足时更可靠。 concluding that CV 是更安全的通用选择。 |
72.5.4.16 代码块 ⑮ — 单元流程图（Mermaid）
72.6 设计中的取舍
在实际项目中，我们常常要在 完整路径可视化 与 仅报告最终模型 之间作出权衡。路径可视化 能够呈现系数随正则化强度的 分段线性结构，帮助我们辨识 特征进入/退出模型的顺序，从而判断模型的 特征选择稳定性 与 求解器差异（例如 Lasso 与 LARS 在小 alpha 区域的微小偏差）。然而，这种细粒度的可视化需要 多次路径求解（每个 alpha 都要计算一次），计算开销为 O(n·p·#α)，在 特征数上万、样本数上万 时会显著增加运行时间。
相对地，仅报告最终模型（如通过交叉验证得到的最佳 alpha）只需要一次模型拟合，代价约为 O(n·p·log p)，非常适合 生产环境的批量预测 或 快速原型迭代。但缺失的路径信息会让我们失去对 模型内部稀疏化过程 的洞见，尤其在高度共线或噪声较大的场景下，可能导致 不稳定的特征选择（相同数据在不同随机种子下选出不同特征）而难以解释。
因此，取舍原则如下：
-
教学、调参或科研阶段，推荐 完整路径可视化，帮助建立对正则化几何的直观认知。
-
大规模生产（如线上推荐、实时风险评估），应优先采用 交叉验证或信息准则 直接得到的单一
alpha，并配合 模型压缩（如特征哈希）以降低推理成本。
72.7 动手练习
练习的步骤如下：首先加载并标准化数据；其次计算 Lasso、ElasticNet 与 Lasso‑LARS 三种正则化路径；然后使用 LassoCV 通过交叉验证得到最佳正则化强度 alpha_cv；最后在对数坐标下绘制三条路径曲线，并用竖虚线标记 alpha_cv 的位置，以直观比较不同求解器的系数轨迹并观察交叉验证选择点在路径中的相对位置。完整代码如下所示。
练习的步骤如下：首先加载并标准化数据；其次计算 Lasso、ElasticNet 与 Lasso‑LARS 三种正则化路径；然后使用 LassoCV 通过交叉验证得到最佳正则化强度 alpha_cv；最后在对数坐标下绘制三条路径曲线，并用竖虚线标记 alpha_cv 的位置，以直观比较不同求解器的系数轨迹并观察交叉验证选择点在路径中的相对位置。完整代码如下所示。
# 第 72 章 —— -*- coding: utf-8 -*-
"""
Exercise 1 – 绘制 Lasso、ElasticNet、Lasso‑LARS 三条路径并标记 CV 选出的 α*
"""
import numpy as np
import matplotlib.pyplot as plt
from sklearn.datasets import load_diabetes
from sklearn.linear_model import lasso_path, lars_path, enet_path, LassoCV
# 第 72 章 —— 1️⃣ 加载并标准化
X, y = load_diabetes(return_X_y=True)
X = X / X.std(axis=0)
# 第 72 章 —— 2️⃣ 计算路径（与单元1相同）
eps = 5e-3
alphas_lasso, coefs_lasso, _ = lasso_path(X, y, eps=eps)
alphas_lars, _, coefs_lars = lars_path(X, y, method="lasso")
alphas_enet, coefs_enet, _ = enet_path(X, y, eps=eps, l1_ratio=0.8)
# 第 72 章 —— 3️⃣ 交叉验证得到最佳 α（使用 LassoCV 为例）
lasso_cv = LassoCV(cv=5, n_alphas=200).fit(X, y)
alpha_cv = lasso_cv.alpha_
# 第 72 章 —— 4️⃣ 绘图
plt.figure(figsize=(8, 6))
colors = plt.cm.tab10(np.linspace(0, 1, 3))
# 第 72 章 —— Lasso
plt.semilogx(alphas_lasso, coefs_lasso.T, color=colors[0], lw=1, label="Lasso")
# 第 72 章 —— ElasticNet
plt.semilogx(alphas_enet, coefs_enet.T, '--', color=colors[1], lw=1, label="ElasticNet")
# 第 72 章 —— LARS
plt.semilogx(alphas_lars, coefs_lars.T, ':', color=colors[2], lw=1, label="LARS")
# 第 72 章 —— 标记 CV 选出的 α*
plt.axvline(alpha_cv, color='k', linestyle='-.', linewidth=2,
label=f'CV α* = {alpha_cv:.3e}')
plt.xlabel('alpha (log scale)')
plt.ylabel('coefficients')
plt.title('Lasso / ElasticNet / LARS 正则化路径\n并标记 CV 选出的 α*')
plt.legend(loc='lower left')
plt.tight_layout()
plt.show()
练习的步骤如下：首先参考单元 4 的代码框架；其次在外层循环中加入噪声强度参数（例如从 0.1 到 2.0 的等间距值）；第三步对每个噪声水平重复实验 50 次，每次使用不同随机种子重新生成数据并训练 LassoCV 与 LassoLarsIC 模型；第四步记录每个模型在每次实验中选择的特征索引；最后统计所有实验中每个特征被选中的频率，绘制条形图或热图以比较两种模型在不同噪声下的特征选择稳定性。完整代码框架可参见单元 4，仅需在数据生成部分引入噪声参数并增加重复循环。
72.8 单元 5‑24（简要概览）
对于剩余 单元 5‑24（涉及 RANSAC、Theil‑Sen、Huber、GLM、贝叶斯模型、SGD 系列、Logistic、稀疏 Logistic、Polynomial、Multi‑Task Lasso、PCA/FA/ICA、核 PCA、增量 PCA、流形学习），本章采用 统一模板：
-
文件列表（已在 71.3 表中给出）。
-
核心代码块（从
# %%注释开始的每个独立块）。 -
逐行解释（采用表格形式，尽可能涵盖变量意义、函数作用、关键超参数的几何/统计解释）。
-
小结：该单元的几何意义、数值取舍、适用场景。
-
流程图（Mermaid）展示 数据→模型→评估→可视化 的完整路径。
下面给出 单元 12（Huber vs Ridge） 的示例结构，后续各单元可按此模板自行展开。
72.8.1 单元 12 — HuberRegressor vs Ridge（离群点鲁棒性）
| 步骤 | 代码 | 解释 |
|------|------|------|
| 1️⃣ | rng = np.random.RandomState(0) | 固定随机种子。 |
| 2️⃣ | X, y = make_regression(..., noise=4.0, bias=100.0) | 生成 线性 + 高噪声 的回归数据。 |
| 3️⃣ | 手动添加四个强离群点（X_outliers、y_outliers） | 将离群点 放在特征空间的极端，以检验模型鲁棒性。 |
| 4️⃣ | huber = HuberRegressor(alpha=0.0, epsilon=ε) | alpha=0.0 关闭 L2 正则，epsilon 控制 二次→线性 的转折点；epsilon 越大，模型越接近 Ridge。 |
| 5️⃣ | ridge = Ridge(alpha=0.0) | 纯 最小二乘（无正则），对离群点极为敏感。 |
| 6️⃣ | plt.plot(x, coef_, ...) | 在相同的 x 轴范围 上绘制两条回归直线，便于肉眼比较斜率、截距的差异。 |
| 7️⃣ | plt.legend, plt.title | 添加图例与标题，说明实验意图。 |
小结：HuberRegressor 在 |残差| < ε 区间使用 二次损失，在 |残差| ≥ ε 区间使用 线性损失，从而对离群点的影响仅为 线性（而不是二次），显著降低了异常值对模型的拉伸作用。
流程图
其余单元（如 Theil‑Sen、GLM（Poisson、Gamma、Tweedie）、SGD 系列、Logistic、稀疏 Logistic、Polynomial / Spline、Multi‑Task Lasso、PCA / FA / ICA、核/增量 PCA、流形学习）均遵循上述 “代码块 → 逐行解释 → 小结 → 流程图” 的结构，确保读者可以系统、统一地掌握每个示例的 数学背景、实现细节 与 可视化要点。
72.9 本章小结（表格前引导）
在本章中，我们系统地梳理了 线性模型族（稀疏、鲁棒、贝叶斯、SGD）以及 降维/流形学习（PCA、FA、ICA、核/增量 PCA、LLE、Isomap、MDS、Spectral、t‑SNE）对应的 源码实现、几何直觉 与 工程取舍。下面的对照表帮助快速回顾每类技术的核心属性与适用场景。
| 方法族 | 关键特性 | 典型超参数 | 适用场景 |
|--------|----------|-------------|----------|
| Lasso / ElasticNet / Ridge | 稀疏 vs 稳定 | alpha, l1_ratio | 高维特征选择、共线性抑制 |
| Huber / RANSAC / Theil‑Sen | 鲁棒性（对离群点） | epsilon, residual_threshold | 噪声异常数据的回归 |
| GLM (Poisson / Gamma / Tweedie) | 非正态、计数/保险数据 | link, power (Tweedie) | 频率、严重度、复合分布建模 |
| ARD / Bayesian Ridge | 贝叶斯稀疏、模型不确定性 | alpha_1, lambda_1 | 自动相关性判定、后验解释 |
| SGDClassifier / SGDRegressor | 大规模在线学习 | loss, penalty, learning_rate | 实时推荐、流式数据 |
| Logistic (L1/L2/ElasticNet) | 分类稀疏、概率估计 | C, l1_ratio, solver | 二/多分类、文本稀疏特征 |
| Multi‑Task Lasso (L₂₁) | 任务间稀疏共享 | alpha | 多任务回归、时间序列共享特征 |
| PCA / Probabilistic PCA | 方差最大化、概率解释 | n_components, svd_solver | 可视化、降噪、特征压缩 |
| Factor Analysis / Varimax | 异方差建模、旋转解释 | n_components, rotation | 潜在因子解释、心理测量 |
| ICA / FastICA | 非高斯独立分量 | n_components, whiten | 盲源分离、信号处理 |
| NMF / SparsePCA | 非负/稀疏分解 | n_components, alpha | 图像分解、主题建模 |
| Kernel PCA | 非线性映射 | kernel, gamma | 径向基分离、非线性特征抽取 |
| Incremental PCA | 大数据增量分解 | batch_size | 内存受限的大规模数据 |
| LLE / Isomap / LTSA / Hessian LLE | 保持局部/全局几何 | n_neighbors, n_components | 流形展开、非线性结构发现 |
| MDS / Classical MDS | 距离保持（度量/非度量） | metric, init | 低维距离可视化 |
| Spectral Embedding | 拉普拉斯特征映射 | n_neighbors, affinity | 社区检测、图嵌入 |
| t‑SNE | 局部聚类保留 | perplexity, learning_rate | 高维聚类可视化、数据探索 |
阅读建议：先从 稀疏正则化（单元 1‑4）入手，掌握 系数路径 与 模型选择；再逐步深入 鲁棒回归、贝叶斯模型 与 SGD；最后学习 降维/流形学习（单元 35‑51），形成从 线性到非线性、全局到局部 的完整知识链。
End of Chapter 71.
72.10 源码地图（按单元拆分）
examples/linear_model/plot_lasso_and_elasticnet.py
├── __main__ # Lasso vs ElasticNet 稀疏模式对比
examples/linear_model/plot_lasso_dense_vs_sparse_data.py
├── __main__ # 稀疏/稠密数据一致性验证
examples/linear_model/plot_lasso_lars_ic.py
├── __main__ # AIC/BIC 自动选择 alpha
examples/linear_model/plot_lasso_lasso_lars_elasticnet_path.py
├── __main__ # 正则化路径可视化对比
examples/linear_model/plot_lasso_model_selection.py
├── __main__ # CV 与信息准则模型选择对比
examples/linear_model/plot_elastic_net_precomputed_gram_matrix_with_weighted_samples.py
├── __main__ # 预计算 Gram 矩阵加速与加权样本
examples/linear_model/plot_ols_ridge.py
├── __main__ # OLS vs Ridge 共线性数据对比
examples/linear_model/plot_nnls.py
├── __main__ # 非负最小二乘物理约束
examples/linear_model/plot_omp.py
├── __main__ # OMP 贪心稀疏编码演示
examples/linear_model/plot_ridge_coeffs.py
├── __main__ # 岭回归系数随正则化收缩轨迹
examples/linear_model/plot_ridge_path.py
├── __main__ # 岭回归正则化路径可视化
examples/linear_model/plot_huber_vs_ridge.py
├── __main__ # Huber vs Ridge 离群点鲁棒性对比
examples/linear_model/plot_ransac.py
├── __main__ # RANSAC 随机共识剔除离群点
examples/linear_model/plot_theilsen.py
├── __main__ # Theil-Sen 中位数斜率高崩溃点回归
examples/linear_model/plot_robust_fit.py
├── __main__ # 多鲁棒估计器 MSE 基准测试
examples/linear_model/plot_quantile_regression.py
├── __main__ # Pinball 损失分位数回归建模尾部风险
examples/linear_model/plot_poisson_regression_non_normal_loss.py
├── __main__ # Poisson 回归对数链接建模计数数据
examples/linear_model/plot_tweedie_regression_insurance_claims.py
├── __main__ # Tweedie 复合分布保险理赔建模
examples/linear_model/plot_ard.py
├── __main__ # ARD 自动相关性确定稀疏贝叶斯
examples/linear_model/plot_bayesian_ridge_curvefit.py
├── __main__ # 贝叶斯岭回归证据最大化与预测区间
examples/linear_model/plot_sgd_early_stopping.py
├── __main__ # SGD 验证集早停防过拟合
examples/linear_model/plot_sgd_iris.py
├── __main__ # 多类 SGD 分类器收敛演示
examples/linear_model/plot_sgd_loss_functions.py
├── __main__ # 损失函数景观几何对比
examples/linear_model/plot_sgd_penalties.py
├── __main__ # L1/L2/ElasticNet 正则化几何
examples/linear_model/plot_sgd_separating_hyperplane.py
├── __main__ # 感知机更新决策边界演化可视化
examples/linear_model/plot_sgd_weighted_samples.py
├── __main__ # 样本权重拉拽决策边界效应
examples/linear_model/plot_sgdocsvm_vs_ocsvm.py
├── __main__ # SGDOneClassSVM vs LibSVM 流式异常检测
examples/linear_model/plot_logistic_l1_l2_sparsity.py
├── __main__ # L1 稀疏 vs L2 稠密系数对比
examples/linear_model/plot_logistic_multinomial.py
├── __main__ # 多项式 vs OvR 多分类决策边界
examples/linear_model/plot_logistic_path.py
├── __main__ # 逻辑回归正则化路径追踪
examples/linear_model/plot_sparse_logistic_regression_20newsgroups.py
├── __main__ # 20 Newsgroups 稀疏文本分类实战
examples/linear_model/plot_sparse_logistic_regression_mnist.py
├── __main__ # MNIST 像素稀疏逻辑回归基准
examples/linear_model/plot_polynomial_interpolation.py
├── __main__ # 多项式特征高次展开与 Runge 现象
examples/linear_model/plot_multi_task_lasso_support.py
├── __main__ # 多任务 Lasso 联合特征选择 L21 范数
examples/decomposition/plot_pca_iris.py
├── __main__ # Iris 数据集上 PCA 的 3D 可视化
examples/decomposition/plot_pca_vs_fa_model_selection.py
├── __main__ # Probabilistic PCA 与 Factor Analysis 的交叉验证模型选择对比
examples/decomposition/plot_pca_vs_lda.py
├── __main__ # LDA 与 PCA 在 Iris 数据集上的 2D 投影对比
examples/decomposition/plot_faces_decomposition.py
├── __main__ # Olivetti 人脸数据集上的分解方法可视化
├── plot_gallery # 组件可视化通用网格布局
examples/decomposition/plot_varimax_fa.py
├── __main__ # 带 Varimax 旋转的 Factor Analysis 可视化
examples/decomposition/plot_kernel_pca.py
├── __main__ # Kernel PCA 在 make_circles 数据上的非线性分离
examples/decomposition/plot_incremental_pca.py
├── __main__ # Incremental PCA 在 Iris 数据集上的近似验证
examples/decomposition/plot_ica_blind_source_separation.py
├── __main__ # FastICA 盲源分离与噪声鲁棒性验证
examples/decomposition/plot_ica_vs_pca.py
├── __main__ # ICA 与 PCA 在 2D 点云上的几何对比
├── plot_samples # 几何 ICA 可视化核心：源分布散点+投影方向箭头
examples/decomposition/plot_sparse_coding.py
├── __main__ # SparseCoder 在 Ricker 小波字典上的稀疏编码比较
├── ricker_function # 离散 Ricker 小波生成
├── ricker_matrix # 多宽度字典矩阵构造
examples/decomposition/plot_image_denoising.py
├── __main__ # 在线字典学习图像去噪完整流程
├── show_with_diff # 重构图像与差值图并排显示，Frobenius 范数量化误差
examples/manifold/plot_compare_methods.py
├── __main__ # 各流形学习方法在 S-curve 上的嵌入对比
├── plot_3d # 3D 散点图渲染管线
├── plot_2d # 2D 散点图渲染管线
├── add_2d_scatter # 2D 散点绘制原语
examples/manifold/plot_lle_digits.py
├── __main__ # 手写数字数据上各嵌入方法投射与计时对比
├── plot_embedding # 嵌入可视化全流程：归一化+标记+缩略图
examples/manifold/plot_manifold_sphere.py
├── __main__ # 各流形学习方法在 severed sphere 数据上的计时与嵌入对比
examples/manifold/plot_mds.py
├── __main__ # metric、non-metric、classical MDS 在噪声数据上的投影对比
examples/manifold/plot_swissroll.py
├── __main__ # LLE 与 t-SNE 在 Swiss Roll 与 Swiss-Hole 数据上的嵌入行为
examples/manifold/plot_t_sne_perplexity.py
├── __main__ # 不同 perplexity 下 t-SNE 在 circles、S-curve、uniform grid 上的聚类行为
第 73 章 —— 降维与流形学习 —— 探索高维数据的“几何密码”
73.1 学习目标
-
难度：★★★☆☆（3/5）
-
预备知识：Python 基础、面向对象编程与 Markdown/代码阅读基础
理解PCA的正交投影原理及其在方差最大化中的应用。PCA通过特征分解或奇异值分解（SVD）寻找方差最大的正交投影方向，投影后各特征互不相关。这一过程就像将多维数据投射到方差最大的轴上，从而捕获数据的主要变化方向。
掌握因子分析与独立成分分析在潜变量建模中的区别与联系。因子分析假设观测数据是潜在因子的线性组合加上噪声，侧重解释潜在协方差结构；独立成分分析则通过最大化非高斯性来估计统计独立的源信号，关注更高阶统计的独立性，二者在潜变量建模中各有侧重。
理解核PCA如何通过核技巧实现非线性降维。核PCA利用核函数将数据映射到高维特征空间，在该空间中执行线性PCA，从而实现非线性降维，能够处理如同心圆等线性方法无法分离的数据结构。
掌握t‑SNE在高维数据可视化中的拓扑展开机制。t‑SNE将高维相似度转化为低维条件概率，通过最小化KL散度来保持局部邻居关系，perplexity参数控制局部与全局结构的平衡，使得相近的点在低维空间仍保持相近。
理解流形学习算法（如LLE、Isomap）如何保持局部几何结构。LLE通过在每个点的局部邻域内进行线性重构来保持局部线性关系；Isomap则构建K‑近邻图并计算测地线距离，再用MDS保持这些全局距离，适合展开有洞的流形如瑞士卷。
73.2 生活类比
想象降维技术是一位“维度雕塑师”。PCA就像一台“正方体投影机”，它把多维数据投射到方差最大的正交轴上，好比用光从正方体的最宽面投出最长的影子，让我们一眼看出数据的主要方向。因子分析是“隐藏结构探测器”，它不只看表面的波动，而是聆听数据背后的“心跳”，寻找潜在因子，就像医生通过脉搏判断健康。独立成分分析像“信号源分离师”，把混合的音轨拆成独立的人声和乐器，就像调音台把每根音线单独抽出来。核PCA是“弯曲空间投影仪”，先用核函数把数据扭曲成可以直线分割的形状，再做线性投影，类似用透镜看扭曲的画作。t‑SNE则是“邻居保持制图师”，它尽量让距离近的点在低维空间仍保持相近，就像在社交网络图上把朋友们聚在一起。LLE是“局部几何保康师”，用小三角形贴合球面，保持局部平面结构；Isomap则是“测地线保持者”，沿着流形的最短路径（测地线）展开数据，像蚂蚁在苹果表面爬行而不走捷径。MDS是“距离保持制图师”，它根据高维数据之间的距离来重建低维空间，保持全局几何关系，就像根据城市间的实际距离绘制地图。SpectralEmbedding像“振动模式分析师”，通过数据构建的图的振动特征来寻找嵌入，捕捉数据的全局结构，类似分析鼓面的振动模式来理解其形状。
73.3 主成分与因子分析 —— 方差最大化的“正交投影仪”
核心概念
PCA 通过特征分解或奇异值分解（SVD）寻找方差最大的正交投影方向，投影后各特征互不相关。IncrementalPCA 采用增量学习，适用于无法一次性装入内存的大规模数据。KernelPCA 利用核技巧将数据映射到高维特征空间，在该空间执行线性 PCA，从而实现非线性降维。因子分析（FA）则假设观测数据是潜在因子线性组合加噪声的生成过程，侧重解释潜在协方差结构，并可通过 Varimax 旋转提升可解释性。
73.3.1 架构图
73.3.2 核心类型定义
源码路径：examples/decomposition/plot_pca_iris.py - PCA（第 30-110 行）
# 第 73 章 —— ① 导入必要的模块
from sklearn.datasets import load_iris
# 第 73 章 —— ② 导入 PCA 类
from sklearn.decomposition import PCA
# 第 73 章 —— ③ 用于 3D 可视化的导入（尽管未直接使用，但为 3D 投影所必需）
import mpl_toolkits.mplot3d # noqa: F401
# 第 73 章 —— ④ 导入绘图库
import matplotlib.pyplot as plt
# 第 73 章 —— ⑤ 加载鸢尾花数据集
iris = load_iris(as_frame=True)
# 第 73 章 —— ⑥ 打印数据集的键以查看其结构
print(iris.keys())
# 第 73 章 —— ⑦ 使用 seaborn 绘制特征两两配对图，按类别着色
import seaborn as sns
# 第 73 章 —— ⑧ 将目标值映射为类别名称，便于图例显示
iris.frame["target"] = iris.target_names[iris.target]
# 第 73 章 —— ⑨ 绘制成对特征散点图矩阵
_ = sns.pairplot(iris.frame, hue="target")
# 第 73 章 —— ⑩ 实例化 PCA，保留 3 个主成分
pca = PCA(n_components=3)
# 第 73 章 —— ⑪ 拟合并转换原始特征矩阵，将 4 维数据降到 3 维
X_reduced = pca.fit_transform(iris.data) # 输出形状为 (150, 3)
# 第 73 章 —— ⑫ 创建图形对象，准备进行 3D 散点图绘制
fig = plt.figure(1, figsize=(8, 6))
# 第 73 章 —— ⑫ 添加 3D 子图，设定俯视角度
ax = fig.add_subplot(111, projection="3d", elev=-150, azim=110)
# 第 73 章 —— ⑬ 在 3D 空间中绘制降维后的数据点，颜色对应鸢尾花类别
scatter = ax.scatter(
X_reduced[:, 0],
X_reduced[:, 1],
X_reduced[:, 2],
c=iris.target,
s=40,
)
# 第 73 章 —— ⑭ 设置坐标轴标签和标题
ax.set(
title="First three principal components",
xlabel="1st Principal Component",
ylabel="2nd Principal Component",
zlabel="3rd Principal Component",
)
# 第 73 章 —— ⑮ 隐藏坐标轴刻度标签，使图像更清晰
ax.xaxis.set_ticklabels([])
ax.yaxis.set_ticklabels([])
ax.zaxis.set_ticklabels([])
# 第 73 章 —— ⑯ 创建图例，将颜色映射到类别名称
legend1 = ax.legend(
scatter.legend_elements()[0],
iris.target_names.tolist(),
loc="upper right",
title="Classes",
)
# 第 73 章 —— ⑰ 将图例添加到图形中
ax.add_artist(legend1)
# 第 73 章 —— ⑱ 显示图形
plt.show()
# 第 73 章 —— ⑲ 打印每个主成分解释的方差比例，评估信息保留程度
print(pca.explained_variance_ratio_)
这段代码展示了如何使用 PCA 将四维鸢尾花特征压缩到三维，并输出每个主成分捕获的方差比例。
源码路径：examples/decomposition/plot_incremental_pca.py - IncrementalPCA（第 31-100 行）
# 第 73 章 —— ① 导入必要的模块
import matplotlib.pyplot as plt
import numpy as np
from sklearn.datasets import load_iris
from sklearn.decomposition import PCA, IncrementalPCA
# 第 73 章 —— ② 加载鸢尾花数据集
iris = load_iris()
X = iris.data
y = iris.target
# 第 73 章 —— ③ 设置目标降维维度
n_components = 2
# 第 73 章 —— ④ 创建 IncrementalPCA 实例，设定目标维度和批大小
ipca = IncrementalPCA(n_components=n_components, batch_size=10)
# 第 73 章 —— ⑤ 用整个数据集逐批拟合并转换，增量式地学习主成分
X_ipca = ipca.fit_transform(X)
# 第 73 章 —— ⑥ 创建标准 PCA 实例用于比较
pca = PCA(n_components=n_components)
# 第 73 章 —— ⑦ 拟合并转换数据，得到基准投影结果
X_pca = pca.fit_transform(X)
# 第 73 章 —— ⑧ 定义颜色列表，用于按类别着色散点图
colors = ["navy", "turquoise", "darkorange"]
# 第 73 章 —— ⑨ 遍历两种方法的结果，分别绘制散点图
for X_transformed, title in [(X_ipca, "Incremental PCA"), (X_pca, "PCA")]:
# ⑩ 创建新图形，设定大小
plt.figure(figsize=(8, 8))
# ⑪ 遍历每个类别，绘制散点图
for color, i, target_name in zip(colors, [0, 1, 2], iris.target_names):
plt.scatter(
X_transformed[y == i, 0],
X_transformed[y == i, 1],
color=color,
lw=2,
label=target_name,
)
# ⑫ 如果是 Incremental PCA，则计算并显示与标准 PCA 的平均绝对误差
if "Incremental" in title:
err = np.abs(np.abs(X_pca) - np.abs(X_ipca)).mean()
plt.title(title + " of iris dataset\nMean absolute unsigned error %.6f" % err)
else:
# ⑬ 否则，仅显示方法名称作为标题
plt.title(title + " of iris dataset")
# ⑭ 显示图例
plt.legend(loc="best", shadow=False, scatterpoints=1)
# ⑮ 设置坐标轴范围以确保一致比较
plt.axis([-4, 4, -1.5, 1.5])
# 第 73 章 —— ⑯ 显示所有图形
plt.show()
这段代码演示了在内存受限场景下，IncrementalPCA 如何以小批量方式逼近传统 PCA 的投影结果。
源码路径：examples/decomposition/plot_kernel_pca.py - KernelPCA（第 46-96 行）
# 第 73 章 —— ① 导入必要的模块
from sklearn.datasets import make_circles
from sklearn.model_selection import train_test_split
# 第 73 章 —— ② 生成两个同心圆数据集，用于演示非线性降维
X, y = make_circles(n_samples=1_000, factor=0.3, noise=0.05, random_state=0)
# 第 73 章 —— ③ 划分训练集和测试集，保持类别比例
X_train, X_test, y_train, y_test = train_test_split(X, y, stratify=y, random_state=0)
# 第 73 章 —— ④ 可视化训练集和测试集的原始分布
import matplotlib.pyplot as plt
# 第 73 章 —— ⑤ 创建包含两个子图的图形，共享 x 和 y 轴
_, (train_ax, test_ax) = plt.subplots(ncols=2, sharex=True, sharey=True, figsize=(8, 4))
# 第 73 章 —— ⑥ 绘制训练集散点图，颜色对应类别
train_ax.scatter(X_train[:, 0], X_train[:, 1], c=y_train)
# 第 73 章 —— ⑦ 设置训练集图的标签和标题
train_ax.set_ylabel("Feature #1")
train_ax.set_xlabel("Feature #0")
train_ax.set_title("Training data")
# 第 73 章 —— ⑧ 绘制测试集散点图
test_ax.scatter(X_test[:, 0], X_test[:, 1], c=y_test)
# 第 73 章 —— ⑨ 设置测试集图的标签和标题
test_ax.set_xlabel("Feature #0")
# 第 73 章 —— ⑩ 设置测试集图的标题
_ = test_ax.set_title("Testing data")
# 第 73 章 —— ⑪ 导入 PCA 和 KernelPCA 类
from sklearn.decomposition import PCA, KernelPCA
# 第 73 章 —— ⑫ 实例化标准 PCA，保留 2 个主成分
pca = PCA(n_components=2)
# 第 73 章 —— ⑬ 实例化 KernelPCA，使用 RBF 核并启用逆变换功能
kernel_pca = KernelPCA(
n_components=None, # 保留所有特征（等价于原始维度）
kernel="rbf",
gamma=10,
fit_inverse_transform=True,
alpha=0.1,
)
# 第 73 章 —— ⑭ 在训练集上拟合 PCA 模型
pca.fit(X_train)
# 第 73 章 —— ⑮ 在测试集上应用 PCA 投影
X_test_pca = pca.transform(X_test)
# 第 73 章 —— ⑯ 在训练集上拟合 KernelPCA 模型
kernel_pca.fit(X_train)
# 第 73 章 —— ⑰ 在测试集上应用 KernelPCA 投影
X_test_kernel_pca = kernel_pca.transform(X_test)
# 第 73 章 —— ⑱ 创建包含三个子图的图形，用于比较原始数据、PCA 投影和 KernelPCA 投影
fig, (orig_data_ax, pca_proj_ax, kernel_pca_proj_ax) = plt.subplots(
ncols=3, figsize=(14, 4)
)
# 第 73 章 —— ⑲ 在第一个子图中绘制原始测试数据
orig_data_ax.scatter(X_test[:, 0], X_test[:, 1], c=y_test)
orig_data_ax.set_ylabel("Feature #1")
orig_data_ax.set_xlabel("Feature #0")
orig_data_ax.set_title("Testing data")
# 第 73 章 —— ⑳ 在第二个子图中绘制 PCA 投影结果
pca_proj_ax.scatter(X_test_pca[:, 0], X_test_pca[:, 1], c=y_test)
pca_proj_ax.set_ylabel("Principal component #1")
pca_proj_ax.set_xlabel("Principal component #0")
pca_proj_ax.set_title("Projection of testing data\n using PCA")
# 第 73 章 —— ㉑ 在第三个子图中绘制 KernelPCA 投影结果
kernel_pca_proj_ax.scatter(X_test_kernel_pca[:, 0], X_test_kernel_pca[:, 1], c=y_test)
kernel_pca_proj_ax.set_ylabel("Principal component #1")
kernel_pca_proj_ax.set_xlabel("Principal component #0")
# 第 73 章 —— ㉒ 设置第三个子图的标题
_ = kernel_pca_proj_ax.set_title("Projection of testing data\n using KernelPCA")
# 第 73 章 —— ㉓ 使用 PCA 的逆变换将投影数据近似还原到原始特征空间
X_reconstructed_pca = pca.inverse_transform(pca.transform(X_test))
# 第 73 章 —— ㉔ 使用 KernelPCA 的逆变换将投影数据近似还原到原始特征空间（结果为近似值）
X_reconstructed_kernel_pca = kernel_pca.inverse_transform(kernel_pca.transform(X_test))
# 第 73 章 —— ㉕ 创建包含三个子图的图形，用于比较原始数据、PCA 重构和 KernelPCA 重构
fig, (orig_data_ax, pca_back_proj_ax, kernel_pca_back_proj_ax) = plt.subplots(
ncols=3, sharex=True, sharey=True, figsize=(13, 4)
)
# 第 73 章 —— ㉖ 在第一个子图中绘制原始测试数据
orig_data_ax.scatter(X_test[:, 0], X_test[:, 1], c=y_test)
orig_data_ax.set_ylabel("Feature #1")
orig_data_ax.set_xlabel("Feature #0")
orig_data_ax.set_title("Original test data")
# 第 73 章 —— ㉗ 在第二个子图中绘制 PCA 重构结果
pca_back_proj_ax.scatter(X_reconstructed_pca[:, 0], X_reconstructed_pca[:, 1], c=y_test)
pca_back_proj_ax.set_xlabel("Feature #0")
pca_back_proj_ax.set_title("Reconstruction via PCA")
# 第 73 章 —— ㉘ 在第三个子图中绘制 KernelPCA 重构结果
kernel_pca_back_proj_ax.scatter(
X_reconstructed_kernel_pca[:, 0], X_reconstructed_kernel_pca[:, 1], c=y_test
)
kernel_pca_back_proj_ax.set_xlabel("Feature #0")
# 第 73 章 —— ㉙ 设置第三个子图的标题
_ = kernel_pca_back_proj_ax.set_title("Reconstruction via KernelPCA")
这段代码说明了通过 RBF 核将同心圆数据映射到可线性分割的空间，并展示了逆变换的近似性质。
源码路径：examples/decomposition/plot_pca_vs_fa_model_selection.py - FactorAnalysis（第 63-120 行）
# 第 73 章 —— ① 导入必要的模块
import numpy as np
from scipy import linalg
# 第 73 章 —— ② 生成低秩矩阵作为纯净数据，再添加噪声模拟真实场景
n_samples, n_features, rank = 500, 25, 5
sigma = 1.0
rng = np.random.RandomState(42)
U, _, _ = linalg.svd(rng.randn(n_features, n_features))
X = np.dot(rng.randn(n_samples, rank), U[:, :rank].T)
# 第 73 章 —— ③ 添加同方差噪声（所有特征噪声方差相同）
X_homo = X + sigma * rng.randn(n_samples, n_features)
# 第 73 章 —— ④ 添加异方差噪声（每个特征噪声方差不同）
sigmas = sigma * rng.rand(n_features) + sigma / 2.0
X_hetero = X + rng.randn(n_samples, n_features) * sigmas
# 第 73 章 —— ⑤ 导入用于模型评估和比较的模块
import matplotlib.pyplot as plt
from sklearn.covariance import LedoitWolf, ShrunkCovariance
from sklearn.decomposition import PCA, FactorAnalysis
from sklearn.model_selection import GridSearchCV, cross_val_score
# 第 73 章 —— ⑥ 定义要评估的成分数范围
n_components = np.arange(0, n_features, 5) # options for n_components
# 第 73 章 —— ⑦ 定义一个函数，用于通过交叉验证计算 PCA 和 FA 在给定数据上的平均得分
def compute_scores(X):
pca = PCA(svd_solver="full")
fa = FactorAnalysis()
pca_scores, fa_scores = [], []
for n in n_components:
pca.n_components = n
fa.n_components = n
pca_scores.append(np.mean(cross_val_score(pca, X)))
fa_scores.append(np.mean(cross_val_score(fa, X)))
return pca_scores, fa_scores
# 第 73 章 —— ⑧ 定义一个函数，用于评估 shrinkage 协方差估计器的性能
def shrunk_cov_score(X):
shrinkages = np.logspace(-2, 0, 30)
cv = GridSearchCV(ShrunkCovariance(), {"shrinkage": shrinkages})
return np.mean(cross_val_score(cv.fit(X).best_estimator_, X))
# 第 73 章 —— ⑨ 定义一个函数，用于评估 LedoitWolf 协方差估计器的性能
def lw_score(X):
return np.mean(cross_val_score(LedoitWolf(), X))
# 第 73 章 —— ⑩ 分别对同方差噪声和异方差噪声数据进行模型比较
for X, title in [(X_homo, "Homoscedastic Noise"), (X_hetero, "Heteroscedastic Noise")]:
pca_scores, fa_scores = compute_scores(X)
n_components_pca = n_components[np.argmax(pca_scores)]
n_components_fa = n_components[np.argmax(fa_scores)]
# ⑪ 使用 MLE 方法自动选择 PCA 的成分数
pca = PCA(svd_solver="full", n_components="mle")
pca.fit(X)
n_components_pca_mle = pca.n_components_
# ⑫ 打印不同方法选择的最佳成分数
print("best n_components by PCA CV = %d" % n_components_pca)
print("best n_components by FactorAnalysis CV = %d" % n_components_fa)
print("best n_components by PCA MLE = %d" % n_components_pca_mle)
# ⑬ 创建图形，绘制 PCA 和 FA 的交叉验证得分曲线
plt.figure()
plt.plot(n_components, pca_scores, "b", label="PCA scores")
plt.plot(n_components, fa_scores, "r", label="FA scores")
# ⑭ 添加垂直线表示真实秩和各方法选择的成分数
plt.axvline(rank, color="g", label="TRUTH: %d" % rank, linestyle="-")
plt.axvline(
n_components_pca,
color="b",
label="PCA CV: %d" % n_components_pca,
linestyle="--",
)
plt.axvline(
n_components_fa,
color="r",
label="FactorAnalysis CV: %d" % n_components_fa,
linestyle="--",
)
plt.axvline(
n_components_pca_mle,
color="k",
label="PCA MLE: %d" % n_components_pca_mle,
linestyle="--",
)
# ⑮ 添加水平线表示 shrinkage 和 LedoitWolf 估计器的得分
plt.axhline(
shrunk_cov_score(X),
color="violet",
label="Shrunk Covariance MLE",
linestyle="-.",
)
plt.axhline(
lw_score(X),
color="orange",
label="LedoitWolf MLE" % n_components_pca_mle,
linestyle="-.",
)
# ⑯ 设置坐标轴标签、图例和标题
plt.xlabel("nb of components")
plt.ylabel("CV scores")
plt.legend(loc="lower right")
plt.title(title)
# 第 73 章 —— ⑰ 显示所有图形
plt.show()
这段代码展示了如何在同一数据集上使用交叉验证比较 PCA 与 FA 的模型选择性能。
源码路径：examples/decomposition/plot_pca_vs_lda.py - PCA 和 LDA 对比
# 第 73 章 —— ① 导入必要的模块
import matplotlib.pyplot as plt
from sklearn import datasets
from sklearn.decomposition import PCA
from sklearn.discriminant_analysis import LinearDiscriminantAnalysis
# 第 73 章 —— ② 加载鸢尾花数据集
iris = datasets.load_iris()
# 第 73 章 —— ③ 提取特征矩阵和目标向量
X = iris.data
y = iris.target
target_names = iris.target_names
# 第 73 章 —— ④ 实例化 PCA，保留 2 个主成分
pca = PCA(n_components=2)
# 第 73 章 —— ⑤ 拟合 PCA 模型并转换数据到主成分空间
X_r = pca.fit(X).transform(X)
# 第 73 章 —— ⑥ 实例化 LDA，保留 2 个判别分量（监督方法，需使用标签）
lda = LinearDiscriminantAnalysis(n_components=2)
# 第 73 章 —— ⑦ 拟合 LDA 模型并转换数据
X_r2 = lda.fit(X, y).transform(X)
# 第 73 章 —— ⑧ 打印前两个主成分解释的方差比例
print(
"explained variance ratio (first two components): %s"
% str(pca.explained_variance_ratio_)
)
# 第 73 章 —— ⑨ 创建第一个图形，用于可视化 PCA 投影结果
plt.figure()
# 第 73 章 —— ⑩ 定义颜色列表和线宽，用于按类别着色散点图
colors = ["navy", "turquoise", "darkorange"]
lw = 2
# 第 73 章 —— ⑪ 遍历每个类别，在 PCA 投影空间中绘制散点图
for color, i, target_name in zip(colors, [0, 1, 2], target_names):
plt.scatter(
X_r[y == i, 0], X_r[y == i, 1], color=color, alpha=0.8, lw=lw, label=target_name
)
# 第 73 章 —— ⑫ 显示图例，设置位置和样式
plt.legend(loc="best", shadow=False, scatterpoints=1)
# 第 73 章 —— ⑬ 设置图形标题
plt.title("PCA of IRIS dataset")
# 第 73 章 —— ⑭ 创建第二个图形，用于可视化 LDA 投影结果
plt.figure()
# 第 73 章 —— ⑮ 遍历每个类别，在 LDA 投影空间中绘制散点图
for color, i, target_name in zip(colors, [0, 1, 2], target_names):
plt.scatter(
X_r2[y == i, 0], X_r2[y == i, 1], alpha=0.8, color=color, label=target_name
)
# 第 73 章 —— ⑯ 显示图例
plt.legend(loc="best", shadow=False, scatterpoints=1)
# 第 73 章 —— ⑰ 设置图形标题
plt.title("LDA of IRIS dataset")
# 第 73 章 —— ⑱ 显示所有图形
plt.show()
这段代码展示了如何使用 PCA 和 LDA 对鸢尾花数据集进行二维投影，并通过可视化比较两种方法在类别分离上的差异。PCA 关注总方差最大化，而 LDA 监督地最大化类间方差。
源码路径：examples/decomposition/plot_faces_decomposition.py - FactorAnalysis 在人脸数据上的应用
# 第 73 章 —— ① 导入必要的模块
import logging
import matplotlib.pyplot as plt
from numpy.random import RandomState
from sklearn import cluster, decomposition
from sklearn.datasets import fetch_olivetti_faces
# 第 73 章 —— ② 设置随机状态以确保结果可复现
rng = RandomState(0)
# 第 73 章 —— ③ 配置日志输出格式和级别
logging.basicConfig(level=logging.INFO, format="%(asctime)s %(levelname)s %(message)s")
# 第 73 章 —— ④ 加载 Olivetti 人脸数据集，返回特征矩阵和目标（此处仅使用特征）
faces, _ = fetch_olivetti_faces(return_X_y=True, shuffle=True, random_state=rng)
n_samples, n_features = faces.shape
# 第 73 章 —— ⑤ 全局中心化：减去每个特征的均值，使每列均值为零
faces_centered = faces - faces.mean(axis=0)
# 第 73 章 —— ⑥ 局部中心化：进一步减去每个样本的均值，使每行均值为零
faces_centered -= faces_centered.mean(axis=1).reshape(n_samples, -1)
# 第 73 章 —— ⑦ 打印数据集中的人脸数量
print("Dataset consists of %d faces" % n_samples)
# 第 73 章 —— ⑧ 定义一个基础函数，用于以画廊形式可视化图像（如特征向量）
n_row, n_col = 2, 3
n_components = n_row * n_col
image_shape = (64, 64)
def plot_gallery(title, images, n_col=n_col, n_row=n_row, cmap=plt.cm.gray):
fig, axs = plt.subplots(
nrows=n_row,
ncols=n_col,
figsize=(2.0 * n_col, 2.3 * n_row),
facecolor="white",
constrained_layout=True,
)
fig.get_layout_engine().set(w_pad=0.01, h_pad=0.02, hspace=0, wspace=0)
fig.set_edgecolor("black")
fig.suptitle(title, size=16)
for ax, vec in zip(axs.flat, images):
vmax = max(vec.max(), -vec.min())
im = ax.imshow(
vec.reshape(image_shape),
cmap=cmap,
interpolation="nearest",
vmin=-vmax,
vmax=vmax,
)
ax.axis("off")
fig.colorbar(im, ax=axs, orientation="horizontal", shrink=0.99, aspect=40, pad=0.01)
plt.show()
# 第 73 章 —— ⑨ 使用定义好的函数，可视化居中后的人脸数据（前 n_components 张）
plot_gallery("Faces from dataset", faces_centered[:n_components])
# 第 73 章 —— ⑩ 实例化 PCA 估计器，使用随机 SVD 并启用白化
pca_estimator = decomposition.PCA(
n_components=n_components, svd_solver="randomized", whiten=True
)
# 第 73 章 —— ⑪ 拟合 PCA 模型到居中后的人脸数据
pca_estimator.fit(faces_centered)
# 第 73 章 —— ⑫ 可视化 PCA 学习到的特征（特征向量），即“特征脸”
plot_gallery(
"Eigenfaces - PCA using randomized SVD", pca_estimator.components_[:n_components]
)
# 第 73 章 —— ⑬ 实例化 NMF 估计器，设定容忍度
nmf_estimator = decomposition.NMF(n_components=n_components, tol=5e-3)
# 第 73 章 —— ⑪ 拟合 NMF 模型到原始非负人脸数据
nmf_estimator.fit(faces) # original non- negative dataset
# 第 73 章 —— ⑫ 可视化 NMF 学习到的非负成分
plot_gallery("Non-negative components - NMF", nmf_estimator.components_[:n_components])
# 第 73 章 —— ⑬ 实例化 FastICA 估计器，设定最大迭代次数、白化方式和容忍度
ica_estimator = decomposition.FastICA(
n_components=n_components, max_iter=400, whiten="arbitrary-variance", tol=15e-5
)
# 第 73 章 —— ⑭ 拟合 FastICA 模型到居中后的人脸数据
ica_estimator.fit(faces_centered)
# 第 73 章 —— ⑮ 可视化 FastICA 学习到的独立成分
plot_gallery(
"Independent components - FastICA", ica_estimator.components_[:n_components]
)
# 第 73 章 —— ⑯ 实例化 MiniBatchSparsePCA 估计器，设定稀疏性参数、最大迭代次数、批大小和随机状态
batch_pca_estimator = decomposition.MiniBatchSparsePCA(
n_components=n_components, alpha=0.1, max_iter=100, batch_size=3, random_state=rng
)
# 第 73 章 —— ⑰ 拟合 MiniBatchSparsePCA 模型到居中后的人脸数据
batch_pca_estimator.fit(faces_centered)
# 第 73 章 —— ⑱ 可视化 MiniBatchSparsePCA 学习到的稀疏成分
plot_gallery(
"Sparse components - MiniBatchSparsePCA",
batch_pca_estimator.components_[:n_components],
)
# 第 73 章 —— ⑲ 实例化 MiniBatchDictionaryLearning 估计器，设定成分数、正则化参数、最大迭代次数和批大小
batch_dict_estimator = decomposition.MiniBatchDictionaryLearning(
n_components=n_components, alpha=0.1, max_iter=50, batch_size=3, random_state=rng
)
# 第 73 章 —— ⑳ 拟合字典学习模型到居中后的人脸数据
batch_dict_estimator.fit(faces_centered)
# 第 73 章 —— ㉑ 可视化学习到的字典（每列为一个基向量）
plot_gallery("Dictionary learning", batch_dict_estimator.components_[:n_components])
# 第 73 章 —— ㉒ 实例化 MiniBatchKMeans 估计器，设定聚类数、容忍度、批大小、最大迭代次数和随机状态
kmeans_estimator = cluster.MiniBatchKMeans(
n_clusters=n_components,
tol=1e-3,
batch_size=20,
max_iter=50,
random_state=rng,
)
# 第 73 章 —— ㉓ 拟合 KMeans 模型到居中后的人脸数据
kmeans_estimator.fit(faces_centered)
# 第 73 章 —— ㉔ 可视化 KMeans 学习到的聚类中心
plot_gallery(
"Cluster centers - MiniBatchKMeans",
kmeans_estimator.cluster_centers_[:n_components],
)
# 第 73 章 —— ㉕ 实例化 FactorAnalysis 估计器，设定成分数和最大迭代次数
fa_estimator = decomposition.FactorAnalysis(n_components=n_components, max_iter=20)
# 第 73 章 —— ㉖ 拟合 FactorAnalysis 模型到居中后的人脸数据
fa_estimator.fit(faces_centered)
# 第 73 章 —— ㉗ 可视化 FactorAnalysis 学习到的因子（成分向量）
plot_gallery("Factor Analysis (FA)", fa_estimator.components_[:n_components])
# 第 73 章 —— ㉘ --- Pixelwise variance
# 第 73 章 —— ㉙ 创建图形，用于可视化 FactorAnalysis 估计的像素级噪声方差
plt.figure(figsize=(3.2, 3.6), facecolor="white", tight_layout=True)
# 第 73 章 —— ㉚ 从 FA 模型中提取噪声方差向量
vec = fa_estimator.noise_variance_
# 第 73 章 —— ㉛ 计算用于颜色映射的最大绝对值，确保正负值对称显示
vmax = max(vec.max(), -vec.min())
# 第 73 章 —— ㉜ 将噪声方差重塑为图像形状并以灰度图显示
plt.imshow(
vec.reshape(image_shape),
cmap=plt.cm.gray,
interpolation="nearest",
vmin=-vmax,
vmax=vmax,
)
# 第 73 章 —— ㉝ 隐藏坐标轴
plt.axis("off")
# 第 73 章 —— ㉞ 设置图形标题，解释这是像素级噪声方差图
plt.title("Pixelwise variance from \n Factor Analysis (FA)", size=16, wrap=True)
# 第 73 章 —— ㉟ 添加水平颜色条，说明噪声方差的大小
plt.colorbar(orientation="horizontal", shrink=0.8, pad=0.03)
# 第 73 章 —— ㊱ 显示图形
plt.show()
这段代码展示了如何使用 FactorAnalysis 对 Olivetti 人脸数据进行分解，提取潜在因子并可视化，同时估计每个像素的噪声方差，以评估模型在各方向上的假设噪声水平。
源码路径：examples/decomposition/plot_varimax_fa.py
# 第 73 章 —— ① 导入必要的模块
import matplotlib.pyplot as plt
import numpy as np
from sklearn.datasets import load_iris
from sklearn.decomposition import PCA, FactorAnalysis
from sklearn.preprocessing import StandardScaler
# 第 73 章 —— ② 加载鸢尾花数据集
data = load_iris()
# 第 73 章 —— ③ 对数据进行标准化（均值为0，方差为1），确保各特征贡献相当
X = StandardScaler().fit_transform(data["data"])
feature_names = data["feature_names"]
# 第 73 章 —— ④ 计算特征间的相关系数矩阵并以热图形式可视化
ax = plt.axes()
im = ax.imshow(np.corrcoef(X.T), cmap="RdBu_r", vmin=-1, vmax=1)
ax.set_xticks([0, 1, 2, 3])
ax.set_xticklabels(list(feature_names), rotation=90)
ax.set_yticks([0, 1, 2, 3])
ax.set_yticklabels(list(feature_names))
plt.colorbar(im).ax.set_ylabel("$r$", rotation=0)
ax.set_title("Iris feature correlation matrix")
plt.tight_layout()
plt.show()
# 第 73 章 —— ⑤ 设置要提取的成分数
n_comps = 2
# 第 73 章 —— ⑥ 定义要比较的方法列表：标准 PCA、未旋转的因子分析和 Varimax 旋转的因子分析
methods = [
("PCA", PCA()),
("Unrotated FA", FactorAnalysis()),
("Varimax FA", FactorAnalysis(rotation="varimax")),
]
# 第 73 章 —— ⑦ 创建包含多个子图的图形，用于并排比较不同方法的结果
fig, axes = plt.subplots(ncols=len(methods), figsize=(10, 8), sharey=True)
# 第 73 章 —— ⑧ 遍历每个方法
for ax, (method, fa) in zip(axes, methods):
# ⑨ 设置成分数
fa.set_params(n_components=n_comps)
# ⑩ 拟合模型到标准化后的数据
fa.fit(X)
# ⑪ 提取成分矩阵（转置以使特征在行上，成分在列上）
components = fa.components_.T
# ⑫ 打印当前方法的成分矩阵
print("\n\n %s :\n" % method)
print(components)
# ⑬ 计算成分矩阵中的最大绝对值，用于设置颜色映射范围
vmax = np.abs(components).max()
# ⑭ 以热图形式可视化成分矩阵，正负值用不同颜色表示
ax.imshow(components, cmap="RdBu_r", vmax=vmax, vmin=-vmax)
# ⑮ 设置 y 轴刻度标签为特征名称
ax.set_yticks(np.arange(len(feature_names)))
ax.set_yticklabels(feature_names)
# ⑯ 设置 x 轴刻度标签为成分名称
ax.set_title(str(method))
ax.set_xticks([0, 1])
ax.set_xticklabels(["Comp. 1", "Comp. 2"])
# 第 73 章 —— ⑰ 设置整个图形的总标题
fig.suptitle("Factors")
plt.tight_layout()
plt.show()
这段代码展示了如何在鸢尾花数据集上应用因子分析及其 Varimax 旋转，通过可视化成分载荷矩阵来帮助解释潜在结构，并与 PCA 结果进行对比。
73.4 独立成分与稀疏编码 —— 信号分离的“盲源拆解术”
核心概念
FastICA 通过最大化负熵（或最小化非高斯性）来估计统计独立的源信号，要求原始混合信号是独立非高斯源的线性组合。与 PCA 仅关注二阶统计（协方差）不同，ICA 追求更高阶统计的独立性。稀疏编码（SparseCoder）在给定字典的情况下，以 L1 正则化寻找最稀疏的系数向量；DictionaryLearning 则在同一过程里学习字典和稀疏系数，MiniBatchDictionaryLearning 通过小批量迭代提升大规模数据的学习效率。
73.4.1 时序图
73.4.2 核心类型定义
源码路径：examples/decomposition/plot_ica_blind_source_separation.py - FastICA（第 38-95 行）
# 第 73 章 —— ① 导入必要的模块
import numpy as np
from scipy import signal
# 第 73 章 —— ② 设置随机种子以确保结果可复现
np.random.seed(0)
# 第 73 章 —— ③ 定义时间点和样本数量
n_samples = 2000
time = np.linspace(0, 8, n_samples)
# 第 73 章 —— ④ 生成三个独立的非高斯源信号
s1 = np.sin(2 * time) # Signal 1 : sinusoidal signal
s2 = np.sign(np.sin(3 * time)) # Signal 2 : square signal
s3 = signal.sawtooth(2 * np.pi * time) # Signal 3: saw tooth signal
# 第 73 章 —— ⑤ 将信号合并为矩阵 S，每列为一个源信号
S = np.c_[s1, s2, s3]
# 第 73 章 —— ⑥ 向源信号添加高斯噪声以模拟真实测量
S += 0.2 * np.random.normal(size=S.shape) # Add noise
# 第 73 章 —— ⑦ 对每个源信号进行标准化（零均值，单位方差）
S /= S.std(axis=0) # Standardize data
# 第 73 章 —— ⑧ 定义混合矩阵 A，描述源信号如何线性组合形成观测值
A = np.array([[1, 1, 1], [0.5, 2, 1.0], [1.5, 1.0, 2.0]]) # Mixing matrix
# 第 73 章 —— ⑨ 生成观测数据 X，即混合后的信号
X = np.dot(S, A.T) # Generate observations
# 第 73 章 —— ⑩ 导入 PCA 和 FastICA 类
from sklearn.decomposition import PCA, FastICA
# 第 73 章 —— ⑪ 实例化 FastICA，设定成分数和白化方式
ica = FastICA(n_components=3, whiten="arbitrary-variance")
# 第 73 章 —— ⑫ 拟合 FastICA 模型到混合信号 X 并获取估计的源信号
S_ = ica.fit_transform(X) # Reconstruct signals
# 第 73 章 —— ⑬ 获取估计的混合矩阵
A_ = ica.mixing_ # Get estimated mixing matrix
# 第 73 章 —— ⑭ 验证重构误差：近似等于原始混合信号（考虑均值偏移）
assert np.allclose(X, np.dot(S_, A_.T) + ica.mean_)
# 第 73 章 —— ⑮ 为比较，计算 PCA 重构结果
pca = PCA(n_components=3)
H = pca.fit_transform(X) # Reconstruct signals based on orthogonal components
# 第 73 章 —— ⑯ 导入绘图库
import matplotlib.pyplot as plt
# 第 73 章 —— ⑰ 定义要绘制的模型和名称
models = [X, S, S_, H]
names = [
"Observations (mixed signal)",
"True Sources",
"ICA recovered signals",
"PCA recovered signals",
]
colors = ["red", "steelblue", "orange"]
# 第 73 章 —— ⑱ 创建多行单列的图形，依次绘制每个信号
plt.figure()
for ii, (model, name) in enumerate(zip(models, names), 1):
plt.subplot(4, 1, ii)
plt.title(name)
# ⑲ 遍历每个信号的通道，绘制时间序列
for sig, color in zip(model.T, colors):
plt.plot(sig, color=color)
# 第 73 章 —— ⑳ 调整子图间距并显示图形
plt.tight_layout()
plt.show()
这段代码演示了从混合信号到独立源的完整恢复过程，并验证了重构的正确性。
源码路径：examples/decomposition/plot_sparse_coding.py - SparseCoder（第 22-78 行）
# 第 73 章 —— ① 导入必要的模块
import matplotlib.pyplot as plt
import numpy as np
from sklearn.decomposition import SparseCoder
# 第 73 章 —— ② 定义 Ricker 小波函数（亦称 Mexican hat 小波）
def ricker_function(resolution, center, width):
"""Discrete sub-sampled Ricker (Mexican hat) wavelet"""
x = np.linspace(0, resolution - 1, resolution)
x = (
(2 / (np.sqrt(3 * width) * np.pi**0.25))
* (1 - (x - center) ** 2 / width**2)
* np.exp(-((x - center) ** 2) / (2 * width**2))
)
return x
# 第 73 章 —— ③ 定义生成 Ricker 小波字典的函数
def ricker_matrix(width, resolution, n_components):
"""Dictionary of Ricker (Mexican hat) wavelets"""
centers = np.linspace(0, resolution - 1, n_components)
D = np.empty((n_components, resolution))
for i, center in enumerate(centers):
D[i] = ricker_function(resolution, center, width)
D /= np.sqrt(np.sum(D**2, axis=1))[:, np.newaxis]
return D
# 第 73 章 —— ④ 设置信号分辨率、子采样因子和基波宽度
resolution = 1024
subsampling = 3 # subsampling factor
width = 100
n_components = resolution // subsampling
# 第 73 章 —— ⑤ 生成固定宽度和多宽度的 Ricker 小波字典
D_fixed = ricker_matrix(width=width, resolution=resolution, n_components=n_components)
D_multi = np.r_[
tuple(
ricker_matrix(width=w, resolution=resolution, n_components=n_components // 5)
for w in (10, 50, 100, 500, 1000)
)
]
# 第 73 章 —— ⑥ 生成测试信号：前四分之一为 3.0，其余为 -1.0 的阶跃信号
y = np.linspace(0, resolution - 1, resolution)
first_quarter = y < resolution / 4
y[first_quarter] = 3.0
y[np.logical_not(first_quarter)] = -1.0
# 第 73 章 —— ⑦ 定义要比较的稀疏编码方法列表（标题、算法、alpha、非零系数数、颜色）
estimators = [
("OMP", "omp", None, 15, "navy"),
("Lasso", "lasso_lars", 2, None, "turquoise"),
]
lw = 2
# 第 73 章 —— ⑧ 创建包含两个子图的图形，用于比较在两种字典下的稀疏编码效果
plt.figure(figsize=(13, 6))
for subplot, (D, title) in enumerate(
zip((D_fixed, D_multi), ("fixed width", "multiple widths"))
):
# ⑨ 设置子图标题
plt.subplot(1, 2, subplot + 1)
plt.title("Sparse coding against %s dictionary" % title)
# ⑩ 绘制原始信号作为参考（虚线）
plt.plot(y, lw=lw, linestyle="--", label="Original signal")
# ⑪ 遍历每种稀疏编码方法
for title, algo, alpha, n_nonzero, color in estimators:
# ⑫ 实例化 SparseCoder，使用指定字典和参数
coder = SparseCoder(
dictionary=D,
transform_n_nonzero_coefs=n_nonzero,
transform_alpha=alpha,
transform_algorithm=algo,
)
# ⑬ 对输入信号进行稀疏编码，获取稀疏系数向量
x = coder.transform(y.reshape(1, -1))
# ⑭ 计算非零系数的数量
density = len(np.flatnonzero(x))
# ⑮ 通过字典重构信号
x = np.ravel(np.dot(x, D))
# ⑯ 计算重构误差（平方和）
squared_error = np.sum((y - x) ** 2)
# ⑰ 绘制重构信号并添加图例标签
plt.plot(
x,
color=color,
lw=lw,
label="%s: %s nonzero coefs,\n%.2f error" % (title, density, squared_error),
)
# ⑱ 应用软阈值去偏方法作为基准
coder = SparseCoder(
dictionary=D, transform_algorithm="threshold", transform_alpha=20
)
x = coder.transform(y.reshape(1, -1))
_, idx = (x != 0).nonzero()
x[0, idx], _, _, _ = np.linalg.lstsq(D[idx, :].T, y, rcond=None)
x = np.ravel(np.dot(x, D))
squared_error = np.sum((y - x) ** 2)
# ⑲ 绘制去偏后的阈值结果
plt.plot(
x,
color="darkorange",
lw=lw,
label="Thresholding w/ debiasing:\n%d nonzero coefs, %.2f error"
% (len(idx), squared_error),
)
# ⑳ 设置坐标轴范围并显示图例
plt.axis("tight")
plt.legend(shadow=False, loc="best")
# 第 73 章 —— ㉑ 调整子图间距并显示图形
plt.subplots_adjust(0.04, 0.07, 0.97, 0.90, 0.09, 0.2)
plt.show()
这段代码展示了如何使用预定义的 Ricker 小波字典对一维信号进行稀疏分解，并重构回原始信号。
源码路径：examples/decomposition/plot_image_denoising.py - MiniBatchDictionaryLearning（第 88-138 行）
# 第 73 章 —— ① 导入必要的模块
import numpy as np
from scipy.datasets import face
# 第 73 章 —— ② 加载浣熊面部图像并转换为浮点格式（归一化到 [0,1]）
raccoon_face = face(gray=True)
raccoon_face = raccoon_face / 255.0
# 第 73 章 —— ③ 通过四邻域平均降采样以提高处理速度
raccoon_face = (
raccoon_face[::4, ::4]
+ raccoon_face[1::4, ::4]
+ raccoon_face[::4, 1::4]
+ raccoon_face[1::4, 1::4]
)
raccoon_face /= 4.0
height, width = raccoon_face.shape
# 第 73 章 —— ④ 在图像右半部分添加高斯噪声以模拟失真
print("Distorting image...")
distorted = raccoon_face.copy()
distorted[:, width // 2 :] += 0.075 * np.random.randn(height, width // 2)
# 第 73 章 —— ⑤ 导入绘图库
import matplotlib.pyplot as plt
# 第 73 章 —— ⑥ 定义辅助函数：显示图像及其与原始图像的差异
def show_with_diff(image, reference, title):
"""Helper function to display denoising"""
plt.figure(figsize=(5, 3.3))
plt.subplot(1, 2, 1)
plt.title("Image")
plt.imshow(image, vmin=0, vmax=1, cmap=plt.cm.gray, interpolation="nearest")
plt.xticks(())
plt.yticks(())
plt.subplot(1, 2, 2)
difference = image - reference
plt.title("Difference (norm: %.2f)" % np.sqrt(np.sum(difference**2)))
plt.imshow(
difference, vmin=-0.5, vmax=0.5, cmap=plt.cm.PuOr, interpolation="nearest"
)
plt.xticks(())
plt.yticks(())
plt.suptitle(title, size=16)
plt.subplots_adjust(0.02, 0.02, 0.98, 0.79, 0.02, 0.2)
# 第 73 章 —— ⑦ 显示失真图像及其与原始图像的差异
show_with_diff(distorted, raccoon_face, "Distorted image")
# 第 73 章 —— ⑧ 导入用于提取图像块的函数
from time import time
from sklearn.feature_extraction.image import extract_patches_2d
# 第 73 章 —— ⑨ 提取图像左半部分的所有参考块（无噪声区域）
print("Extracting reference patches...")
t0 = time()
patch_size = (7, 7)
data = extract_patches_2d(distorted[:, : width // 2], patch_size)
data = data.reshape(data.shape[0], -1)
data -= np.mean(data, axis=0)
data /= np.std(data, axis=0)
print(f"{data.shape[0]} patches extracted in %.2fs." % (time() - t0))
# 第 73 章 —— ⑩ 导入 MiniBatchDictionaryLearning 并学习字典
from sklearn.decomposition import MiniBatchDictionaryLearning
print("Learning the dictionary...")
t0 = time()
dico = MiniBatchDictionaryLearning(
# 增加迭代次数可提高质量但会增加训练时间
n_components=50,
batch_size=200,
alpha=1.0,
max_iter=10,
)
V = dico.fit(data).components_
dt = time() - t0
print(f"{dico.n_iter_} iterations / {dico.n_steps_} steps in {dt:.2f}.")
# 第 73 章 —— ⑪ 可视化学习到的字典（前 100 个原子）
plt.figure(figsize=(4.2, 4))
for i, comp in enumerate(V[:100]):
plt.subplot(10, 10, i + 1)
plt.imshow(comp.reshape(patch_size), cmap=plt.cm.gray_r, interpolation="nearest")
plt.xticks(())
plt.yticks(())
plt.suptitle(
"Dictionary learned from face patches\n"
+ "Train time %.1fs on %d patches" % (dt, len(data)),
fontsize=16,
)
plt.subplots_adjust(0.08, 0.02, 0.92, 0.85, 0.08, 0.23)
# 第 73 章 —— ⑫ 提取图像右半部分的噪声块（待重构区域）
print("Extracting noisy patches... ")
t0 = time()
data = extract_patches_2d(distorted[:, width // 2 :], patch_size)
data = data.reshape(data.shape[0], -1)
intercept = np.mean(data, axis=0)
data -= intercept
print("done in %.2fs." % (time() - t0))
# 第 73 章 —— ⑩ 定义要比较的变换算法列表（标题、算法、参数）
transform_algorithms = [
("Orthogonal Matching Pursuit\n1 atom", "omp", {"transform_n_nonzero_coefs": 1}),
("Orthogonal Matching Pursuit\n2 atoms", "omp", {"transform_n_nonzero_coefs": 2}),
("Least-angle regression\n4 atoms", "lars", {"transform_n_nonzero_coefs": 4}),
("Thresholding\n alpha=0.1", "threshold", {"transform_alpha": 0.1}),
]
reconstructions = {}
# 第 73 章 —— ⑪ 遍历每种变换算法
for title, transform_algorithm, kwargs in transform_algorithms:
print(title + "...")
reconstructions[title] = raccoon_face.copy()
t0 = time()
# ⑫ 设置当前变换参数
dico.set_params(transform_algorithm=transform_algorithm, **kwargs)
# ⑬ 对噪声块进行稀疏编码
code = dico.transform(data)
# ⑭ 用字典重构图像块
patches = np.dot(code, V)
# ⑮ 添加截距以校正均值偏移
patches += intercept
# ⑯ 将重构块 reshape 为图像格式
patches = patches.reshape(len(data), *patch_size)
# ⑩ 如果是阈值方法，则进行最小-最大归一化以避免负值
if transform_algorithm == "threshold":
patches -= patches.min()
patches /= patches.max()
# ⑪ 将重构后的右半部分合并回完整图像
reconstructions[title][:, width // 2 :] = reconstruct_from_patches_2d(
patches, (height, width // 2)
)
dt = time() - t0
print("done in %.2fs." % dt)
# ⑫ 显示重构结果及其与原始图像的差异
show_with_diff(reconstructions[title], raccoon_face, title + " (time: %.1fs)" % dt)
# 第 73 章 —— ⑬ 显示所有结果
plt.show()
这段代码演示了在图像去噪任务中，如何先学习稀疏字典，再通过稀疏编码重建噪声图像的右半部分。
73.5 流形学习探索 —— 瑞士卷与球面的“拓扑展开图”
核心概念
流形学习假设高维数据采样自低维流形，目标是保持流形内部的几何结构。
- LLE 通过在每个点的局部邻域内进行线性重构，保持局部线性关系；对噪声敏感但能完整保留拓扑。
- Isomap 首先构建 K‑近邻图并在图上计算测地线距离，再用 MDS 保持这些全局距离，适合有洞的流形（如瑞士卷）。
- MDS 仅依据成对距离进行嵌入，分为经典（保欧氏）和非度量（保序）两种；不涉及局部线性假设。
- SpectralEmbedding 使用拉普拉斯特征映射（图的振动模式），对噪声更鲁棒。
- t‑SNE 将高维相似度转化为低维条件概率，通过最小化 KL 散度保持局部邻居关系，perplexity 控制“局部 vs 全局”平衡。
73.5.1 数据流图
73.5.2 核心类型定义
源码路径：examples/manifold/plot_compare_methods.py - LocallyLinearEmbedding（第 57-110 行）
# 第 73 章 —— ① 导入必要的模块
import matplotlib.pyplot as plt
# 第 73 章 —— ② 用于 3D 可视化的导入（尽管未直接使用，但为 3D 投影所必需）
import mpl_toolkits.mplot3d # noqa: F401
from matplotlib import ticker
from sklearn import datasets, manifold
# 第 73 章 —— ③ 生成 S 曲线数据集
n_samples = 1500
S_points, S_color = datasets.make_s_curve(n_samples, random_state=0)
# 第 73 章 —— ④ 定义辅助函数：绘制 3D 散点图
def plot_3d(points, points_color, title):
x, y, z = points.T
fig, ax = plt.subplots(
figsize=(6, 6),
facecolor="white",
tight_layout=True,
subplot_kw={"projection": "3d"},
)
fig.suptitle(title, size=16)
col = ax.scatter(x, y, z, c=points_color, s=50, alpha=0.8)
ax.view_init(azim=-60, elev=9)
ax.xaxis.set_major_locator(ticker.MultipleLocator(1))
ax.yaxis.set_major_locator(ticker.MultipleLocator(1))
ax.zaxis.set_major_locator(ticker.MultipleLocator(1))
fig.colorbar(col, ax=ax, orientation="horizontal", shrink=0.6, aspect=60, pad=0.01)
plt.show()
# 第 73 章 —— ⑤ 定义辅助函数：绘制 2D 散点图
def plot_2d(points, points_color, title):
fig, ax = plt.subplots(figsize=(3, 3), facecolor="white", constrained_layout=True)
fig.suptitle(title, size=16)
add_2d_scatter(ax, points, points_color)
plt.show()
# 第 73 章 —— ⑥ 定义辅助函数：在 2D 坐标系上添加散点图
def add_2d_scatter(ax, points, points_color, title=None):
x, y = points.T
ax.scatter(x, y, c=points_color, s=50, alpha=0.8)
ax.set_title(title)
ax.xaxis.set_major_formatter(ticker.NullFormatter())
ax.yaxis.set_major_formatter(ticker.NullFormatter())
# 第 73 章 —— ⑦ 显示原始 S 曲线数据的 3D 可视化
plot_3d(S_points, S_color, "Original S-curve samples")
# 第 73 章 —— ⑧ 设置邻居数和目标维度
n_neighbors = 12 # neighborhood which is used to recover the locally linear structure
n_components = 2 # number of coordinates for the manifold
# 第 73 章 —— ⑨ 定义 LLE 的通用参数
params = {
"n_neighbors": n_neighbors,
"n_components": n_components,
"eigen_solver": "auto",
"random_state": 0,
}
# 第 73 章 —— ⑩ 实例化标准 LLE 方法
lle_standard = manifold.LocallyLinearEmbedding(method="standard", **params)
# 第 73 章 —— ⑪ 拟合并转换数据
S_standard = lle_standard.fit_transform(S_points)
# 第 73 章 —— ⑫ 实例化 LTSA LLE 方法
lle_ltsa = manifold.LocallyLinearEmbedding(method="ltsa", **params)
# 第 73 章 —— ⑬ 拟合并转换数据
S_ltsa = lle_ltsa.fit_transform(S_points)
# 第 73 章 —— ⑭ 实例化 Hessian LLE 方法
lle_hessian = manifold.LocallyLinearEmbedding(method="hessian", **params)
# 第 73 章 —— ⑮ 拟合并转换数据
S_hessian = lle_hessian.fit_transform(S_points)
# 第 73 章 —— ⑯ 实例化 modified LLE 方法
lle_mod = manifold.LocallyLinearEmbedding(method="modified", **params)
# 第 73 章 —— ⑰ 拟合并转换数据
S_mod = lle_mod.fit_transform(S_points)
# 第 73 章 —— ⑱ 创建包含四个子图的图形，用于比较四种 LLE 变体
fig, axs = plt.subplots(
nrows=2, ncols=2, figsize=(7, 7), facecolor="white", constrained_layout=True
)
fig.suptitle("Locally Linear Embeddings", size=16)
# 第 73 章 —— ⑲ 定义要比较的 LLE 方法及其显示名称
lle_methods = [
("Standard locally linear embedding", S_standard),
("Local tangent space alignment", S_ltsa),
("Hessian eigenmap", S_hessian),
("Modified locally linear embedding", S_mod),
]
# 第 73 章 —— ⑳ 遍历每个方法，在对应子图中绘制结果
for ax, method in zip(axs.flat, lle_methods):
name, points = method
add_2d_scatter(ax, points, S_color, name)
# 第 73 章 —— ㉑ 显示图形
plt.show()
这段代码演示了四种 LLE 方法在 S‑curve 数据上的应用，分别展示了不同局部线性假设的效果。
源码路径：examples/manifold/plot_compare_methods.py - Isomap（第 118-130 行）
# 第 73 章 —— ① 导入必要的模块（已在前文导入，此处假设可用）
# 第 73 章 —— ② 构造 Isomap 实例，设定邻居数、目标维度和闵可夫斯基距离参数
isomap = manifold.Isomap(n_neighbors=n_neighbors, n_components=n_components, p=1)
# 第 73 章 —— ③ 在 S 曲线数据上拟合并转换，利用测地线距离保持全局拓扑结构
S_isomap = isomap.fit_transform(S_points)
# 第 73 章 —— ④ 导入绘图库并绘制二维嵌入结果（假设 plot_2d 和 S_color 已定义）
plot_2d(S_isomap, S_color, "Isomap Embedding")
这段代码展示了 Isomap 如何通过测地线距离保持全局拓扑，将 S‑curve “展开” 成平面。
源码路径：examples/manifold/plot_compare_methods.py - MDS（第 132-152 行）
# 第 73 章 —— ① 导入必要的模块（已在前文导入，此处假设可用）
# 第 73 章 —— ② 实例化度量 MDS，设定目标维度、最大迭代次数、初始化方式等参数
md_scaling = manifold.MDS(
n_components=n_components,
max_iter=50,
n_init=1,
random_state=0,
init="classical_mds",
normalized_stress=False,
)
# 第 73 章 —— ③ 在 S 曲线数据上拟合并转换，保持高维空间中的成对欧氏距离
S_scaling_metric = md_scaling.fit_transform(S_points)
# 第 73 章 —— ④ 实例化非度量 MDS，关闭度量约束以保持非顺序距离关系
md_scaling_nonmetric = manifold.MDS(
n_components=n_components,
max_iter=50,
n_init=1,
random_state=0,
init="classical_mds",
metric_mds=False,
normalized_stress=False,
)
# 第 73 章 —— ⑤ 在 S 曲线数据上拟合并转换
S_scaling_nonmetric = md_scaling_nonmetric.fit_transform(S_points)
# 第 73 章 —— ⑥ 实例化古典 MDS，通过特征分解直接计算嵌入
md_scaling_classical = manifold.ClassicalMDS(n_components=n_components)
# 第 73 章 —— ⑦ 在 S 曲线数据上拟合并转换
S_scaling_classical = md_scaling_classical.fit_transform(S_points)
# 第 73 章 —— ⑧ 导入绘图库并创建包含三个子图的图形（假设所需函数和变量已定义）
fig, axs = plt.subplots(
nrows=1, ncols=3, figsize=(7, 3.5), facecolor="white", constrained_layout=True
)
fig.suptitle("Multidimensional scaling", size=16)
# 第 73 章 —— ⑨ 定义要比较的 MDS 方法及其显示名称
mds_methods = [
("Metric MDS", S_scaling_metric),
("Non-metric MDS", S_scaling_nonmetric),
("Classical MDS", S_scaling_classical),
]
# 第 73 章 —— ⑩ 遍历每个方法，在对应子图中绘制结果
for ax, method in zip(axs.flat, mds_methods):
name, points = method
add_2d_scatter(ax, points, S_color, name)
# 第 73 章 —— ⑪ 显示图形
plt.show()
这段代码比较了度量 MDS 与非度量 MDS 对同一流形的不同保距策略。
源码路径：examples/manifold/plot_compare_methods.py - SpectralEmbedding（第 154-164 行）
# 第 73 章 —— ① 导入必要的模块（已在前文导入，此处假设可用）
# 第 73 章 —— ② 实例化 SpectralEmbedding，设定目标维度、邻居数和随机状态
spectral = manifold.SpectralEmbedding(
n_components=n_components, n_neighbors=n_neighbors, random_state=42
)
# 第 73 章 —— ③ 在 S 曲线数据上拟合并转换，利用拉普拉斯特征映射捕捉数据的全局结构
S_spectral = spectral.fit_transform(S_points)
# 第 73 章 —— ④ 导入绘图库并绘制二维嵌入结果（假设 plot_2d 和 S_color 已定义）
plot_2d(S_spectral, S_color, "Spectral Embedding")
这段代码展示了拉普拉斯特征映射在保持局部相似性方面的表现。
源码路径：examples/manifold/plot_compare_methods.py - TSNE（第 166-176 行）
# 第 73 章 —— ① 导入必要的模块（已在前文导入，此处假设可用）
# 第 73 章 —— ② 实例化 t-SNE，设定目标维度、困惑度、初始化方式、最大迭代次数和随机状态
t_sne = manifold.TSNE(
n_components=n_components,
perplexity=30,
init="random",
max_iter=250,
random_state=0,
)
# 第 73 章 —— ③ 在 S 曲线数据上拟合并转换，通过最小化 KL 散度保持局部邻居结构
S_t_sne = t_sne.fit_transform(S_points)
# 第 73 章 —— ④ 导入绘图库并绘制二维嵌入结果（假设 plot_2d 和 S_color 已定义）
plot_2d(S_t_sne, S_color, "T-distributed Stochastic \n Neighbor Embedding")
这段代码体现了 t‑SNE 如何通过 KL 散度最小化在二维平面上保留局部邻居结构。
73.6 设计中的取舍
在实际应用中，我们常需要根据数据特性、计算资源和可解释性需求来选择合适的降维方法。当数据呈现明显的线性结构且解释性重要时，PCA 是首选，因其具备强理论基础、计算高效且结果易于解释。然而，当数据蕴藏非线性模式（如同心圆、瑞士卷）时，线性 PCA 无法有效捕捉其内在几何结构，此时应考虑核PCA或流形学习方法。核PCA通过将数据映射到高维特征空间实现非线性降维，能够处理复杂流形，但其计算复杂度为 O(n²)，对大规模数据不友好，且缺乏显式映射函数，限制了 out-of-sample 扩展。流形学习方法如 LLE 和 Isomap 能够更好地保持数据的局部或全局几何结构，适合揭示非凸流形，但它们对噪声较为敏感，且超参数（如邻居数）的调优较为困难；更重要的是，这些方法通常缺少显式的映射函数，使得对新样本进行 out‑of‑sample 投影变得困难。相比之下，t‑SNE 能够生成极具可视化吸引力的局部结构保持嵌入，但它不保留全局距离，其结果高度依赖于 perplexity 参数的经验调优，并且由于随机初始化，多次运行可能得到不同的结果，缺乏可重现性。因此，在实际项目中，若需构建可推广的模型，应优先考虑具有显式映射的方法（如 PCA、KernelPCA）；若仅关注可视化且数据规模适中，t‑SNE 是强有力的选择；若数据流形具有明显的几何结构（如孔洞）且可接受一定的计算成本，Isomap 提供了良好的全局结构保持；若主要关心局部几何保持且能容忍对噪声的敏感性，则 LLE 是合适的选择。
73.7 动手练习
-
比较 PCA、KernelPCA 与因子分析在人脸数据上的分解效果
-
使用
fetch_olivetti_faces加载 Olivetti 人脸数据。 -
分别训练
PCA(n_components=6),KernelPCA(kernel='rbf', gamma=15, n_components=6)与FactorAnalysis(n_components=6)。 -
可视化前 6 个主成分/因子/核主成分的重构图像；讨论 KernelPCA 在捕捉光照变化时的优势，FactorAnalysis 的因子与 PCA 主成分的差异。
-
-
探索 t‑SNE 在瑞士卷数据上的 perplexity 敏感性
-
生成瑞士卷数据
make_swiss_roll(n_samples=1500)。 -
分别使用 perplexity=5、30、50 运行
TSNE。 -
观察不同 perplexity 下嵌入的展开程度与簇分离度；思考 perplexity 越大，t‑SNE 越倾向全局结构，太小会导致局部碎片化。
-
-
对比 LLE、Isomap 与 MDS 在 S 曲线数据上的嵌入行为
-
生成 S 曲线
datasets.make_s_curve(n_samples=1000)。 -
使用
LocallyLinearEmbedding,Isomap(n_neighbors=10)与MDS（metric）进行二维嵌入。 -
可视化原始 3D 曲线与三种 2D 嵌入，检查是否出现撕裂或折叠。
-
分析：LLE 通过局部线性保持拓扑而不撕裂；Isomap 通过测地线保持全局结构；MDS 在仅靠距离矩阵时，若曲线具有非欧氏距离，可能产生折叠。
-
73.8 本章小结
本章我们围绕降维与流形学习展开，先从线性投影的 PCA 与 IncrementalPCA 入手，理解了方差最大化的正交投影机制以及大数据增量学习的实现方式；随后引入 KernelPCA，展示了核技巧如何把非线性结构映射到高维空间实现线性降维；接着深入 因子分析 与 FastICA，分别从生成模型与独立性角度解释潜变量的建模方式，并通过 SparseCoder 与 MiniBatchDictionaryLearning 说明稀疏编码在信号分离与图像去噪中的实际作用；随后我们系统比较了 LLE、Isomap、MDS、SpectralEmbedding 与 t‑SNE 在不同流形（S‑curve、瑞士卷、球面）上的表现，阐明了各算法保持局部或全局几何的原理以及它们的参数敏感性；最后通过设计取舍的讨论，帮助我们在实际项目中根据数据特性、计算资源与可解释性需求选择合适的降维方案。
以下表格总结了本章介绍的关键概念及其解释：
| 概念 | 解释 |
|------|------|
| PCA / IncrementalPCA | 线性降维，通过特征分解或 SVD 寻找方差最大的正交投影方向，IncrementalPCA 适用于大规模数据的增量学习 |
| KernelPCA | 非线性降维，利用核函数将数据映射至高维特征空间后执行线性 PCA，实现复杂流形的展开 |
| FactorAnalysis | 潜变量模型，假设观测数据由少数潜在因子线性组合加噪声，侧重解释协方差结构而非方差最大化 |
| FastICA | 盲源分离，通过最大化非高斯性估计统计独立的成分，适用于信号分离等场景 |
| SparseCoder / DictionaryLearning | 稀疏表示框架，学习或利用字典将信号表示为少数原子组合，用于去噪、特征提取等 |
| LLE | 局部线性嵌入，保持每个点的局部线性重构关系，适合保留流形拓扑 |
| Isomap | 测地线保持，将局部邻域图的最短路径（测地线）用于全局 MDS，能够展开有孔洞的流形 |
| MDS | 多维缩放，仅依据成对距离重构低维空间，分为经典（保欧氏）与非度量（保序） |
| SpectralEmbedding | 拉普拉斯特征映射，通过图的谱分解捕捉全局结构，对噪声更鲁棒 |
| t‑SNE | 概率流Manifold 嵌入，用 KL 散度最小化保留局部邻居，perplexity 决定局部 vs 全局关注度 |
在深入理解这些降维技术的原理与适用场景之后，我们将进入下一章的学习：预处理与特征工程 —— 数据转化的“炼金工坊”。在那一章中，我们将深入探讨特征缩放、离散化、稀疏编码与特征选择等技术，学习如何把原始数据转化为模型友好的高质量特征，为后续的建模工作奠定坚实基础。
73.9 模块地图/架构图
examples/decomposition/plot_pca_iris.py
├── __main__
examples/decomposition/plot_pca_vs_fa_model_selection.py
├── __main__
examples/decomposition/plot_pca_vs_lda.py
├── __main__
examples/decomposition/plot_faces_decomposition.py
├── __main__
examples/decomposition/plot_varimax_fa.py
├── __main__
examples/decomposition/plot_kernel_pca.py
├── __main__
examples/decomposition/plot_incremental_pca.py
├── __main__
examples/decomposition/plot_ica_blind_source_separation.py
├── __main__
examples/decomposition/plot_ica_vs_pca.py
├── __main__
examples/decomposition/plot_sparse_coding.py
├── __main__
examples/decomposition/plot_image_denoising.py
├── __main__
examples/manifold/plot_compare_methods.py
├── __main__
examples/manifold/plot_lle_digits.py
├── __main__
examples/manifold/plot_manifold_sphere.py
├── __main__
examples/manifold/plot_mds.py
├── __main__
examples/manifold/plot_swissroll.py
├── __main__
examples/manifold/plot_t_sne_perplexity.py
├── __main__
以上地图列出本章源码模块及其职责，后文将按数据流逐一解析。
第 74 章 —— 预处理与特征工程 —— 数据转化的“炼金工坊”
74.1 学习目标
-
难度：★★★☆☆（3/5）
-
预备知识：Python 基础、面向对象编程与 Markdown/代码阅读基础
-
理解特征缩放的目的和不同缩放器的适用场景
-
掌握非线性变换器（PowerTransformer、QuantileTransformer）的工作原理
浙公网安备 33010602011771号