2022数模C题玻璃成分分析实战：成分数据与小样本建模避坑指南

简介：面向数学建模竞赛参赛者和数据分析学习者，提供2022年全国大学生数学建模C题完整解题方案。内容围绕古代玻璃文物成分分析与鉴别，系统涵盖数据预处理与探索、特征降维、因子分析构建风化前后成分转换、亚类划分、灰色关联度分析及敏感性检验，并配有对应实现代码，可直接对照论文理解模型思路与计算细节。资源为单份PDF文件，容量3.55MB，兼具论文阐述与代码展示，便于按章节研读和复现。目前已有927人浏览学习，适合正在备战国赛或练习分类预测与降维建模的选手作为范例参考。

1. 2022全国大学生数学建模C题在做什么：一道“成分数据”+“小样本”的复合题

2022年全国大学生数学建模竞赛C题，赛题名称是《古代玻璃制品的成分分析与鉴别》。表面上是考古学场景，实际是一个典型的小样本成分数据挖掘问题：给你一张几十行的玻璃成分表（SiO₂、PbO、BaO等十几种氧化物的百分比），先按类型分成铅钡玻璃和高钾玻璃，再在每类里找亚类，最后还要预测“风化前”的化学成分。拿奖的关键不在用多深的模型，而在数据清洗、成分数据变换和防数据泄漏这三件事做到位。适合谁？想冲国赛奖的建模组，以及需要快速上手“化学组分表格”这类结构化数据的同学。

2. 拿到C题数据先做三件事：表对齐、缺失值处理和成分特征构造

2.1 附件1和附件2怎么关联：以“采样编号”为主键的唯一合并

C题的数据由两部分组成。一部分是成分实测表，记录了每个采样点的氧化物百分比和“是否风化”标记；另一部分是样品外观信息表，包含颜色、纹饰等字段。两表之间靠“采样编号”关联。拿到数据的第一步不是建模，而是确认两个表里有多少编号能对齐、有多少孤儿行。常见翻车是直接用pandas默认的inner join，合并完才发现样本数少了一截，后面的聚类和回归全在残缺数据上跑。

import pandas as pd # 附件1: 成分实测, 附件2: 颜色/纹饰等外观信息 comp = pd.read_excel('附件1.xlsx') info = pd.read_excel('附件2.xlsx') # 先看两个表各自的形状和主键 print("成分表:", comp.shape, comp['采样编号'].nunique()) print("信息表:", info.shape, info['采样编号'].nunique()) # 左连接保留成分表的全部样本, 检查哪些编号没匹配上 df = comp.merge(info, on='采样编号', how='left') not_matched = df[df['颜色'].isna()]['采样编号'].tolist() print("未匹配到附件2的编号:", not_matched)

这里用how='left'而不是how='inner'，目的是暴露对齐问题而不是隐藏它。左连接后如果发现某行颜色全为空，就该回头查附件1里是否存在重复编号或不可见字符。两个表的主键都是“采样编号”，但实际数据里出现过编号带空格或全角半角混用导致merge失败的情况，排查时先对编号做strip()和str.upper()统一样式。对齐确认后，再进入清洗流程。

2.2 缺失值不是闹着玩的：成分数据缺失的两种处理姿势

成分表的缺失一般分两种。一种是某个成分在所有样本中都没测到，这种直接删列；另一种是部分样本缺失，常见原因是检测仪器能力下限或样品量不够。千万不要对所有缺失列一律填0，因为化学成分百分比在后续计算几何平均时会直接导致对数变换失败。我自己打比赛时的习惯是：先看缺失率矩阵，再做逐列处理。

chem_cols = ['SiO2', 'Na2O', 'K2O', 'CaO', 'MgO', 'Al2O3', 'Fe2O3', 'CuO', 'PbO', 'BaO', 'P2O5', 'SrO', 'SnO2', 'SO2'] # 缺失率统计 miss = df[chem_cols].isna().mean().sort_values(ascending=False) print(miss[miss > 0]) # 缺失率超过50%的列直接删, 低于20%的用该类型(高钾/铅钡)组内中位数填充 drop_cols = miss[miss > 0.5].index.tolist() fill_cols = miss[(miss > 0) & (miss <= 0.5)].index.tolist() print("删列:", drop_cols, "填充列:", fill_cols) for c in fill_cols: df[c] = df.groupby('类型')[c].transform(lambda s: s.fillna(s.median()))

