简介面向备战数学建模竞赛、尤其是研究2012年A题葡萄酒质量评价的读者这份一等奖论文PDF完整呈现从评分信度检验到葡萄分级、理化指标关联与质量预测的建模思路。压缩包仅含1个PDF文件大小约2.41MB便于直接查阅与打印精读。论文以K-S检验、Wilcoxon符号秩检验和肯德尔和谐系数判断两组评酒员差异与可信度用主成分分析构建葡萄分级体系并以典型相关分析和多元线性回归讨论理化指标对葡萄酒质量的影响同时交代MATLAB、SPSS、SAS与EXCEL的实现路径。该文已有2908人学习适合参赛者借鉴论文结构、模型选择、假设检验表述与结果论证方式也可作为统计建模课程案例参考。1. 一套被反复下载的论文真正值钱的是那条数据处理链路2012 年数学建模 A 题给的是同一批葡萄的两套数据人喝出来的分和仪器测出来的理化值。一等奖论文《葡萄酒的评价》把这两套数据接了起来——单样本 K-S 检验定分布Wilcoxon 符号秩检验判两组评酒员有没有系统性差异肯德尔和谐系数比谁更可信主成分加 AHP 把二十多个葡萄指标压成一个综合分做分级最后用典型相关和多元回归回答「理化指标能不能替掉品酒员」。它不是一篇读完就完的范文而是一条可以照着跑的数据流水线三份附件、四个问题、五种统计模型每一步都有明确的输入、输出和判定阈值。适合正在做统计建模课设的人也适合手里攥着一批评分数据、想搞清楚评分者信度到底怎么算的从业者。2. 三份附件先过一遍字段口径、缺失填补与异常修正动手之前得先把数据摸清楚。这套题最容易劝退人的地方不是模型是附件 1 到附件 3 的字段结构完全不是一个形态一个是人打分两个是仪器测值一个按评酒员排一个按样本排。不先把它们对齐成「样本 × 指标」的矩阵后面任何一行代码都跑不通。2.1 附件 1/2/3 的字段结构和量纲差异三份附件的组织方式差别很大先把结构列清楚再谈建模附件观测主体主要字段规模取值与量纲附件 1评酒员 × 酒样外观澄清度、色调、香气纯正度、浓度、质量、口感纯正度、浓度、持久性、质量、整体评价红葡萄酒 27 个样本、白葡萄酒 28 个样本两组各 10 名评酒员整数等级分各分项上限不一致附件 2酿酒葡萄样本氨基酸总量、蛋白质、VC 含量、花色苷、酒石酸、苹果酸、柠檬酸、多酚氧化酶活力、褐变度、DPPH 自由基 1/IC50、总酚、单宁、总糖、可溶性固形物、干物质、果皮质量、果穗质量、色差 a*/b*、可滴定酸、固酸比、黄酮醇等红葡萄 27 行单位混杂g/L、百分比、色差为无量纲附件 3葡萄酒样本花色苷、单宁、总酚、DPPH 半抑制体积、酒精度、pH 等与附件 1 一一对应部分为体积分数部分为浓度这张表要盯住两点。第一附件 1 的行是评酒员、列是样本而附件 2/3 的行就是样本转置和对齐是第一步第二附件 2 的指标量纲跨度极大总糖在 g/L 量级色差 a* 在个位数量级直接丢进主成分方差全被总糖吃掉标准化不是可选项而是必做项。常见做法是先写一张字段映射表把「样本编号」作为唯一主键三份数据全部左连接到这张主键上。别小看这一步红葡萄 27 个样本、白葡萄 28 个样本编号只要有两位数字不一致01 和 1后面的配对检验就会静默出错——不报错但结果全歪。提示白葡萄酒比红葡萄酒多一个样本所有按样本循环的代码都要写成按组分别处理别用一个写死的 27。2.2 均值替换的两个前提组内取均值、结果取整原文在附件 1 里发现了两类脏数据处理方式不一样。第一类是缺失红葡萄酒样品 20 的色调分4 号评酒员那一格是空的其余 9 人打的是 6、6、4、6、6、8、6、6、8剔掉后取均值再取整填 6。第二类是异常第一组白葡萄酒样品 3 的持久性分7 号评酒员打了 77而相邻评酒员都在 5 到 7 之间明显是手工录入多敲了一个 7。填补代码本身很短关键在于「按什么范围取均值」import pandas as pd import numpy as np def fill_by_item_mean(df, as_intTrue): df: 行评酒员列样本值某一分项的评分 只在该样本的其余评酒员范围内取均值不做全表平均 filled df.copy() for col in filled.columns: s filled[col] if s.isna().any(): m s.mean(skipnaTrue) # 同一列 同一样本 filled[col] s.fillna(round(m) if as_int else m) return filled def flag_out_of_range(df, low, high): 评分分项有明确上下界越界直接判异常返回待人工复核的位置 mask (df low) | (df high) return df.where(mask) # 非异常位置为 NaN方便定位fill_by_item_mean里最要紧的是s.mean()的作用域。这里的s是某一分项下所有样本的列如果数据结构是「行样本、列评酒员」那s.mean()就是跨评酒员取均值方向对了如果行列反了就变成跨样本取均值等于把不同酒样的水平混在一起方差被抹平后面 Wilcoxon 的配对差值分布会跟着变。as_intTrue不是洁癖。评酒员打的是整数等级分如果插出一个 6.4Wilcoxon 编秩时这一格会和别的样本产生额外的并列秩秩和统计量变化z 值可能和原文差零点几。原文报的红葡萄酒 z 约 -2.53、白葡萄酒 z 约 -2.23都在 1.96 的临界线外侧差一点点不影响结论方向但换一批数据就未必了。2.3 异常值走「先剔除后替换」别用 3σ 一把梭原文对 77 这种数据的处理是「先剔除后替换」先置为缺失再走 2.2 的均值填补。这个顺序比直接改数要稳因为它保留了「这里原来有问题」的痕迹。我一般不会对评分数据用 3σ 准则。原因很直接评分样本量小一个分项就 10 个评酒员标准差本身估不准3σ 的边界会随样本漂而且评分有硬上下界用领域规则比统计规则可靠得多。分项的上限是已知的超上限就是错低于下限也是错这种判断不需要统计。真正需要小心的是「看起来不合理但其实合理」的数据。比如某位评酒员对某个样品的整体评价给了很低的分而同组其他人给得高。这不一定是异常可能是这位评酒员对这款酒真的不买账——评酒员之间的评价尺度差异本来就是这套题要研究的东西把它当异常清掉等于人为制造了一致性。判断标准可以这样定越界超过分项上限、低于 0、出现两位数的判异常组内偏离中位数超过该分项极差一半的标出来人工看一眼别自动改。2.4 把评分表压成样本级配对矩阵问题一要比较「两组评酒员的评价结果」比较对象是样本级的。所以中间要做一次聚合10 名评酒员的分项分先求和成每人的总分再取组内均值作为该样本在该组的代表值。item_cols [澄清度, 色调, 纯正度_香, 浓度_香, 质量_香, 纯正度_口, 浓度_口, 持久性, 质量_口, 整体评价] # 1) 每人总分等权求和和原文口径一致 score raw.groupby([组别, 样本, 评酒员])[item_cols].sum().reset_index() score score.rename(columns{0: 总分}) if 总分 not in score else score # 2) 样本级每组每个样本一个代表值红白分开 red score.query(酒种 红).pivot(index样本, columns组别, values总分).mean(axis1)第一段是分项等权求和。原文的口径就是这样虽然香气质量、口感质量这些分项的分值上限本来就比澄清度高等权求和实际上已经隐含了权重。如果自己再乘一套权重问题一的 z 值会变但方向通常不变。第二段容易写错。pivot之后如果直接把 10 名评酒员的行铺开得到的是 10 行 × 27 列的宽表Wilcoxon 需要的是「同一批样本、两组各一个值」的配对序列长度必须是样本数。所以这里先用groupby把评酒员维度压掉再 pivot 出两组最后两列一一对应。注意压掉评酒员维度会丢失组内离散程度。问题一如果还想讨论「组内分歧」得另外保留标准差别在这一步扔掉。3. 评分差异与葡萄分级Wilcoxon 编秩到 AHP 赋权的完整链路问题一和问题二看着是两件事实际上是一条链先用假设检验确认两组评分不可混用选出一组相对可信的分数作为「质量真值」再把这个真值当作评价体系里的一级指标和葡萄的理化指标一起做综合评价最后落到分级。链断在哪一环后面的结论都站不住。3.1 单样本 K-S 检验p0.525 说明的是「不拒绝」原文对四组评分数据各做了一遍正态性检验用的是 P-P 图和单样本 K-S第一组红葡萄酒得到 p0.525大于 0.05于是认为评分服从正态分布进而说明评分数据是科学的、合理的。这个推理链条要打折看。K-S 检验的原假设是「样本来自指定分布」p 大于 0.05 只能说不拒绝原假设不等于证明了正态。而且这里的样本是 27 个样本均值样本量小检验功效低本来就很难拒绝。from scipy import stats x red[第一组].values # 27 个样本的组内均分 z (x - x.mean()) / x.std(ddof1) d, p stats.kstest(z, norm) # 参数由样本估计严格说要用 Lilliefors 修正 print(fD{d:.4f}, p{p:.4f}, n{len(x)})kstest这里传的是标准化后的数据和标准正态等价于检验正态性但临界值表是按「分布参数已知」算的。参数是从样本估出来的p 值会偏大偏保守专业的做法是用 Lilliefors 修正或者直接上 Shapiro-Wilk。原文用 SPSS 的 K-S 选项同样存在这个偏差。不过这一步在整条链路里不是关键。问题一后面用的是 Wilcoxon 符号秩检验本身就是非参数方法不要求正态。K-S 在这里更多是交代「为什么选这个模型」的叙事作用。复现时对不上 p 值不用纠结对上了反而是巧合。3.2 Wilcoxon 符号秩检验配对、编秩与 z 值的算法两组评酒员评的是同一批酒样这是天然的配对数据所以用两配对样本的非参数检验。原文的步骤写得很规范算成对差值 D、按绝对值编秩、恢复符号、分别求正秩和 T 与负秩和 T−、取较小者为统计量 T、再用正态近似算 z。以红葡萄酒为例n27E(T)n(n1)/4189D(T)n(n1)(2n1)/24≈1732.5标准差约 41.6原文报的 z≈-2.53反推 T 在 84 左右落在拒绝域里。白葡萄酒 z≈-2.23同样拒绝。结论是两组对红、白葡萄酒的评价都存在显著性差异。Python 里一行就能跑from scipy.stats import wilcoxon stat, p wilcoxon(red[第一组], red[第二组], zero_methodwilcox, # 差值为 0 的样本如何处理 correctionFalse, # 是否做连续性校正 alternativetwo-sided) print(fT{stat}, p{p:.4f})三个参数值得逐个说。zero_method决定差值为 0 的配对怎么算wilcox直接丢弃、pratt保留进秩、zsplit把零差值平分到正负两侧。27 个样本里出现一两个零差值不稀奇选哪个会改变有效 n进而改变 z。原文没有明确交代这一点复现时先统计一下零差值个数写进方法说明里。correction是连续性校正样本量小的时候打开会让 z 往 0 缩一点更保守。27 对数据不算大我一般会开着跑一遍和不开的结果对照看结论是否翻转。alternative用双侧因为原文的备择假设就是「两组评价结果不同」没预设方向。还要提醒一句wilcoxon返回的stat是较小的秩和不是 z。想要和原文的 -2.53 对上得自己按 n 换算z (T - n*(n1)/4) / sqrt(n*(n1)*(2*n1)/24)。第一次复现时对不上数很正常不是代码写错了。3.3 肯德尔和谐系数 W四组信度比较与卡方边界差异检验只能说明两组不一样回答不了「哪组更可信」。原文换了个角度把评分转成秩看组内 10 名评酒员对 27 个样本的排序有多一致。一致性强说明这组评酒员的评价尺度稳定可信度更高。公式是 W 12Q / [m²(n³−n)]其中 Q 是各样本秩和与平均秩和之差的平方和m 是评酒员数10n 是样本数27 或 28。检验统计量 χ² m(n−1)W。import numpy as np from scipy.stats import rankdata def kendall_w(mat): mat: 行评酒员列样本值评分 返回和谐系数 W 与卡方统计量 m, n mat.shape R rankdata(mat, axis0).sum(axis0) # 每个样本的秩和 Q ((R - R.mean()) ** 2).sum() W 12 * Q / (m ** 2 * (n ** 3 - n)) return W, m * (n - 1) * Wrankdata(mat, axis0)是按列样本在评酒员之间编秩方向别搞反。默认methodaverage遇到并列分数取平均秩和原文的做法一致。原文算出的四组 W 分别是 0.044366788、0.043314652、0.040039039、0.041112055对应 χ² 为 11.535、11.262、10.811、11.100。红葡萄酒第一组的 W 和 χ² 都略高白葡萄酒第二组略高于是得出「红葡萄酒第一组更可信、白葡萄酒第二组更可信」并在后续分析中分别采用这两组分数。这里有一个必须说清楚的边界。χ²11.535自由度按 n−126 算0.05 水平的临界值约 38.9按 27 算也在 40.1 左右。四组数据的 χ² 都远没有越过这道线也就是说没有一组达到「评分者显著一致」的统计门槛。W0.044 这个量级本身就说明 10 名评酒员之间的排序一致性很弱。原文的用法是横向比大小——四个值里谁高谁可信。这种相对排序在方案比较里是成立的但不能顺着写成「评分具有显著一致性」。做工程判断时这反而是有用信息两组评酒员的内部分歧都很大「谁更可信」只是矮子里拔将军。复现时把这句话改准能省掉答辩时一轮追问。3.4 主成分降维27 个理化指标压到 8 个的阈值问题二要拿葡萄的理化指标做分级但指标太多且彼此相关——总酚和单宁相关、总糖和可溶性固形物相关直接加权会重复计权、放大某类信息。原文用主成分分析先把 27 个指标压成 8 个互不相关的主成分累计方差贡献率 83.044%。from sklearn.preprocessing import StandardScaler from sklearn.decomposition import PCA X StandardScaler().fit_transform(grapes) # 行样本列理化指标 pca PCA(n_components0.83044, svd_solverfull) Z pca.fit_transform(X) print(pca.n_components_, pca.explained_variance_ratio_.cumsum()[-1])n_components传浮点表示「累计方差比例达到这个值就停」传整数表示固定个数。原文是按 83.044% 截的正好落在第 8 个主成分上所以浮点阈值能对上如果传 8 也等价但换一批数据就不一定了浮点写法更通用。StandardScaler必须放在 PCA 前面。附件 2 里总糖是 g/L 量级、色差是无量纲个位数不标准化的话第一主成分基本就是总糖其他指标全被压成噪声。主成分解读还有个坑。原文报的是旋转后的载荷方差最大化旋转所以能说出「主成分 1 是总黄酮、总酚、DPPH 自由基和单宁的组合」「主成分 2 是总糖、可溶性固形物和干物质的组合」这种话。sklearn 默认不旋转载荷矩阵和 SPSS 的输出长得不一样同一个主成分对应的指标组合可能变。复现时要么补一步 varimax 旋转要么只报载荷绝对值排前几位的指标别硬套原文的主成分命名。3.5 AHP 判断矩阵与阶梯分级CR0.029 之后还剩什么指标压完了接下来要定权重。原文用层次分析法构造了一个 8×8 的判断矩阵最大特征根为 8CI≈0.013CR≈0.029小于 0.1通过一致性检验。归一化后的权重分布在 0.0253 到 0.0987本文还有配套的精品资源点击获取
