鸢尾花数据集PCA、LDA、GMM与留一法实战

简介：这份数理统计大作业文档面向修读数理统计、模式识别或机器学习入门课程的学生，以经典鸢尾花数据集为主线，串联马氏距离、混合高斯模型、主成分分析、线性判别分析和刀切法等知识点，适合用作课程作业参考、方法复现或复习提纲。压缩包仅含1个docx文件，约511KB，篇幅紧凑但结构完整，从符号说明、样本整理到算法原理、数据处理与结果比对均有覆盖，可对照目录快速定位公式推导与实验结论。其中主成分分析与线性判别分析的降维去噪、高斯混合模型的极大似然聚类、刀切法对降维方式和模型的评估构成核心脉络，读者可借此理解各方法的适用条件、参数估计思路及彼此差异。文档已有1429人学习，兼具课程规范性与案例可读性，便于在短时间内搭建从理论到鸢尾花数据分析的完整框架。

1. 数理统计大作业挑鸢尾花：150 条样本怎么撑起 PCA、LDA、GMM 和刀切法

一份数理统计大作业的题目落到鸢尾花数据集上，很多人的第一反应是 import sklearn、三行 API、截张图交差。但带这门课的导师真正想看的，是把马氏距离、主成分分析、线性判别分析、高斯混合模型、刀切法这几块拼成一条从数据到结论的完整链路。这套资源就是一份完成度很高的参照：一份带公式推导和结果分析的文档，配一套纯 numpy 实现的分析代码，数据是 150 条样本、4 个数值属性、3 个类别的鸢尾花数据集。它适合两类人——正在赶数理统计或模式识别课程大作业、想要一份能跑通又能讲清原理的底稿；以及平时调包调惯了，想把降维、聚类、评估这条链路自己写一遍的从业者。下面按数据、降维、聚类、评估的顺序拆，代码可以直接抄。

2. 数据摸底与马氏距离基线：4 个属性的尺度差和一条 0.94 起跳的线

动手之前先做两件事：把数据读进来核对形状，用一个不看标签的距离基线把上限探出来。很多人跳过这一步直接上 PCA 和 GMM，最后准确率不对也说不清是算法问题还是数据问题。这一节的两个脚本，一个负责把数据摆平，一个负责给出“闭着眼睛也能拿到多少分”的参照线。

2.1 读入、核对与量纲检查

公开版本的鸢尾花数据一般是逗号分隔的文本，前四列是花萼长度、花萼宽度、花瓣长度、花瓣宽度，第五列是类别字符串。看着简单，实际读的时候有三个地方会翻车：文件末尾的空行、类别字符串两端的空格、以及列顺序和你以为的不一样（有些整理版把花萼宽度放在第一列）。所以读入之后必须断言形状和类别计数。

import numpy as np LABEL_MAP = {"Iris-setosa": 0, "Iris-versicolor": 1, "Iris-virginica": 2} def load_iris(path="iris.data"): """返回 X:(150,4) float64, y:(150,) int，标签 0/1/2""" raw = np.genfromtxt(path, delimiter=",", dtype=str) raw = raw[raw[:, 0] != ""] # 末尾空行不剔除，astype(float) 会直接抛 ValueError X = raw[:, :4].astype(float) y = np.array([LABEL_MAP[s.strip()] for s in raw[:, 4]], dtype=int) assert X.shape == (150, 4), X.shape assert np.bincount(y).tolist() == [50, 50, 50], np.bincount(y) return X, y X, y = load_iris() print(np.round(X.mean(0), 3)) # [5.843 3.057 3.758 1.199] print(np.round(X.std(0, ddof=0), 3)) # [0.825 0.434 1.759 0.763]

assert那两行不是摆设。标签映射写错一位，比如把 versicolor 和 virginica 的键写反，后面所有准确率都会掉到 0.6 以下，而 PCA 和 GMM 的中间结果看起来完全正常，这种错误能让人查一晚上。np.bincount是最快的体检手段，三个 50 一出来，标签对齐就确认了。

四个属性的单位都是厘米，但标准差差了四倍：花瓣长度 1.759，花萼宽度只有 0.434。这意味着直接对协方差矩阵做 PCA，第一主成分会被花瓣长度主导，花萼宽度那点信息基本被淹没。文档里选择用相关矩阵代替协方差矩阵，本质就是做了一次 z-score 标准化，把量纲抹平——这一步在数据性质不同的时候是必须的，同一量纲下则看你想不想保留方差作为权重。

