简介一份围绕鸢尾花数据集展开的数理统计大作业完整文档适合正在学习多元统计、模式识别或机器学习入门课程的高校学生参考。资源包含1个docx文档压缩包约511KB内容系统演示了马氏距离、混合高斯模型、主成分分析、线性判别分析和刀切法在花萼、花瓣共四个属性数据上的综合应用。文档从符号说明与样本数据选取整理入手依次给出各类算法的原理推导与公式说明并通过PCA与LDA降维、GMM聚类以及留一法模型评估形成完整的数据处理分析和结果对比其中还讨论了不同降维程度与降维方式的差异能够帮助读者理解统计模型在实际数据上的建模、评估与选择思路。目前已有1426人学习下载可作为撰写数理统计实验报告、复现经典分类聚类方法或准备课程答辩时的参考样板。1. 五十年经典数据集遇上四把统计利器这套组合拳到底在验证什么一份数理统计大作业把鸢尾花数据集作为试验台同时用上马氏距离、高斯混合模型、主成分分析、线性判别分析和刀切法很容易被看作在罗列算法。但如果只盯着算法清单就会错过真正值得搞懂的东西这五样东西按「先降维、再建模、后判类、最后验证」的顺序排列组成了一条完整的统计决策链。任何一个环节换掉结果的含义都会变。鸢尾花数据集很小150 个样本、4 个数值属性、3 个类别却是数理统计里最合适的最小实验场。它属性之间存在线性相关量纲又不一致于是有了降维的必要属性基本符合高斯分布假设于是高斯混合模型能派上用场类别之间又不完全可分所以需要一个不偏不倚的评价方式来比较不同降维方案。读懂这份作业等于把「为什么数据决定方法方法又如何在数据上兑现」这条主线完整看了一遍比单独记住某个算法的公式更有复用价值。2. 数据与马氏距离为什么欧氏距离在四维特征上直接失灵2.1 鸢尾花数据的真实结构相关性与量纲如何影响距离鸢尾花数据集包含 150 个观察值分为 3 类每类 50 个每个观察值有花萼长度、花萼宽度、花瓣长度、花瓣宽度 4 个属性。数据从公开机器学习数据库中下载列的顺序是花萼长度、花萼宽度、花瓣长度、花瓣宽度最后一列是类别标签。符号上称自变量为 x它是 4 维输入向量n 表示样本个数取 150μ 表示均值向量Σ 表示协方差矩阵R 表示相关矩阵。这份数据的结构藏着一个陷阱花萼长度分布在 4.3 到 7.9 厘米之间花瓣宽度则只有 0.1 到 2.5 厘米的量级。如果用欧氏距离去度量两个样本的相似度数值范围更大的变量会在平方和中占绝对主导花萼长度和花瓣长度几乎决定距离而花瓣宽度的差异被淹没。更隐蔽的是属性之间的相关性。花萼长度和花瓣长度、花瓣宽度之间存在明显的线性关联这些变量反映的信息有重叠。欧氏距离假设各维度相互独立、等权重这个假设在鸢尾花数据上从头到尾都不成立。从类别结构看Setosa 这一类和另外两类在花瓣属性上有非常清晰的边界几乎可以线性切分Versicolour 与 Virginica 则在花瓣长度和花瓣宽度上有较大重叠。这意味着后续方法要解决的核心问题不是「能不能分」而是「在不完全可分的地方如何度量相似度、如何建模密度、如何评价降维后信息保留了多少」。这三层问题环环相扣第一环就落在距离度量上。2.2 马氏距离的本质旋转坐标系再按方差缩放马氏距离的定义是一个均值为 μ、协方差矩阵为 Σ 的多变量向量 x 到该分布中心的距离形式为 d(x, μ) sqrt((x − μ)^T Σ^(−1) (x − μ))。表面看只是把欧氏距离里的单位矩阵换成了 Σ^(−1)但这个替换改变了度量的几何含义它先把坐标系按协方差矩阵旋转到各主轴上再按每个方向上的方差做缩放最后才计算欧氏距离。用一个类比来说欧氏距离看的是「直线距离」马氏距离看的是「在这个分布里有多异常」。马氏距离有三个性质值得展开。第一当 Σ I 时各项独立且方差相同马氏距离退化为欧氏距离所以欧氏距离只是马氏距离的特例。第二马氏距离等价于对 x 和中心 μ 做标准化之后再求欧氏距离标准化的尺度由协方差矩阵决定因此它天然消除了量纲影响。第三它对变量的线性变换具有不变性也就是尺度无关无论把厘米换成英寸还是把数据整体做线性旋转距离排序不变。这三个性质合起来说明一件事当类别分布不是正圆而是旋转椭球时欧氏距离会同时忽略方向和展布马氏距离则把椭球的形状直接编码进度量里。在后续的判别流程里高斯混合模型会对每个类别拟合出一个高斯成分每个成分有自己的 μ 和 Σ。此时样本到底属于哪个类别最合理的度量就是到该成分中心的马氏距离因为 Σ 就是这个椭圆簇的形状。把马氏距离换回欧氏距离等于强行把所有簇当成正圆重叠区域里的样本会被系统性地误判。2.3 落成代码一个最小可用的马氏距离函数下面是实现马氏距离的 Python 函数输入一组样本、一个中心向量和一个协方差矩阵import numpy as np def mahalanobis_distance(x, mu, cov): # x: (n, d) 每个样本一行mu: (d,) 中心cov: (d, d) 协方差 # 用伪逆而不是逆防止协方差矩阵奇异时直接报错 cov_inv np.linalg.pinv(cov) diff x - mu.reshape(1, -1) left diff cov_inv return np.sqrt(np.einsum(ij,ij-i, left, diff))参数说明x 的形状是 (n, d)鸢尾花场景里通常是 150 行 4 列但降维后可能变成 150 行 2 列mu 是某个高斯成分的均值向量长度 dcov 是该成分的协方差矩阵形状 (d, d)。einsum 这一步计算的是 left 和 diff 逐行求内积相当于对角化地算出行向量 diff cov_inv diff.T 的结果避免用矩阵乘法再取对角线的写法。返回的是每个样本到该中心的马氏距离数组。提示用伪逆 pinv 只是一种兜底手段。当协方差矩阵奇异时它给出最小二乘意义下的距离避免程序崩溃但也会掩盖协方差估计本身退化的问题。真正要紧的是在建模阶段通过正则化防止协方差矩阵病态而不是在距离计算阶段补救。这套思路在后面 GMM 的 EM 迭代里会再次出现。3. 主成分分析与线性判别分析降维的两条路线差在标签用没用3.1 两种降维的分水岭方差最大化和可分性最大化主成分分析PCA和线性判别分析LDA都是线性降维方法形式上都是找一个线性变换把 p 维数据压到低维空间但目标函数完全不同。PCA 是无监督方法不理会类别标签只要求投影后的数据方差尽量大因为方差被认为是信息量的代理。LDA 是有监督方法由分类任务驱动目标是把同类样本的投影点拉近、把异类样本的投影点推远也就是同时放大类间离差、缩小类内离差。这个差异直接解释了为什么在有标签的情况下 LDA 理论上比 PCA 效果好PCA 找到的方差最大方向不一定能区分类别可能两个类别在该方向上大致重合投影后反而更难分LDA 的优化目标天然就是「判别友好」的它找到的方向直接服务于后续分类。在鸢尾花数据上三个类别的形状类似各类协方差结构差异不大正好落在 LDA 的假设范围内所以 LDA 的降维结果在分类任务上更有优势。需要强调 LDA 并不是在所有场景下都优于 PCA。它的监督信息来自标签如果标签噪声大、或类别数不准确LDA 会把标签里的错误也学进去PCA 则完全不依赖标签在无监督探索、数据可视化、去噪场景下依然是更稳健的起点。两者的选择本质上是「要不要为分类目标牺牲无偏性」的取舍。3.2 PCA 的源码实现用相关矩阵等价于先标准化再分解下面是 PCA 类的最小实现接口设计是 fit 学习投影方向、transform 执行投影class PCA: def __init__(self, n_components2): self.n_components n_components def fit(self, X): # 先做标准化(X - mean) / std self.mean_ X.mean(axis0) self.std_ X.std(axis0) X_std (X - self.mean_) / self.std_ # 标准化后的协方差矩阵就是相关矩阵 cov np.cov(X_std, rowvarFalse) # 形状 (d, d) eigvals, eigvecs np.linalg.eigh(cov) idx np.argsort(eigvals)[::-1] # 特征值降序 self.components_ eigvecs[:, idx[:self.n_components]] return self def transform(self, X): X_std (X - self.mean_) / self.std_ return X_std self.components_代码里的关键逻辑有两点。第一fit 里先标准化再求协方差矩阵得到的其实就是相关矩阵这避免花萼宽度0.1 到 2.5 厘米和花瓣长度1.0 到 6.9 厘米之间的量纲差异主导主成分方向。第二np.linalg.eigh 返回升序排列的特征值和特征向量所以要用 argsort 反转索引再取前 n_components 列构成投影矩阵。投影后的方差正好等于对应特征值因此特征值大小直接反映该主成分保留了多少信息量。提示PCA 的标准化步骤不能省。如果直接对原始协方差矩阵做特征分解数值范围大的变量会主导第一个主成分方向低方差但有用的变量可能被完全忽略。这不是实现细节问题而是 PCA 定义的一部分——用相关矩阵就是默认「每个原始变量先被压缩到单位方差再竞争信息量」。这份大作业里专门强调用相关矩阵代替协方差矩阵来消除量纲影响原因就在这。3.3 LDA 的源码实现广义特征值问题里藏着类别数上限LDA 的求解逻辑比 PCA 多了一步构造类内散度矩阵 Sw 和类间散度矩阵 Sb然后解广义特征值问题。Sw 衡量每个类别内部样本的离散程度Sb 衡量各类中心相对全局中心的离散程度。目标是让投影后的类间散度尽量大、类内散度尽量小数学上等价于最大化瑞利商其极大值由 S_w^(−1) S_b 的最大特征值给出。实现代码如下class LDA: def __init__(self, n_components2): self.n_components n_components def fit(self, X, y): classes np.unique(y) mu X.mean(axis0) Sw np.zeros((X.shape[1], X.shape[1])) Sb np.zeros((X.shape[1], X.shape[1])) for c in classes: Xc X[y c] muc Xc.mean(axis0) # 类内散度各类样本到类中心的离差平方和 Sw (Xc - muc).T (Xc - muc) n_c len(Xc) diff (muc - mu).reshape(-1, 1) # 类间散度类中心到全局中心按样本数加权 Sb n_c * diff diff.T # 广义特征值Sb w lambda Sw w用伪逆转成普通特征问题 eigvals, eigvecs np.linalg.eigh(np.linalg.pinv(Sw) Sb) idx np.argsort(eigvals)[::-1] self.components_ eigvecs[:, idx[:self.n_components]] return self def transform(self, X): return X self.components_代码中 Sw 累加每个类别的 (Xc − muc)^T (Xc − muc)它本质上是每个类别协方差矩阵的按样本数加权和。Sb 的构造用到了每个类中心 muc 与全局中心 mu 的差的外积再乘上该类的样本数 n_c。这里用 np.linalg.pinv(Sw) 而不是直接求逆目的和上一节一致Sw 在特征维度较高时可能不可逆伪逆让程序至少能跑通。但 LDA 有一个结构性的限制类别数为 C 时Sb 的秩最多是 C−1因此 S_w^(−1) S_b 的非零特征值也最多只有 C−1 个。对鸢尾花数据集C3LDA 最多给出 2 个有效投影方向。如果强行设置 n_components3第三个方向对应的是接近零的特征值方向由数值噪声决定投影结果没有任何统计意义。这也是后续分析里 LDA 降维只取 1 维或 2 维的原因——不是不想保留更多信息而是方法本身不允许。相比之下 PCA 可以降到 1 维、2 维或 3 维于是刀切法就有了对比的空间同样降到 1 维、2 维时PCA 和 LDA 哪个误判率更低。4. 高斯混合模型与极大似然软聚类如何把方差和相关性织进密度4.1 从 K 均值到 GMM硬聚类和软聚类的本质差别k-均值是最常见的聚类方法它把每个样本硬性归到距离最近的中心更新中心时只取该类内样本的均值。这个做法的局限在于它只关心中心位置完全忽略簇的形状。当不同类别在各个方向上的方差不一样、维度之间存在相关时k-均值会把旋转椭球形的簇当成球形来切分边界自然不准。高斯混合模型GMM解决了这个问题。它假设数据由 K 个高斯成分混合生成密度函数写为 p(x) Σπ_k N(x|μ_k, Σ_k)其中 π_k 是第 k 个成分的混合权重、满足 Σπ_k 1N(x|μ_k, Σ_k) 是均值为 μ_k、协方差矩阵为 Σ_k 的多元高斯密度。每个成分不仅学自己的中心还学一个完整的协方差矩阵因此可以表达不同方向上展布不同、维度间有相关的椭圆簇。在鸢尾花数据上4 个属性是连续数值、近似正态用高斯成分描述每个类别非常自然。与 k-均值最本质的差别是软分类GMM 不给每个样本只贴一个标签而是计算样本属于每个成分的后验概率称为责任。当数据点落在两个簇的重叠区域时它同时以一定概率归属两边。这种软输出不仅在聚类上更贴近真实分布也为后面的判别提供了后验概率可以和马氏距离无缝衔接。4.2 EM 的 E 步与 M 步参数更新的底层机制GMM 的参数没有解析解因为每个样本究竟由哪个成分生成是隐变量标准做法是极大似然意义下的 EM 迭代。每次迭代分两步E 步在固定参数时计算每个样本属于每个成分的责任M 步用责任作为权重重新估计 π、μ、Σ。具体代码形式如下def gmm_em(X, K3, max_iter200, tol1e-4, seed42): rng np.random.default_rng(seed) n, d X.shape pi np.ones(K) / K # 随机选 K 个样本作为初始中心更稳的做法是用 k-means 初始化 mu X[rng.choice(n, K, replaceFalse)] cov np.tile(np.eye(d), (K, 1, 1)) 1e-6 * np.eye(d) for _ in range(max_iter): resp np.zeros((n, K)) # E 步计算责任即第 k 个成分对样本 x_i 的后验概率 for k in range(K): resp[:, k] pi[k] * multivariate_normal_pdf(X, mu[k], cov[k]) resp / resp.sum(axis1, keepdimsTrue) # M 步按责任加权更新参数 Nk resp.sum(axis0) # 每个成分的等效样本数 pi Nk / n for k in range(K): mu[k] (resp[:, k:k1] * X).sum(axis0) / Nk[k] diff X - mu[k] cov[k] (resp[:, k:k1] * diff).T diff / Nk[k] cov[k] 1e-6 * np.eye(d) # 正则防止协方差奇异 return pi, mu, covE 步的关键在责任矩阵 resp 的归一化先按每个成分的似然加权算出原始值再对每行求和归一保证每个样本对各个成分的归属概率加起来等于 1。M 步里 Nk 是第 k 个成分累积的责任直观理解就是「有多少个等效样本被这个成分认领」。μ 的更新是责任加权的均值Σ 的更新是责任加权的外积和除以 Nk相当于对协方差矩阵做加权最大似然估计。收敛判据通常用相邻两次迭代的对数似然增量小于 tol 就停止。注意协方差矩阵加 1e-6 的对角项不是可有可无的细节。当某个成分的等效样本数很少或者数据在高维空间里稀疏时Σ 会趋向半正定甚至奇异这时候马氏距离里的 Σ^(−1) 会出现极大值分类结果瞬间失真。正则项相当于给协方差设了一个下限让数值计算保持稳定。4.3 降维、GMM、马氏距离如何构成一条完整预测管线整个分析流程是原始 4 维数据先经过 PCA 或 LDA 降到 12 维然后用 GMM 在这个低维空间拟合 3 个高斯成分每个成分对应一个类别。对新样本做预测时把它投影到同样的低维空间分别计算它到 3 个成分中心的马氏距离距离最近的那个成分就是预测类别。这个判断方式和直接用后验概率 p(k|x) 等价因为多元高斯密度的对数恰恰是负的马氏距离平方的一半再减去一个与方向有关的常数项对同一个样本比较不同成分时常数项虽有差异但主要排序由马氏距离决定。降维在这里不是可有可无的预处理而是直接影响 GMM 的训练质量。4 维时每个高斯成分的协方差矩阵是 4×4包含 4 个方差和 6 个协方差参数3 个成分一共 30 个参数要靠 150 条样本估计方差很大。降到 2 维后每个协方差矩阵只有 5 个独立参数3 个成分共 15 个参数估计稳定性明显提升。同时降维丢弃的维度往往对应噪声和冗余信息GMM 在这些方向上不必再拟合无意义的分布。这就是「降维同时去噪」在概率建模层面的真实含义。5. 刀切法与避坑手册留一法如何拆穿降维的自我感觉良好5.1 刀切法的机制每个样本都有一次「没见过」的机会刀切法在这里的本质就是留一法交叉验证。具体步骤是从某一类 Gi 的 n_i 个训练样本中剔除一个用剩下 n_i−1 个样本重新建立判别规则然后用这条规则对剔除的那一个样本判类如此重复直到该类每个样本都被判过一次记录错判个数 k_i所有类别执行相同操作后总误判率等于总错判数除以总样本数。为什么不用回报法回报法用全部样本训练、再对全部样本测试训练集和测试集完全重叠相当于「考自己出过的题」。模型只是记住了数据误判率低得没有意义。刀切法让每个样本都经历一次「模型没在训练时见过我」的考验得到的误判率是对泛化能力更合理的估计。它相比 K 折交叉验证的差别在于每次只留一个样本训练集与完整数据集只差 1 条模型学到的信息损失最小适合 150 条样本这样的小数据集代价是运算量相对大因为整个「降维 GMM 距离判类」流程要重新执行 150 次。5.2 刀切法循环里最容易被忽视的顺序问题刀切法代码结构看似简单但顺序问题会直接毁掉结论。正确做法是把降维、GMM 拟合全部放进循环里每次只用训练集的样本完成这些拟合测试样本只参与最后的距离计算def jackknife_evaluate(fit_predict, X, y): n len(y) errors 0 for i in range(n): train_idx np.delete(np.arange(n), i) test_idx i pred fit_predict( X[train_idx], y[train_idx], X[test_idx:test_idx1] ) errors (pred[0] ! y[i]) return errors / nfit_predict 内部依次完成降维拟合、GMM 拟合、马氏距离判类所有拟合参数都只从训练样本中得到。如果为了省时间先在整个数据集上做 PCA再用降维后的结果循环留一属于信息泄漏——降维阶段已经看过测试样本的分布测试误差会被系统性低估。刀切法比较不同降维方案时必须保证每种方案都在同等的信息条件下被评估否则「PCA 误判率比 LDA 低」这类结论可能是偷看数据造成的假象。5.3 我在鸢尾花任务上踩过的五个坑第一个坑直接用欧氏距离判类错判率明显偏高。现象是 Setosa 与两类边界清晰的样本没问题但 Versicolour 与 Virginica 重叠区域的样本大量分错。原因是欧氏距离忽略特征相关性和量纲把椭圆簇当成球簇。解决方法是改用马氏距离并配合 GMM 学出的协方差矩阵。第二个坑EM 结果每次跑都不一样。现象是多次运行后有时两个高斯成分挤在同一簇里有时某个类别的样本没有任何成分认领。原因是随机初始化让小概率的局部最优反复出现。解决方法是固定随机种子最好用 k-均值的结果初始化 μ再进入 EM 迭代同时保留协方差正则项防止单个成分的等效样本数趋近零。第三个坑不降维直接用 GMM 建模误判率反而波动大。现象是 4 维下交叉验证的结果时好时坏。原因是 4 维协方差矩阵参数多、150 条样本估计不稳过拟合风险高。解决方法是先用 PCA 或 LDA 降到 2 维让每个成分的参数降到 5 个。第四个坑LDA 降维设置 n_components3。现象是第三个投影方向上的特征值几乎为零降维后 GMM 在这条轴上疯狂拟合噪声。原因是三类的类间散度矩阵秩最多为 2不可能存在第三个有效方向。解决方法是把 LDA 的降维维度上限设为类别数减一即 2。第五个坑在整个数据集上先做 PCA 再做刀切法。现象是误判率低得惊人但换成新样本后表现明显变差。原因是降维阶段的均值、标准差和特征向量偷看了所有样本的信息测试样本的信息已经参与训练。解决方法是在刀切法循环内部重新执行完整的降维拟合确保每个测试样本在训练阶段完全不可见。6. 设计取舍的最后一道思考四种方法怎么咬合成一个自洽框架把这四种方法拆开看各自解决不同层面的问题马氏距离是度量决定「怎么算相似」PCA 和 LDA 是特征变换决定「在哪个空间里算相似」GMM 是密度建模决定「用什么样的分布假设描述每个类别」刀切法是验证策略决定「凭什么相信前面的选择」。它们不在同一个层面上却恰好覆盖了统计分析从描述到推断的完整路径。真正有价值的不是某一算法的公式而是这条流程的先后关系先根据数据相关性选度量再用降维控制模型复杂度然后建概率模型最后用留一法检验泛化能力。这个框架也有明显的边界。PCA 不考虑标签对分类任务可能保留无用方向LDA 假设各类协方差相近如果类别形状差异很大线性判别会失真GMM 假设每个簇可被高斯刻画数据严重偏态时拟合效果有限刀切法在小样本上偏差小但方差大且运算成本高。进阶路线也清晰可见用 SVD 直接做 PCA 可以绕开协方差矩阵的条件数问题若怀疑各类协方差不同可以改用二次判别分析QDA它相当于让每个类别保留自己的协方差矩阵这和 GMM 的思想同源若要决定 GMM 的成分个数可以用 BIC 对模型复杂度做权衡而不是固定 3 个成分。把这些变体对照着读原本孤立的公式会连成一张网。从那以后我每次拿到一份新的数值型数据都会强制走一遍「看相关性 → 定度量 → 选降维 → 建概率模型 → 留一验证」这个流程先不管算法多高级把每个环节的假设和数据对齐再谈调参。踩过的坑少了结论也更能站住脚。希望帮到你。本文还有配套的精品资源点击获取