简介这份数理统计大作业文档面向高校统计学、机器学习课程的学习者围绕经典鸢尾花数据集展开完整分析帮助读者理解多重变量数据的建模流程。内容涵盖马氏距离度量样本相似度、混合高斯模型聚类、主成分分析降维去噪、线性判别分析兼顾类间与类内差异以及刀切法评估降维方法与模型泛化能力并附有数据预处理与结果讨论。资源包共1个docx文件约511KB结构完整、论述连贯适合作为课程作业参考或复习数理统计知识点的范本。目前已有1429人学习下载读者可从中获得从算法原理到数据处理分析的完整思路对照梳理各方法的适用场景与实现要点快速搭建自己的分析框架。1. 从一份鸢尾花大作业说起为什么降维和聚类要放在同一套流程里数理统计大作业最容易翻车的地方不是公式不会推而是把 PCA、LDA、GMM 各写一个 demo最后拼不出一条完整的分析链。鸢尾花数据集只有 150 个样本、4 个特征、3 个类别看起来简单但它恰好能把“降维—建模—评估”这条主线完整跑一遍先用主成分分析PCA和线性判别分析LDA把四维数据压到低维再在高斯混合模型GMM下做极大似然估计最后用刀切法留一法比较不同降维程度和方法之间的差异。这套方案适合正在做数理统计课程设计、机器学习毕设开题或者需要交一份可复现分析报告的人。下面按需求、设计、实现、避坑、验收的顺序拆开讲代码和参数都能直接改。2. 需求与数据底座鸢尾花数据集到底要解决什么问题2.1 四个特征、三个类别与一条分析主线鸢尾花数据集由 Fisher 在 1936 年整理1988 年公开包含 150 个观测每类 50 个特征是花萼长度、花萼宽度、花瓣长度、花瓣宽度单位都是厘米标签是 Setosa、Versicolour、Virginica 三类。这个数据集的价值不在于“分类准确率刷到多高”而在于它同时满足几个条件特征都是数值型、量纲一致、近似服从高斯分布、特征之间存在明显线性相关。这几点决定了后面可以放心用马氏距离、GMM 和 PCA。大作业要回答的问题可以拆成三层。第一层是降维四个特征里花瓣长度和花瓣宽度高度相关花萼长度和花萼宽度也有一定相关性直接建模会引入冗余和噪声需要 PCA 或 LDA 压到一维或二维。第二层是建模降维后的数据用 GMM 拟合分布在极大似然意义下估计均值、协方差和混合权重再用样本到各高斯成分均值的马氏距离判定类别。第三层是评估降维到几维、用哪种降维方法不能拍脑袋要用刀切法算判错率来比。我一般会把整份作业的符号先固定下来避免后面公式和代码对不上。样本矩阵记为 $D$样本数 $n150$输入向量 $x(x_1,x_2,x_3,x_4)^T$ 对应四个特征类别 $y\in{0,1,2}$均值向量 $\mu$协方差矩阵 $\Sigma$相关矩阵 $R$线性变换方向 $u$。符号表看着是形式主义但写到第四章比对结果时你会发现没有它根本说不清“第几主成分对应哪个特征”。2.2 数据读取与预处理的工程化写法原始数据从公开机器学习仓库获取格式是逗号分隔、无表头、最后一列是类别。读取时不要手写 150 行直接用 pandas 加列名顺手把类别映射成 0/1/2。预处理有两个关键动作一是检查缺失值和重复值鸢尾花数据一般干净但作业里要写出检查过程二是标准化PCA 用相关矩阵代替协方差矩阵来消除量纲影响代码里体现为对特征做 z-score。import pandas as pd import numpy as np COLUMNS [sepal_len, sepal_wid, petal_len, petal_wid, label] LABEL_MAP {Iris-setosa: 0, Iris-versicolor: 1, Iris-virginica: 2} def load_iris(path): df pd.read_csv(path, headerNone, namesCOLUMNS) # 类别字符串映射为 0/1/2方便后续按整数索引 df[label] df[label].str.strip().map(LABEL_MAP) assert df[label].notna().all(), 存在未识别的类别标签 assert df.shape[0] 150, 样本数应为 150 return df def standardize(X): # 用总体标准差做 z-score等价于用相关矩阵做 PCA mu X.mean(axis0) sigma X.std(axis0, ddof0) return (X - mu) / sigma, mu, sigma这段代码里ddof0是有意为之PCA 用相关矩阵时标准化要除以总体标准差和协方差矩阵的特征分解保持口径一致。如果这里用ddof1特征值会有细微偏差虽然对鸢尾花这种量级影响不大但答辩时被问到“你的相关矩阵怎么算的”会答不圆。assert不是摆设数据行数或标签映射出错时能第一时间拦住比后面 GMM 跑出诡异结果再回头查要省事得多。3. 算法链路设计PCA、LDA、GMM 与刀切法怎么串3.1 马氏距离把协方差结构带进分类判据马氏距离是整条链路的粘合剂。欧氏距离假设各方向方差相同且特征独立鸢尾花显然不满足花瓣长度的方差远大于花萼宽度特征之间还相关。马氏距离的定义是 $d(x,\mu)\sqrt{(x-\mu)^T\Sigma^{-1}(x-\mu)}$当 $\SigmaI$ 时退化为欧氏距离本质上是把数据白化后再算欧氏距离所以它无量纲、不受量纲变化影响。在方案里马氏距离承担两个职责。一是 GMM 内部计算样本到各高斯成分均值的“归属强度”虽然 GMM 的 E 步用的是带协方差的行列式和指数项但判别阶段用马氏距离更直观二是最终分类规则样本归到马氏距离最小的那个成分对应的类别。这里有个容易忽略的点协方差矩阵要加正则项否则某类样本在降维后维度接近样本数时$\Sigma$ 可能奇异求逆直接报错。通常的做法是 $\Sigma\lambda I$$\lambda$ 取 $10^{-6}$ 到 $10^{-3}$ 量级。3.2 高斯混合模型极大似然估计与 EM 迭代GMM 把总体分布写成 $p(x)\sum_{k1}^{K}\pi_k \mathcal{N}(x|\mu_k,\Sigma_k)$$K$ 是成分数$\pi_k$ 是混合权重。对数似然 $\sum_{i}\log p(x_i)$ 没有闭式解用 EM 迭代E 步算后验责任度 $\gamma_{ik}$M 步更新 $\pi_k$、$\mu_k$、$\Sigma_k$。鸢尾花三类$K$ 取 3 是自然选择但作业里值得加一组 $K2$ 和 $K4$ 的对照说明 $K$ 不是拍出来的。from scipy.stats import multivariate_normal def gmm_em(X, K3, max_iter200, tol1e-6, reg1e-6): n, d X.shape # 用分位数初始化均值避免所有成分挤在一起 idx np.linspace(0, n - 1, K).astype(int) mu X[np.argsort(X[:, 0])][idx].copy() sigma [np.cov(X.T) reg * np.eye(d) for _ in range(K)] pi np.ones(K) / K prev_ll -np.inf for _ in range(max_iter): # E 步责任度矩阵行和为 1 resp np.zeros((n, K)) for k in range(K): resp[:, k] pi[k] * multivariate_normal.pdf(X, mu[k], sigma[k]) ll np.log(resp.sum(axis1) 1e-12).sum() resp / resp.sum(axis1, keepdimsTrue) # M 步加权更新权重、均值、协方差 nk resp.sum(axis0) pi nk / n mu (resp.T X) / nk[:, None] for k in range(K): diff X - mu[k] sigma[k] (resp[:, k][:, None] * diff).T diff / nk[k] reg * np.eye(d) if abs(ll - prev_ll) tol: break prev_ll ll return pi, mu, sigma, resp参数说明reg是协方差正则项防止奇异tol控制对数似然收敛阈值鸢尾花这种规模 200 次迭代足够初始化用排序后等间隔取点比纯随机稳能减少 EM 陷入局部最优的概率。这段代码没有调库的GaussianMixture是为了让作业里“极大似然估计”这句话有对应的实现答辩时能指着 E 步 M 步讲。3.3 PCA 与 LDA无监督降维和有监督降维的对照PCA 的目标是找线性变换 $u$使投影方差最大解是协方差矩阵或相关矩阵最大特征值对应的特征向量。鸢尾花四维降到二维通常前两个主成分能解释 95% 以上的方差。LDA 则利用标签目标是最大化类间散度与类内散度的比值解是 $S_w^{-1}S_b$ 的最大特征值对应特征向量。两者放在一起对比正好说明“有标签时 LDA 通常优于 PCA”这个结论而不是空口断言。def pca_fit(X, n_components2): R np.corrcoef(X, rowvarFalse) # 相关矩阵消除量纲 eigvals, eigvecs np.linalg.eigh(R) order np.argsort(eigvals)[::-1] eigvals, eigvecs eigvals[order], eigvecs[:, order] return eigvecs[:, :n_components], eigvals def lda_fit(X, y, n_components2): classes np.unique(y) mu_all X.mean(axis0) Sw np.zeros((X.shape[1], X.shape[1])) Sb np.zeros_like(Sw) for c in classes: Xc X[y c] muc Xc.mean(axis0) Sw (Xc - muc).T (Xc - muc) Sb len(Xc) * np.outer(muc - mu_all, muc - mu_all) eigvals, eigvecs np.linalg.eigh(np.linalg.inv(Sw 1e-6 * np.eye(Sw.shape[0])) Sb) order np.argsort(eigvals)[::-1] return eigvecs[:, order[:n_components]], eigvals[order]pca_fit用np.corrcoef而不是np.cov对应前面说的“用相关矩阵代替协方差矩阵消除量纲”。lda_fit里Sw加正则项同理鸢尾花每类 50 个样本、四维特征Sw本来可逆但降维到一维再算时数值稳定性会变差加正则是个保险。LDA 最多只能降到 $C-12$ 维这一点在作业里要写清楚否则有人试图降到三维会得到无意义的第三个方向。3.4 刀切法用留一法给降维方案打分刀切法留一法的流程是每次从训练集剔除一个样本用剩余样本建立判别规则再预测被剔除的样本重复 150 次统计判错个数判错率 $ \sum k_i / n$。它比简单回报法更客观因为每个样本都当过测试集比 k 折交叉验证更彻底因为训练集几乎用满。代价是运算量大150 次 GMM 拟合在普通笔记本上要跑几分钟但完全可接受。方案里刀切法要比较的对象是“降维方法 × 降维维度”的组合PCA 降到 1 维、2 维LDA 降到 1 维、2 维再加一组不降维的原始四维作为基线。每组跑一遍留一法记录判错率最后画一张对比表。这张表就是整份作业的结论支撑比任何文字描述都有说服力。4. 实现与结果分析从类设计到对比表4.1 模块划分与类接口按作业里的 UML 思路PCA 和 LDA 各设计成一个类fit接收训练样本求变换矩阵transform做投影fit_transform串起来。GMM 单独一个类fit跑 EMpredict用马氏距离判类。刀切法写成函数接收“降维器 分类器”的组合返回判错率。这样拆的好处是刀切法里换降维方法只改一个参数不用复制粘贴五份代码。class PCAReducer: def __init__(self, n_components2): self.n_components n_components def fit(self, X, yNone): self.W_, self.eigvals_ pca_fit(X, self.n_components) return self def transform(self, X): return X self.W_ class LDAReducer: def __init__(self, n_components2): self.n_components n_components def fit(self, X, y): self.W_, self.eigvals_ lda_fit(X, y, self.n_components) return self def transform(self, X): return X self.W_接口统一后刀切法的主循环可以写成通用形式对每个样本用剩余样本fit降维器和 GMM再transform被留出的样本、算马氏距离、取最小距离对应的类别。注意 LDA 的fit需要标签PCA 不需要所以通用循环里要把y传进去PCA 类用yNone兼容。4.2 降维维度与判错率的对比结果跑完五组配置典型结果会呈现这样的规律原始四维基线判错率在 2% 到 4% 之间PCA 降到二维判错率接近基线降到一维明显变差因为第一主成分主要承载花瓣长度方向的信息Setosa 和另外两类能分开但 Versicolour 和 Virginica 有重叠LDA 降到二维判错率最低通常能到 2% 以下降到一维也比 PCA 一维好因为它用了标签投影方向直接冲着类间可分性去。配置维度判错数约判错率原始特征 GMM442.7%PCA GMM253.3%PCA GMM1128.0%LDA GMM221.3%LDA GMM164.0%这张表要配合解释PCA 二维接近基线说明四维里确实有冗余两个主成分保住了大部分判别信息PCA 一维掉得厉害说明“方差最大”不等于“判别性最强”LDA 二维最优印证了有监督降维在分类任务上的优势。数字不必和上表完全一致但趋势要能自圆其说答辩时被追问“为什么你的 PCA 一维这么差”要能答到“第一主成分被花瓣长度主导类间重叠没被拉开”。4.3 结果解读里容易漏掉的两个细节第一个细节是标准化对 PCA 的影响。如果不标准化直接用协方差矩阵花瓣长度的方差远大于其他特征第一主成分几乎就是花瓣长度本身降维效果和标准化后完全不同。作业里要明确写“采用相关矩阵”并在附录放一张标准化前后的特征值对比说明为什么这么选。第二个细节是 GMM 的类别分配。GMM 是无监督的跑出来的成分编号和真实标签没有对应关系判类前要先做一次匹配比如用混淆矩阵找最优排列或者直接用马氏距离到各类均值的最近邻规则绕开编号问题。我一般选后者逻辑简单答辩时也好解释GMM 负责估计分布参数判别用马氏距离两步分开说。5. 避坑与排查这类大作业最容易翻车的五个点5.1 协方差矩阵奇异导致求逆报错现象是 GMM 迭代到一半抛LinAlgError: Singular matrix或者马氏距离算出inf。原因通常是某类样本在降维后维度接近或超过样本数或者某个高斯成分只分到极少数样本协方差估计退化。解决办法是给协方差加正则项 $\lambda I$$\lambda$ 从 $10^{-6}$ 起调同时在 EM 里加一个最小样本数保护某成分责任度之和小于阈值就重新初始化而不是硬算。5.2 标准化口径不一致导致 PCA 和 LDA 结果对不上现象是同一份数据PCA 用相关矩阵、LDA 用原始协方差两者降维后的尺度差一个量级对比表看起来像两组不同数据。原因是 PCA 前做了 z-scoreLDA 前没做。解决方法是整条链路统一先标准化一次PCA 和 LDA 都在标准化后的数据上跑LDA 的Sw、Sb自然也在同一尺度下。这一点在代码里体现为standardize只调用一次结果传给两个降维器。5.3 刀切法循环里重复拟合造成结果不可复现现象是同一份代码跑两次判错数不一样。原因是 GMM 初始化带随机性而刀切法每次循环都重新fit随机种子没固定。解决办法是在刀切法入口设np.random.seed(42)并且 GMM 初始化改成确定性的分位数取点如 3.2 节代码所示。如果坚持随机初始化就在每次fit前重置种子保证 150 次循环的随机序列可复现。5.4 LDA 降维维度超过 C-1 得到无意义方向现象是试图把鸢尾花 LDA 降到三维第三个特征值接近零投影后数据几乎塌缩。原因是 LDA 的类间散度矩阵秩最多为 $C-1$三类数据最多两个有效判别方向。解决办法是在lda_fit里加断言n_components len(classes) - 1作业正文里也写明这个约束避免读者照抄时踩同一个坑。5.5 判错率统计把训练样本和测试样本混在一起现象是判错率异常低接近 0。原因是判别规则用了被预测样本本身的信息比如先对所有样本拟合 GMM 再留一预测等于测试集泄漏进训练。解决办法是严格按刀切法流程每次剔除一个样本用剩余 149 个样本重新fit降维器和 GMM再预测被剔除的那个。运算量是大但这是留一法的定义不能偷懒。6. 验收与答辩怎么证明这套分析链站得住验收这份大作业我一般看四样东西按重要性排序。第一是流程完整性数据读取、标准化、PCA、LDA、GMM、马氏距离判别、刀切法比对七个环节缺一不可缺了哪个答辩时被问到就是硬伤。第二是结果可复现代码里固定随机种子附录给出运行环境和依赖版本别人拿到能跑出同样的判错率表。第三是对比有解释PCA 一维为什么差、LDA 二维为什么好要落到特征相关性和标签信息上不能只说“实验结果符合预期”。第四是符号和公式自洽第二章的符号表、第三章的公式、第四章的代码变量名要对得上这是数理统计作业区别于纯编程作业的地方。答辩时最常被追问的三个问题提前准备好答案。一是“为什么用相关矩阵而不是协方差矩阵做 PCA”答量纲不一致花瓣长度方差远大于花萼宽度不标准化会让第一主成分被单一特征主导。二是“GMM 是无监督的你怎么保证成分和类别对应”答判别阶段用马氏距离到各类均值的最近邻规则GMM 只负责估计分布参数两步解耦。三是“刀切法运算量大有没有考虑 k 折”答留一法是刀切法的标准形式150 个样本规模可接受且训练集利用率最高适合小样本场景。一个具体技巧把五组配置的判错率表放在结论页最显眼的位置旁边配一句“LDA 二维在留一法下判错率最低验证了有监督降维在分类任务上的优势”。评委翻到这一页三秒内能抓住你的核心结论后面的追问也会围绕这张表展开你准备起来更有针对性。从那以后我每次带这类分析型作业都强制要求先跑通“原始基线 两种降维 × 两个维度”这五组留一法再动笔写正文。因为只有这张表出来了你才知道哪些结论能写、哪些话是拍脑袋。希望帮到你。本文还有配套的精品资源点击获取