为什么不全局填中位数？因为铅钡玻璃的PbO动辄百分之三四十，高钾玻璃的PbO可能只有个位数，全局填值会严重扭曲两类样本的判别特征。按“类型”分组填充，保留的是类别内部的典型水平。如果某列缺失率在20%~50%之间，填充后还要在论文的“数据处理”一节明确交代填充规则，评委给分时会看这个细节。对个别样本多个成分同时缺失的情况，另一个做法是直接删行，但C题样本量本来就小，删行要谨慎，优先保留关键成分完整的行。

2.3 百分比成分的“伪相关”：为什么不做标准化而是做CLR变换

成分数据有一个让新手初赛翻车的特性：每行的所有氧化物百分比加起来接近100%，这是一个“闭合”约束，导致变量之间存在天然的负相关。直接对成分列做z-score标准化只能处理量纲，不能消除闭合效应。两个本就可比的变量在PCA、聚类和回归里会产生假信号。常规解法是CLR（中心化对数比）变换，这是处理成分数据的标准动作。

import numpy as np def clr_transform(df_in, cols): """对指定的成分列做中心化对数比变换。 先把缺失和0值替换为一个小正数(0.001), 再对每行除以几何平均后取对数。 """ X = df_in[cols].copy() X = X.replace(0, np.nan).fillna(0.001) gm = np.exp(X.apply(np.log).mean(axis=1)) # 每行几何平均 clr = X.apply(np.log).sub(np.log(gm), axis=0) return clr X_clr = clr_transform(df, chem_cols) print(X_clr.head())

CLR把闭合的成分矩阵变成不受和约束影响的欧氏空间数据。变换后的每个值是“该成分相对该样本所有成分平均水平的偏离”，正负号天然可解释：正代表该样本里此成分偏高。对C题来说，后续无论是聚类找亚类、LDA判别还是XGBoost建模，都在CLR空间里做。唯一需要说明的是，CLR会引入一行中所有列的线性相关（每行之和为0），这会让普通协方差矩阵不可逆，所以搭配正则化方法（Ridge、带收缩的LDA），不要用普通线性回归。

提示：如果某行缺失值太多（比如超过一半成分缺失），填充后CLR会失真，这种样本建议直接删行，并在论文里写明删了几条、为什么删。

3. 分类与亚类划分：先聚类再判别的两步法

3.1 用KMeans找亚类：K值怎么定

C题第一问要求在铅钡玻璃和高钾玻璃内部找出亚类。这里的“亚类”没有标签，只能是聚类。常见做法是分别取两类样本做KMeans，用轮廓系数和手肘图选K。样本量小，K一般取2~4，不要贪多。我一般会固定随机种子跑，否则每次结果不一样，论文里没法交代。

from sklearn.cluster import KMeans from sklearn.metrics import silhouette_score def kmeans_k_search(X, max_k=5): scores = {} for k in range(2, max_k + 1): km = KMeans(n_clusters=k, random_state=42, n_init=10) labels = km.fit_predict(X) scores[k] = silhouette_score(X, labels) return scores # 注意: X必须是CLR变换后的矩阵 pb_clr = X_clr[df['类型'] == '铅钡玻璃'] kj_clr = X_clr[df['类型'] == '高钾玻璃'] print("铅钡玻璃轮廓系数:", kmeans_k_search(pb_clr)) print("高钾玻璃轮廓系数:", kmeans_k_search(kj_clr))

轮廓系数越大说明类内越紧凑、类间越分离。但小样本下轮廓系数会偏高且不平滑，我的做法是把轮廓系数和业务常识结合：比如铅钡玻璃按PbO和BaO的高低分离出一类“高铅钡”、一类“低铅钡”，这比纯统计选K更容易在论文里讲圆。KMeans的n_init=10是克服小样本下初始中心随机导致的波动，random_state=42保证复核时结果一致。聚类到的亚类标签，要合并回原数据表，作为后面的判别模型的Y。

3.2 亚类判别：LDA与XGBoost的小样本参数