提示：如果你在别处拿到的版本均值对不上（比如花萼长度均值不是 5.843），先怀疑列顺序被改过，不要怀疑公式。

2.2 马氏距离：把尺度差和相关性一起吃掉

欧氏距离在两个属性方差差四倍、还彼此相关的时候基本没有意义。马氏距离的形式是 d(x, μ) = sqrt((x-μ)ᵀ Σ⁻¹ (x-μ))，它白化掉了协方差结构，所以尺度无关、也考虑了维度之间的相关性。当 Σ 取单位阵时它退化成欧氏距离，这个特例可以用来验证代码写得对不对。

def mahalanobis(x, mu, Sigma): d = np.asarray(x) - mu return float(np.sqrt(d @ np.linalg.pinv(Sigma) @ d)) mu, S = X.mean(0), np.cov(X, rowvar=False) print(np.round([mahalanobis(x, mu, S) for x in X[:3]], 3)) def loo_maha(X, y): """池化协方差的马氏距离最近类判别，留一法评估""" n, pred = len(y), np.empty(len(y), dtype=int) classes = np.unique(y) for i in range(n): tr = np.arange(n) != i Xtr, ytr = X[tr], y[tr] S_pool = np.cov(Xtr - Xtr.mean(0), rowvar=False) Sinv = np.linalg.pinv(S_pool) # 奇异时 pinv 比 inv 稳 d = [] for c in classes: diff = X[i] - Xtr[ytr == c].mean(0) d.append(diff @ Sinv @ diff) # 比较大小不必开方 pred[i] = classes[int(np.argmin(d))] return (pred == y).mean() print(loo_maha(X, y))

这里有两个刻意的选择。一是用平方马氏距离比较大小，省掉开方，留一法要循环 150 次，能省就省。二是用池化协方差（假设三类共享同一个 Σ）而不是每类各算一个，后者是二次判别，参数多、样本少的时候反而不稳。这条基线跑出来通常在 0.94 以上，如果得到 0.33 或者 0.66，基本可以断定是标签映射或者类别顺序错了，而不是距离公式有问题。

参数上值得动一动的是pinv的容差和协方差的估计方式。样本数 150、维数 4，Σ 不会奇异，但到了降维后的低维空间或者你自己构造的特征上，共线是常事，这时候pinv的默认截断阈值 1e-15 可能太松，手动传rcond=1e-6更稳。另外留一法这里已经是完整形式——每次剔除一个样本，用剩余样本重算均值和协方差——没有把全量数据的统计量泄漏到测试折里，这一点在后面的避坑章节还会重点说。

3. 降维的两条路：PCA 用相关矩阵、LDA 用广义特征值的实现细节

鸢尾花四个属性之间有明显的线性相关，花瓣长度和花瓣宽度几乎是同一条直线上的两个刻度。降维在这个场景下不只是为了画图，它是给后面的高斯混合模型减参数：GMM 的协方差矩阵参数个数随维数平方增长，4 维三成分要估 44 个参数，降到 2 维只要 17 个，150 条样本下这个差别直接决定模型稳不稳。

3.1 PCA：相关矩阵、eigh 和特征向量的符号

PCA 的目标函数是让投影后的方差最大，解就是协方差矩阵（或相关矩阵）最大特征值对应的特征向量。原理一句话说完，实现上真正要小心的是三个点：eigh返回的是升序、用相关矩阵时变换必须同步标准化、特征向量的符号不唯一。

class PCA: def __init__(self, n_components=2, use_corr=True): self.k, self.use_corr = n_components, use_corr def fit(self, X): self.mu_ = X.mean(0) Xc = X - self.mu_ self.sd_ = Xc.std(0, ddof=0) M = np.corrcoef(Xc, rowvar=False) if self.use_corr else np.cov(Xc, rowvar=False) w, V = np.linalg.eigh(M) # eigh 专吃对称阵，返回升序 order = np.argsort(w)[::-1] self.w_, self.V_ = w[order], V[:, order] self.ratio_ = self.w_ / self.w_.sum() # 符号归一：让每个主成分绝对值最大的分量为正，两次运行才可比 for j in range(self.V_.shape[1]): if self.V_[np.argmax(np.abs(self.V_[:, j])), j] < 0: self.V_[:, j] *= -1 return self def transform(self, X): Xs = (X - self.mu_) / (self.sd_ if self.use_corr else 1.0) return Xs @ self.V_[:, :self.k]