拿到亚类标签后，第一问还要验证亚类之间是否真的可分。这时候可以用LDA（线性判别分析）做一个基线判别模型。LDA在样本量小到几十行的场景下非常稳，用它算出来的交叉验证准确率，本身就是“亚类差异显著”的直接证据。

from sklearn.discriminant_analysis import LinearDiscriminantAnalysis from sklearn.model_selection import cross_val_score, StratifiedKFold lda = LinearDiscriminantAnalysis(solver='lsqr', shrinkage='auto') cv = StratifiedKFold(n_splits=5, shuffle=True, random_state=42) scores = cross_val_score(lda, X_clr, y_subclass, cv=cv, scoring='accuracy') print("LDA 5折平均准确率: %.3f (+/- %.3f)" % (scores.mean(), scores.std()))

LDA用lsqr+shrinkage='auto'而不是默认的svd，是因为CLR变换后特征矩阵每行和为0，协方差矩阵奇异，收缩估计能稳定求解。如果LDA准确率已经到90%以上，后续用XGBoost只是锦上添花。XGBoost在小样本上的优势是能捕捉非线性交互，劣势是极易过拟合，需要把树压浅。

import xgboost as xgb xgb_model = xgb.XGBClassifier( n_estimators=120, max_depth=2, learning_rate=0.05, subsample=0.8, colsample_bytree=0.8, reg_lambda=1.0, eval_metric='mlogloss', random_state=42 ) cv_xgb = cross_val_score(xgb_model, X_clr, y_subclass, cv=cv, scoring='accuracy') print("XGBoost 5折平均准确率: %.3f" % cv_xgb.mean())

max_depth=2是关键，小样本下每棵树的深度超过3就一定过拟合。learning_rate=0.05配合n_estimators=120，相当于用更多的浅树做集成。colsample_bytree=0.8每次只用80%的特征，给特征选择留随机性。这些参数不是玄学，而是针对“几十行、十几列”的配置习惯；换成几千行的大数据集，同样的参数反而欠拟合。调参时不要开网格搜索硬扫，小样本下网格搜索容易选到对随机种子敏感的参数组合。

3.3 交叉验证切法：GroupKFold避免同源样本泄漏

C题数据里有“对象ID”这类比“采样编号”更高一级的分组字段，一件古代玻璃器物可能被采样检测多次，生成多行数据。如果按普通分层抽样把同一器物的不同采样拆到训练集和验证集，模型会靠记住这个器物来刷分，跨器物的泛化能力其实很差。这里要用GroupKFold。

from sklearn.model_selection import GroupKFold group_cv = GroupKFold(n_splits=5) xgb_model.set_params(eval_metric='mlogloss') gkf_scores = cross_val_score( xgb_model, X_clr, y_subclass, cv=group_cv, groups=df['对象ID'] ) print("GroupKFold准确率: %.3f" % gkf_scores.mean())

对比普通K折和GroupKFold的结果，你会发现后面这个准确率更低，但这才反映真实赛题场景。评委不会单独跑你的交叉验证，但如果你在论文里写“按对象ID分组防止同源泄漏”，这是明显的加分项。如果Excel里没有“对象ID”列，退而求其次用“采样编号”当group，但同一器物多次采样可能编号不同，这种情况要结合题目文字说明来判断分组字段，不要硬造。

4. 风化前成分预测：多输出回归与误差评估

4.1 问题拆解：哪些列是特征、哪些列是目标

C题后半部分围绕“风化”展开。风化点的成分是风化后检测值，想恢复它风化前的成分；未风化点的成分代表风化前的基准。做法是以未风化样本为训练集，以类型、颜色、纹饰作为特征，以各成分为多输出目标。另外，可以把风化不敏感的成分（如Al2O3、SiO2）也加入特征矩阵，因为这些成分在风化前后变化小，能充当样本的“身份信息”。

# 构造特征矩阵: 类型/颜色/纹饰 one-hot feat = pd.get_dummies(df[['类型', '颜色', '纹饰']], prefix='', prefix_sep='') # 目标列: 风化影响明显的成分 target_cols = ['K2O', 'PbO', 'BaO', 'CuO', 'P2O5'] train_mask = df['是否风化'] == '未风化' X_train_all, X_pred_all = feat[train_mask], feat[~train_mask] y_train_all = df.loc[train_mask, target_cols]

这样拆的前提是“是否风化”字段明确给出。实际操作时，部分样本只有风化点没有未风化点成对数据，所以多输出回归是唯一出路。如果有“同一器物风化前后都测了”的配对数据，那可以退化成逐样本差值回归，但C题没有这种配对，别硬找。特征矩阵里不要混入成分列，否则相当于告诉模型答案。

4.2 三套回归方案对比：多元线性、Ridge、随机森林

多输出回归有三种常用方案，复杂度递增：独立多输出线性回归、Ridge多输出、MultiOutputRegressor包随机森林。我的习惯是先跑Ridge，再加随机森林做对照，最后在论文里报告哪个更稳。

from sklearn.linear_model import Ridge from sklearn.ensemble import RandomForestRegressor from sklearn.multioutput import MultiOutputRegressor from sklearn.metrics import mean_absolute_error, mean_squared_error from sklearn.model_selection import cross_val_predict # 方案1: Ridge, 带L2正则 ridge = Ridge(alpha=10.0) ridge.fit(X_train_all, y_train_all) y_pred_ridge = ridge.predict(X_pred_all) # 方案2: 随机森林多输出 rf = MultiOutputRegressor( RandomForestRegressor(n_estimators=200, max_depth=3, random_state=42) ) rf.fit(X_train_all, y_train_all) y_pred_rf = rf.predict(X_pred_all) # 评估: 用未风化样本做交叉验证, 看整体MAE和每个成分的MAE cv_pred = cross_val_predict(ridge, X_train_all, y_train_all, cv=5) mae_ridge = mean_absolute_error(y_train_all, cv_pred) mae_per_col = (y_train_all - cv_pred).abs().mean(axis=0) print("Ridge整体MAE:", mae_ridge) print("各成分MAE:", mae_per_col)

alpha=10.0看起来很大，但在CLR空间里各列经过对数化，正则强度要比默认的1.0大一个量级才能压住方差。RandomForestRegressor(max_depth=3)同样是为了限制单棵树复杂度，防止对小训练集的记忆。交叉验证里报告整体MAE远远不够，按成分列看误差能定位哪个成分最难预测——比如PbO的MAE可能很大，因为它在铅钡玻璃中的绝对值本来就高。写论文时，每个成分的MAE要分开列，评委才信你做过误差分析。

4.3 预测结果怎么写成论文里的统计描述

回归模型的输出只是中间产物，C题最后要落到“某类型风化前的化学成分均值与标准差”。这是论文写作的关键一步：预测出每个风化样本的原始成分后，按亚类合并统计。这里先把CLR预测值逆变换回原始百分比，再分组聚合。

# 逆变换: 把CLR预测值还原为成分百分比 log_gm = np.log((X_clr.iloc[train_mask].apply(np.exp)).mean(axis=1)) # 训练集每行log几何均值 pred_clr = ridge.predict(X_pred_all) # 示例用Ridge pred_orig = np.exp(pred_clr + log_gm.values[:, None]) pred_df = df[~train_mask].copy() pred_df[target_cols] = pred_orig summary = pred_df.groupby(['类型', '亚类'])[target_cols].agg(['mean', 'std']) print(summary.round(2))

按“类型+亚类”算出的均值表，可以直接作为论文中“风化前成分规律”的支撑数据。做这步时注意逆变换细节：CLR只对训练时用过的列有效，逆变换时要保留当时的log_gm向量，不然数值对不上。有参赛队把CLR的预测结果直接当成分百分比写进论文，均值表变成对数尺度，评委一眼就能看出来。

表格里每个数字都要能追溯到代码输出的原始统计，不要手填。这个“均值±标准差”表就是论文结论的核心，后面所有文字分析都得围绕它展开。

5. C题避坑清单：5个让参赛队翻车的细节

5.1 坑1：附件1和附件2的ID对不上

现象：merge之后发现大量行没有颜色和纹饰，或者形状从80行变成60行。原因：两个附件里的“采样编号”字符编码不一致，比如一个全角一个半角，或者末尾带不可见空格。解决：合并前给编号列做astype(str).str.strip()，再用pd.unique对比两边的编号集合，把差异编号列出来人工核对。用inner join默默丢弃数据是最大忌讳，要把对齐过程写进论文附录。