transform里那个除法是最容易漏的一步。用相关矩阵求出的特征向量，是在标准化空间里的方向，你拿原始数据直接乘它，等于用一把刻度和坐标系不匹配的尺子去量。表现是：PCA 降维后的分类准确率比不降维还低一截，而且降的维数越少掉得越狠。要么统一用相关矩阵加标准化，要么统一用协方差矩阵加原始数据，两条路选一条走到底。

eigh和eig的区别也值得说一句。协方差矩阵和相关矩阵都是实对称阵，eigh走的是专门算法，结果一定是实数、特征向量正交，还快；eig是通用解，可能返回复数，还得自己取实部。代价是eigh升序排列，np.argsort(w)[::-1]这步不能省，否则你取到的“第一主成分”其实是方差最小的那个方向。

符号归一那段代码看起来多余，实际很有用。特征向量的方向可以整体乘 -1 而不改变任何数学性质，但两次运行（比如留一法里循环 150 次）如果符号随机翻转，你比较投影坐标、画散点图、算成分载荷的时候会得到一堆自相矛盾的结果，很容易误判成代码有 bug。

累计贡献率用np.cumsum(pca.ratio_)看，相关矩阵下前两个主成分一般能到 95% 上下，第三个之后基本是噪声。这也解释了文档里为什么把降维叫“去噪过程”——被丢掉的那部分方差，在这个数据上确实主要是类内抖动。

3.2 LDA：类间类内散度矩阵与广义特征值

LDA 和 PCA 的差别只有一句话：它看标签。目标函数是最大化 uᵀ S_b u / uᵀ S_w u，其中 S_w 是类内散度、S_b 是类间散度，解是 S_w⁻¹ S_b 的特征向量。多出来的信息量让它在同样的维数下通常比 PCA 分得更开，代价是它只能用有标签的数据，而且判别方向最多只有 C-1 个。

def lda_fit(X, y, reg=1e-6): n, p = X.shape mu = X.mean(0) Sw = np.zeros((p, p)) Sb = np.zeros((p, p)) for c in np.unique(y): Xc, muc = X[y == c], X[y == c].mean(0) Sw += (Xc - muc).T @ (Xc - muc) # 类内散度，不除 n_c Sb += len(Xc) * np.outer(muc - mu, muc - mu) # 类间散度按类样本数加权 Sw += reg * np.eye(p) # 近奇异时兜底 w, V = np.linalg.eigh(np.linalg.pinv(Sw) @ Sb) order = np.argsort(w)[::-1] W = V[:, order][:, :len(np.unique(y)) - 1] # 最多 C-1 个方向 return mu, W, w[order] mu, W, eigvals = lda_fit(X, y) Z = (X - mu) @ W print(Z.shape, np.round(eigvals[:4], 3)) # (150, 2) [32.19 0.285 0. 0.]

Sw和Sb的写法里有两个细节。类内散度不除以样本数，因为后面的广义特征值问题里分母的常数会被约掉；类间散度要乘类样本数，因为均值离整体均值越远的类，对类间距离的贡献应该越大，权重就是它的样本量。这份数据三类各 50 条，权重相同，但如果换成不平衡数据，不乘这个权重结果会明显偏。

特征值的分布本身就是个信号。上面打印出来第一个广义特征值三十多，第二个只有 0.28，第三个直接是 0——这正好印证了 LDA 的方向数上限是 C-1=2，第三维在数学上不存在。所以W只能切两列，你写n_components=3也切不出第三列来，eigh会给你一个数值上接近零的垃圾方向。有人在这里加了一维然后发现 GMM 报协方差奇异，根子就在这。

reg那个正则项也不是可选项。S_w 是 4×4 矩阵、150 个样本，正常不会奇异，但留一法循环里每次少一个样本、如果又恰好做了某种筛选，就可能撞上共线。加 1e-6 的对角线对结果几乎无影响，但能避免eigh返回 NaN。用pinv而不是inv是同一个道理。

3.3 维数怎么定：贡献率、C-1 上限和一组对照

降维维数不能拍脑袋。PCA 这边看累计贡献率，一般卡到 95% 或者固定 2 维方便可视化；LDA 这边被 C-1 卡死，最多 2 维。为了让后面的刀切法比较有意义，配置要铺开成一组，而不是只跑一个最优值。

这张表是后面实验矩阵的骨架。PCA-1 和 LDA-1 的对比最能说明问题：同样是一维，LDA 因为用了标签，通常能把其中一类完全分开，而 PCA 的一维投影里三类还会糊在一起。这不是算法谁高谁低，是有没有用到监督信息的差别，写进报告里正好是一段有说服力的分析。

4. 高斯混合模型：EM 的四个参数和“用马氏距离判类”的边界

降维之后要做的是给数据的分布建模。聚类的入门做法是 k-means，它假设每个簇是一个各向同性的球，且硬分配。鸢尾花的三个类在二维平面上是两个细长的椭圆加一个偏圆的簇，k-means 用球去套椭圆，边界必然是歪的。高斯混合模型给每个成分一套完整的协方差矩阵，能表达不同方向上方差不同、维度之间相关，还能给出后验概率，这就是选它的理由。

4.1 参数账：GMM 比 k-means 多花在哪

一个 p 维、K 个成分的 GMM，参数个数是 K × (p 个均值 + p(p+1)/2 个协方差上三角 + 1 个混合权重) − 1。减一是因为权重和为 1 的约束。k-means 只需要 K × p 个中心坐标。这个差距在低维看不出，在高维是灾难。

150 条样本要估 44 个参数，平均每个参数不到 4 个样本支撑，极大似然估计的方差会很大，某个成分的协方差矩阵很容易估成接近奇异。这就是为什么要先降维再做 GMM——不是因为降维“高级”，是因为参数账算不过来。同样这也解释了文档里为什么强调“在极大似然估计意义下求模型参数”：EM 就是这套参数在隐变量存在时的极大似然求解方法。

4.2 EM 实现与三个必须调的参数

EM 的骨架是 E 步算责任度 γ、M 步用 γ 加权重估参数，循环到对数似然不再明显增长。写出能跑的版本不难，写出不会在某次留一循环里突然崩掉的版本，靠的是三个参数：协方差正则、收敛阈值、初始化的方式。

def gmm_logpdf(X, mu, Sigma, pi): """每个成分的对数密度，用 Cholesky 避免 det/inv 的数值问题""" n, p = X.shape K = len(pi) out = np.zeros((n, K)) for k in range(K): L = np.linalg.cholesky(Sigma[k]) sol = np.linalg.solve(L, (X - mu[k]).T) out[:, k] = (-0.5 * np.sum(sol ** 2, 0) - np.sum(np.log(np.diag(L))) - 0.5 * p * np.log(2 * np.pi) + np.log(pi[k])) return out def gmm_em(X, K=3, iters=200, tol=1e-6, reg=1e-6, seed=0): rng = np.random.default_rng(seed) n, p = X.shape mu = X[rng.choice(n, K, replace=False)] # 随机挑 K 个样本点做初值 Sigma = np.array([np.cov(X, rowvar=False) + reg * np.eye(p) for _ in range(K)]) pi = np.full(K, 1.0 / K) prev = -np.inf for it in range(iters): logp = gmm_logpdf(X, mu, Sigma, pi) # E 步，全程在对数域 m = logp.max(1, keepdims=True) ll = m[:, 0] + np.log(np.exp(logp - m).sum(1)) # 逐样本对数似然 gamma = np.exp(logp - ll[:, None]) # 责任度，每行和为 1 Nk = gamma.sum(0) # M 步 mu = (gamma.T @ X) / Nk[:, None] for k in range(K): d = X - mu[k] Sigma[k] = (gamma[:, k][:, None] * d).T @ d / Nk[k] + reg * np.eye(p) pi = Nk / n if ll.mean() - prev < tol: break prev = ll.mean() return mu, Sigma, pi, gamma, ll

reg是保命的。留一法里某次训练折如果让某个成分只分到一两个点，协方差矩阵直接退化成秩 1，cholesky立刻抛LinAlgError。加 1e-6 的对角线相当于给协方差加一个各向同性的先验，结果几乎不变，但矩阵一定正定。这个值不要设太大，1e-3 以上会明显把椭圆压成圆，聚类效果反而下降。

tol判的是平均对数似然的增量，不是参数的变化量。用参数变化做判据的问题是：协方差矩阵的元素量级和均值差好几个数量级，你得给每个参数单独设阈值，非常麻烦。对数似然是标量，增量小于 1e-6 就停，简单可靠。注意ll.mean()每次循环都要重算，别写成累计值。

初始化用随机挑样本点，简单但对seed敏感。更稳的做法是先跑一遍 k-means，用它收敛的中心做初值，EM 通常 20 轮内就稳定。这两种初始化在报告里可以各跑几遍做对比，正好说明“EM 只保证收敛到局部最优”这句话不是空话。

4.3 判决：后验概率还是马氏距离