5.2 坑2：全样本准确率100%——数据泄漏

现象：亚类判别模型在训练集上准确率100%，交叉验证也接近满分，但论文里写不出有意义的规律。原因：把聚类用的特征和判别用的特征用了同一份，或者把“类型”字段本身当成聚类的输入。类型是给定的标签，聚类应该在“类型内部”做；如果整表直接聚类，聚类结果会复制“高钾/铅钡”二分，而不是亚类。解决：严格按“先按类型分组→组内CLR→组内聚类→组内判别”的流程，判别模型的特征里不许出现用于定义亚类的同一批成分，或者用嵌套交叉验证放行。

5.3 坑3：成分百分比总和不是100%

现象：某几行成分加总只有92%或者104%。原因：原始检测数据存在误差，题目附件没有规范化到100%。如果数据处理时强行按行归一化到100%，会破坏CLR的几何均值和后续回归误差。解决：保留原始百分数，只做CLR；在论文“数据预处理”里注明“部分样本成分总和偏离100%，可能包含未检测出的微量组分，故不强制归一化”。这比假装数据完美要诚实得多。

5.4 坑4：把风化状态直接当标签，忽略了“亚类”才是第一问

现象：第一问写成了“用分类模型区分风化/未风化”，跑了半天发现准确率95%，但评阅标准里要求的是亚类划分。原因：没读懂问题。风化状态不是亚类标签，它是样本的一个属性；第一问要求的是在铅钡玻璃和高钾玻璃内部找亚类，并分析亚类差异。解决：动笔前先列一张“问题→数据集→模型输出”的映射表，第一问输出亚类标签，第二问输出亚类判别规则，第三问输出风化前成分。写完论文再检查一遍这张映射，比调参重要得多。

5.5 坑5：回归预测后成分归一化

现象：预测出的风化前成分加总偏离100%，于是有人又做了一次按行归一化。原因：多输出回归是逐成分独立预测的，没有显式约束加和为100%，预测值的总和天然会漂。解决：不要在回归后强制归一化。加和漂移本身说明模型误差，论文里报告MAE和平均偏差即可。如果非要保证加和为100%，应该在逆变换后用等距对数比变换（ILR）或加一个成分总和误差项做微调，而不是简单除以总和。

6. 提交前代码自查清单：评委复现你的结果只需三步

代码不是论文的附赠品，在国赛评审里代码质量会影响奖项走向。我给一份自己每次提交前都会过的自查清单，按顺序执行，五分钟内能完成。

第一，确认代码能从“读Excel”到“输出结果表”一步跑通，没有隐藏的绝对路径依赖。把pd.read_excel('/Users/xxx/...')改成pd.read_excel('附件1.xlsx')，文件放在和代码同一目录。第二，确认所有随机过程固定了种子。KMeans的random_state、XGBoost的random_state、交叉验证的random_state，一个都不能少。第三，确认输出文件命名规范：建议result1.csv、result2.csv、result3.csv，对应题目三个问的答案表。评委没有时间在你的脚本里翻哪个变量对应哪个问题。

tree # 建议的提交目录结构 . ├── 代码/ │ ├── main.py │ └── utils.py ├── 数据/ │ ├── 附件1.xlsx │ └── 附件2.xlsx ├── 结果/ │ ├── result1.csv │ ├── result2.csv │ └── result3.csv └── 论文.pdf

这个树状结构不是摆设，评委会按这个路径找你的结果文件。目录里不要放训练过程的中间产物或调试脚本，保持干净。最后，在代码开头写一段10行以内的中文注释，说明“运行环境：Python 3.9 + pandas + xgboost + scikit-learn”，并交代主体流程。这段注释是给评委看的，也是给三天后重新打开代码的你自己的后悔药。

做完这场题，我最深的一个教训是：数学建模竞赛比的不是模型多前沿，而是对数据的诚实程度。CLR变换、GroupKFold、按类型分组填缺失值，这些都不炫，但它们决定了论文能不能被复现、结论能不能站住。希望这套以2022 C题为入口的完整处理链能帮到你，下次再碰上成分数据或者小样本分类题，照着这个框架走一遍，至少不会在坑里过夜。

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