模型训完，怎么把一个新样本判到某个类？两条路。第一条是直接取责任度最大的成分，等价于最大后验概率，因为它用了混合权重 π。第二条是文档里写的做法：算样本到各成分均值的马氏距离，取最小的那个。

这两种判法在成分权重相等、协方差矩阵相同时完全等价。不等的时候差别就出来了：马氏距离忽略先验 π，会偏向那些协方差小、分布紧凑的成分；后验概率则会考虑“这个成分本来就没几个样本”这件事。三类各 50 条的均衡数据下两者差别很小，但如果你在别的数据上跑，某一类样本只有十几条，忽略先验会让少数类被吞掉。

def predict_by_post(X, mu, Sigma, pi): return np.argmax(gmm_logpdf(X, mu, Sigma, pi), axis=1) def predict_by_maha(X, mu, Sigma): d = np.empty((len(X), len(mu))) for k in range(len(mu)): diff = X - mu[k] d[:, k] = np.einsum("ij,jk,ik->i", diff, np.linalg.pinv(Sigma[k]), diff) return np.argmax(-d, axis=1) # 距离最小即分数最大

用einsum一次算完所有样本到某个均值的平方马氏距离，比 Python 循环快两个数量级，留一法循环 150 次的时候这个差别很实在。另外注意predict_by_maha用的是每个成分自己的协方差矩阵，这与 3.1 节池化协方差的判别是两回事——前者是生成式模型下的类条件密度比较，后者是最近类中心。

注意：用 GMM 做分类和用 GMM 做聚类是两个任务。聚类时 K 是未知的，需要拿 BIC 或者轮廓系数去选；分类时 K 应该等于类别数，并且每个成分要用带标签的样本单独拟合。混着用会得到一个看起来很好看但没有意义的准确率。

5. 避坑清单：留一法结果对不上的五个典型现象

刀切法（留一法）的逻辑很干净：每次留一个样本做测试，其余全部训练，重复 n 次，用判错率评价方法好坏。它比简单的回代法（拿训练集自己测自己）合理得多，也比 k 折交叉验证更充分地利用了样本。但正因为循环次数多、每折都要重新做一遍降维和聚类，任何一个环节写错都会被放大成“这个方法不行”的错觉。下面五条是我在这类作业里见过频率最高的。

5.1 降维与标准化环节的三个坑

现象一：PCA 降维后的准确率低于不降维的原始四维。原因几乎总是标准化不一致——用相关矩阵求了特征向量，transform时却拿原始数据直接乘，或者反过来用协方差矩阵求特征向量、又画蛇添足地除了标准差。这两种情况下投影坐标的尺度完全错位，后面的距离和密度计算全部失真。 解决是把“用相关矩阵”和“投影前除以标准差”绑成一个开关，写成同一个参数控制，代码里不要留两条独立的分支。跑完之后检查一下降维前后样本的类间可分性，如果 PCA-2 的准确率比 raw-4 低超过 5 个百分点，先回来看这一条。

现象二：同样的代码跑两次，LDA 投影出的坐标符号相反。原因是特征向量的符号不唯一，加上 S_w 在留一法里偶尔接近奇异，eigh每次给出的方向在正负之间随机跳。 解决是在lda_fit里做和 PCA 一样的符号归一，同时给 S_w 加 1e-6 的正则。判断有没有归一好的方法很简单：把 150 个样本的投影坐标排序，两次运行的排序结果应该完全一致，只是整体符号可能相反。

现象三：LDA 降到 3 维之后，GMM 报协方差矩阵非正定。LDA 的判别方向上限是 C-1=2，第三列特征向量对应的特征值在数值上是 0 附近的噪声，投出来的那一维基本是常数。一个方差接近零的维度进 GMM，协方差矩阵必然奇异。 解决办法是在lda_fit返回时硬性截断到len(np.unique(y)) - 1列，不要相信调用方传进来的维数。如果确实需要更多维，说明 LDA 这条路走不通，该换回 PCA。

5.2 聚类与评估环节的两个坑

现象四：留一法准确率虚高，接近 100%，但换个数据集就原形毕露。原因是在全量数据上做了标准化、PCA、甚至选了超参数，然后再去跑留一法。测试样本的均值和方差已经参与过训练折的变换，信息泄漏了。150 条样本、4 维、3 类，泄漏带来的虚高能有 3 到 5 个百分点，足以让一个平庸的配置看起来完美。 解决办法是把“fit 变换”这一整段严格关在训练折里面：PCA().fit(Xtr)，然后transform(Xtr)和transform(xte)分别调用，均值和标准差都只能用训练折的。听起来是常识，但sklearn的Pipeline一旦忘了套在外面，这个错误就会悄悄发生。

现象五：EM 的对数似然曲线震荡或者中途下降。理论保证是单调不减，出现下降只有两种可能：E 步在概率域直接算，密度太小导致下溢成 0，log(0)变成-inf；或者 M 步的协方差更新漏了Nk归一化，参数在两步之间来回拉扯。 解决办法是全程在对数域算，用“减去每行最大值再取 exp”这个技巧防下溢，这一点gmm_logpdf和 E 步里都已经写进去了。另外在每次 M 步之后加两个断言，np.allclose(gamma.sum(1), 1)和np.isclose(pi.sum(), 1)，能当场抓住大部分公式笔误。

6. 把留一法跑成实验矩阵：一组配置表和一个可复用脚本

单独跑一个配置说明不了什么，报告里真正有说服力的是“同一套评估协议下，六种配置横向比较”。所以我一般的做法是把降维方法和维数做成配置列表，外面套一层留一法，内层做“降维 + 每类拟合高斯 + 判决”，一次跑完。关键约束只有一个：所有 fit 都必须发生在训练折内部。

def loo_experiment(X, y, configs): n, rows = len(y), [] classes = np.unique(y) for name, red, k, Kc in configs: pred = np.empty(n, dtype=int) for i in range(n): tr = np.arange(n) != i Xtr, ytr, xte = X[tr], y[tr], X[i:i + 1] if red == "pca": m = PCA(n_components=k).fit(Xtr) # 只在训练折上 fit Ztr, zte = m.transform(Xtr), m.transform(xte) elif red == "lda": mu, W, _ = lda_fit(Xtr, ytr) Ztr, zte = (Xtr - mu) @ W, (xte - mu) @ W else: Ztr, zte = Xtr, xte # 生成式分类：每个类单独拟合 Kc 个成分的高斯混合 scores = [] for c in classes: mu_c, S_c, pi_c, _, _ = gmm_em(Ztr[ytr == c], K=Kc, seed=int(c)) scores.append(gmm_logpdf(zte, mu_c, S_c, pi_c).max()) pred[i] = classes[int(np.argmax(scores))] rows.append((name, round(float((pred == y).mean()), 4))) return rows configs = [ ("raw-4", "none", 4, 1), ("PCA-1", "pca", 1, 1), ("PCA-2", "pca", 2, 1), ("PCA-3", "pca", 3, 1), ("LDA-1", "lda", 1, 1), ("LDA-2", "lda", 2, 1), ("LDA-2-GMM2", "lda", 2, 2), ] for name, acc in loo_experiment(X, y, configs): print(f"{name:12s} {acc:.4f}")

配置列表里最后一个Kc=2值得单独说。前面所有配置都是每类一个高斯成分，那其实就是二次判别分析，GMM 在这里没发挥出“多个成分”的能力。把某一类的成分数提到 2，是在假设这个类内部还有子结构——鸢尾花里 versicolor 和 virginica 本身就有重叠，多一个成分有时能把重叠区拆开。这也是这份大作业里 GMM 相对 k-means 的差异能被真正体现出来的地方：k-means 只能给每个簇一个球心，没法在一个类内部表达双峰。

跑完这张表，比较的正确姿势是成对看，不是只看最高分。PCA-2和LDA-2的差反映监督信息值多少；PCA-1和LDA-1的差在同样维数下把这个结论再确认一遍；PCA-3如果比PCA-2低，说明第三主成分确实是噪声，降维起到了滤波作用；raw-4是所有比较的锚点，它要是最高，那这套数据上根本不需要降维，报告里就得诚实写这一点。

计算量上也提醒一句：7 个配置 × 150 折 = 1050 次 EM，每次两百轮上限，纯 Python 大概要跑几分钟。想快就把iters降到 100、tol放到 1e-5，对最终准确率的影响通常在 0.5 个百分点以内，但调试阶段能省掉大量等待。真正出结果的那一遍再调回去，别在没定稿的代码上追求精度。

从那以后我每次做留一法，都会强制在循环体第一行加一句assert not set(tr_idx) & {i}，再用一个np.random.permutation打乱后的标签跑一遍——如果准确率还在 0.9 以上，说明泄漏了，正常应该掉到 0.33 附近。这个自检只要几行代码，能省掉后面所有“为什么我的结果比别人高这么多”的自我怀疑。希望帮到你。

本文还有配套的精品资源，点击获取