1. 项目概述当数学建模遇见文化遗产去年带队参加国赛拿到C题《古代玻璃制品的成分分析与鉴别》时第一感觉是题目出得很有意思。它把数学建模这个看似“硬核”的竞赛和考古学、材料科学这类人文与自然科学交叉的领域结合了起来。题目给了一堆古代玻璃文物的化学成分检测数据比如二氧化硅、氧化钠、氧化钾这些氧化物的含量百分比然后让你去分析这些玻璃的风化情况、判断它们的类型比如是钾玻璃还是铅钡玻璃甚至还要推测文物之间的关联性。这本质上是一个典型的数据分析问题但背景赋予了它独特的挑战和魅力。为什么说它适合用Python来做因为整个解题流程从数据清洗、探索性分析到特征工程、建立分类与预测模型再到结果可视化Python的生态库Pandas, NumPy, Scikit-learn, Matplotlib/Seaborn几乎提供了“一站式”的解决方案。你不需要在多个软件间来回切换一个Jupyter Notebook或者PyCharm项目就能搞定所有。对于参赛队伍来说这意味着更高的效率和更流畅的协作。这篇分享我就结合去年的解题实战把用Python处理这类成分分析鉴别问题的完整思路、核心代码以及那些容易踩的“坑”系统地梳理一遍目标是让你读完就能上手复现并掌握处理类似多维数据分析问题的通用方法。2. 解题核心思路与整体设计面对这道题首要任务是理解数据并拆解问题。题目数据通常是Excel或CSV格式每一行代表一个玻璃样品列包括文物编号、表面风化情况“风化”或“未风化”、玻璃类型“高钾”、“铅钡”等以及十几种化学成分的含量。我们的目标可以分解为几个子问题1) 分析风化对化学成分的影响2) 根据成分对玻璃类型进行鉴别分类3) 对未知类型的文物进行预测4) 分析不同类别文物的化学成分关联。这形成了一个标准的数据挖掘流水线。我的整体设计遵循“数据驱动”和“模型可解释性”并重的原则。对于数学建模竞赛尤其是国赛评委非常看重你分析问题的逻辑性和模型的物理/化学意义而不仅仅是最终的准确率。因此方案设计如下数据预处理与探索性数据分析EDA这是基石。处理缺失值古玻璃数据常有成分缺失进行描述性统计可视化成分分布初步观察风化与类型对成分的影响。风化影响分析将数据按风化状态分组运用统计检验如T检验、Mann-Whitney U检验定量分析各成分在风化前后是否有显著差异并结合箱线图等可视化手段进行阐述。玻璃类型鉴别分类建模这是核心。将已知类型的样品作为训练集构建分类模型。这里的关键是特征化学成分的筛选与处理以及模型的选择与评估。未知类型预测与敏感性分析用训练好的模型对未知类型的样品进行预测并分析模型决策的稳定性例如通过交叉验证或扰动数据观察结果变化。相关性分析与深层次探索利用统计方法如卡方检验、方差分析或机器学习方法如决策树规则提取探究不同因素之间的关联完成题目要求的其他分析。Python在这个流程中的角色是“全能工具箱”。Pandas负责数据操纵SciPy/Statsmodels负责统计检验Scikit-learn负责机器学习建模Matplotlib/Seaborn负责让一切结果变得直观。注意竞赛中一定要先吃透题目背景。例如要知道铅钡玻璃和高钾玻璃是古代中国两种主要体系其成分差异有历史渊源。这能帮你理解数据并在解释模型结果时提供更有说服力的依据避免陷入纯数字游戏。3. 数据预处理清洗、转换与探索拿到数据后的第一步不是急着跑模型而是花足够的时间“认识”你的数据。这一步往往决定了后续分析的成败。3.1 数据读取与初步审视import pandas as pd import numpy as np import matplotlib.pyplot as plt import seaborn as sns # 设置中文显示和图形样式 plt.rcParams[font.sans-serif] [SimHei, DejaVu Sans] plt.rcParams[axes.unicode_minus] False sns.set_style(whitegrid) # 读取数据 df pd.read_excel(古代玻璃制品成分数据.xlsx) print(f数据形状: {df.shape}) print(df.info()) print(df.head())运行df.info()会立刻告诉你是否有缺失值。古玻璃数据中某些成分可能因检测限或风化流失而缺失显示为NaN。df.head()则让你对数据样貌有个直观感受。3.2 缺失值处理与特征工程处理缺失值需要谨慎。对于化学成分直接删除缺失样本可能损失信息简单用均值填充可能引入偏差。我采用的策略是分析缺失模式使用missingno库的矩阵图看缺失是随机分布还是集中在某些特定类型或风化状态的样品上。分组建模填充如果缺失与类别明显相关例如所有风化玻璃的某个成分都缺失则按组如“高钾-风化”、“铅钡-未风化”计算中位数进行填充这比全局均值更合理。标记缺失对于某些关键成分也可以考虑将“是否缺失”作为一个新的布尔特征加入模型有时缺失本身就有信息量。# 示例按玻璃类型和风化状态分组填充某成分的缺失值 df[二氧化硅(SiO2)] df.groupby([类型, 风化])[二氧化硅(SiO2)].transform( lambda x: x.fillna(x.median()) ) # 对于分组后中位数仍为NaN的情况可能该分组数据全缺失用全局中位数兜底 if df[二氧化硅(SiO2)].isnull().any(): df[二氧化硅(SiO2)].fillna(df[二氧化硅(SiO2)].median(), inplaceTrue)此外特征工程可能包括创建总和特征所有氧化物百分比之和应接近100%考虑误差。计算总和可以检查数据一致性但通常不作为特征加入模型因为会导致多重共线性。创建比率特征例如K2O/Na2O钾钠比对于区分高钾和钠钙玻璃可能有化学意义。类型编码将“类型”、“风化”等分类变量转换为数值标签LabelEncoder或独热编码OneHotEncoder。3.3 探索性数据分析EDAEDA的目标是形成假设指导后续建模。我会做以下几类可视化成分分布直方图看每个化学成分的分布是否正态有无异常值。component_cols [col for col in df.columns if col not in [文物编号, 类型, 风化]] df[component_cols].hist(bins30, figsize(20, 15)) plt.suptitle(化学成分含量分布直方图) plt.show()风化与类型对成分影响的箱线图这是最直观的对比。fig, axes plt.subplots(4, 4, figsize(20, 20)) axes axes.ravel() for idx, col in enumerate(component_cols[:16]): # 假设有16种成分 sns.boxplot(x类型, ycol, hue风化, datadf, axaxes[idx]) axes[idx].set_title(f{col} by 类型 and 风化) axes[idx].tick_params(axisx, rotation45) plt.tight_layout() plt.show()从箱线图可以初步判断哪些成分在“高钾”和“铅钡”玻璃间差异明显风化是否使得某种成分如Na2O、K2O含量系统性降低。相关性热力图查看化学成分之间的相关性。高相关性共线性可能影响某些线性模型的稳定性。plt.figure(figsize(16, 12)) corr_matrix df[component_cols].corr() sns.heatmap(corr_matrix, annotTrue, fmt.2f, cmapcoolwarm, center0) plt.title(化学成分相关性热力图) plt.show()实操心得EDA阶段要多看图多思考化学解释。比如你可能会发现PbO和BaO高度正相关这很合理因为它们常共同出现在铅钡玻璃中。这个发现可以支持你在后续特征选择时考虑使用主成分分析PCA来降维或者直接选取代表性成分。4. 风化影响的定量统计分析题目要求分析风化对化学成分的影响这需要用统计检验来提供量化证据。我们不能仅凭“看起来不一样”就下结论。4.1 统计检验方法选择正态性与方差齐性检验首先对每个成分按“风化”和“未风化”分组用Shapiro-Wilk检验正态性用Levene检验方差齐性。参数检验 vs. 非参数检验如果数据满足正态分布且方差齐使用独立样本T检验否则使用Mann-Whitney U检验秩和检验。对于小样本数据非参数检验更稳健。4.2 Python实现与结果解读from scipy import stats def compare_weathering(component): 比较某成分在风化与未风化组间的差异 weathered df[df[风化]风化][component].dropna() unweathered df[df[风化]未风化][component].dropna() # 1. 正态性检验 (Shapiro-Wilk) _, p_normal_w stats.shapiro(weathered) _, p_normal_uw stats.shapiro(unweathered) normal (p_normal_w 0.05) and (p_normal_uw 0.05) # 原假设数据来自正态分布 # 2. 方差齐性检验 (Levene) _, p_var stats.levene(weathered, unweathered) equal_var p_var 0.05 # 原假设方差相等 # 3. 选择检验方法 if normal and equal_var: # T检验 t_stat, p_val stats.ttest_ind(weathered, unweathered, equal_varTrue) test_used 独立样本T检验 else: # Mann-Whitney U检验 u_stat, p_val stats.mannwhitneyu(weathered, unweathered, alternativetwo-sided) test_used Mann-Whitney U检验 # 计算效应量 (Cohens d for T, Rank-biserial correlation for U) if test_used 独立样本T检验: pooled_std np.sqrt(((len(weathered)-1)*weathered.std()**2 (len(unweathered)-1)*unweathered.std()**2) / (len(weathered)len(unweathered)-2)) d (weathered.mean() - unweathered.mean()) / pooled_std effect_size fCohens d {d:.3f} else: # 秩双列相关 r 1 - (2 * u_stat) / (len(weathered) * len(unweathered)) effect_size fRank-biserial r {r:.3f} return { 成分: component, 检验方法: test_used, p值: p_val, 效应量: effect_size, 风化组均值: weathered.mean(), 未风化组均值: unweathered.mean(), 差异显著(p0.05): p_val 0.05 } # 对每个成分进行检验 results [] for col in component_cols: results.append(compare_weathering(col)) results_df pd.DataFrame(results) print(results_df.sort_values(p值))4.3 结果分析与表述得到检验结果后如何写在论文里列出表格将成分、检验方法、p值、效应量、均值等制成清晰表格。重点解读挑选p值显著如0.05且效应量较大的成分进行重点说明。例如结果可能显示Na2O、K2O在风化后含量显著降低p0.001效应量大这与化学常识一致因为碱金属氧化物易溶于水流失。而SiO2、Al2O3等可能变化不显著因为它们化学性质稳定。结合可视化在论文中附上这些关键成分的箱线图或小提琴图让统计结论有直观支撑。注意事项不要滥用p值。p0.05只意味着“在统计上显著”但不一定代表“在实际化学意义上影响巨大”。一定要结合效应量如Cohen‘s d和专业知识共同判断。一个成分变化了0.1%即使p值显著其实际意义也可能有限。5. 玻璃类型鉴别分类模型构建与评估这是题目的核心即建立一个能根据化学成分准确预测玻璃类型的分类器。我们以已知类型的样品为训练集。5.1 特征选择与数据准备首先分离特征X和标签y。特征就是化学成分含量标签是玻璃类型。# 假设我们已经处理好了缺失值并创建了特征DataFrame X和标签Series y from sklearn.model_selection import train_test_split # 划分训练集和测试集用于评估模型泛化能力 X_train, X_test, y_train, y_test train_test_split(X, y, test_size0.2, random_state42, stratifyy) # stratifyy 确保训练集和测试集中各类别比例与原数据集一致直接使用所有成分作为特征可能存在问题1) 特征间可能存在多重共线性2) 某些成分对分类贡献很小是噪声。因此特征选择很重要过滤法计算每个特征与标签之间的相关性对于数值标签用相关系数对于分类标签可以用ANOVA F值。SelectKBest可以帮我们选择排名靠前的K个特征。包裹法如递归特征消除RFE它使用一个基模型如逻辑回归、SVM来迭代选择特征。效果更好但计算量更大。嵌入法使用带正则化的模型如L1正则化的逻辑回归、基于特征重要度的树模型训练然后根据模型系数或重要性选择特征。对于此题我推荐结合树模型的特征重要性和领域知识进行选择。先跑一个随机森林看看哪些成分最重要。from sklearn.ensemble import RandomForestClassifier rf RandomForestClassifier(n_estimators100, random_state42) rf.fit(X_train, y_train) # 获取特征重要性 importances rf.feature_importances_ feature_names X_train.columns feat_imp_df pd.DataFrame({特征: feature_names, 重要性: importances}).sort_values(重要性, ascendingFalse) plt.figure(figsize(10,6)) sns.barplot(x重要性, y特征, datafeat_imp_df) plt.title(随机森林特征重要性排序) plt.tight_layout() plt.show()5.2 模型选择与训练没有哪个模型永远最好需要尝试和比较。对于成分数据通常特征数~15不算多样本量~100也有限以下模型值得尝试逻辑回归简单、可解释性强。但前提是特征间线性可分且需要处理共线性。支持向量机在高维空间表现好尤其是使用RBF核时。但对参数和尺度敏感。随机森林对非线性关系、特征交互捕捉能力强不易过拟合还能给出特征重要性。是此类问题的“常胜将军”。XGBoost/LightGBM更强大的梯度提升树精度可能更高但需要更多调参。建议构建一个简单的模型流水线进行比较from sklearn.linear_model import LogisticRegression from sklearn.svm import SVC from sklearn.ensemble import RandomForestClassifier from sklearn.preprocessing import StandardScaler from sklearn.pipeline import Pipeline from sklearn.model_selection import cross_val_score # 定义模型和管道 models { LR: Pipeline([(scaler, StandardScaler()), (clf, LogisticRegression(max_iter1000, random_state42))]), SVM: Pipeline([(scaler, StandardScaler()), (clf, SVC(kernelrbf, random_state42))]), RF: RandomForestClassifier(n_estimators100, random_state42) } # 使用交叉验证评估 cv_scores {} for name, model in models.items(): scores cross_val_score(model, X_train, y_train, cv5, scoringaccuracy) cv_scores[name] scores.mean() print(f{name} 平均交叉验证准确率: {scores.mean():.4f} (/- {scores.std():.4f}))5.3 模型评估与优化准确率不是唯一指标尤其是当类别不平衡时。要查看混淆矩阵、精确率、召回率和F1-score。from sklearn.metrics import classification_report, confusion_matrix, ConfusionMatrixDisplay # 选择表现最好的模型例如随机森林在测试集上最终评估 best_model RandomForestClassifier(n_estimators100, random_state42) best_model.fit(X_train, y_train) y_pred best_model.predict(X_test) print(测试集分类报告:) print(classification_report(y_test, y_pred)) # 绘制混淆矩阵 cm confusion_matrix(y_test, y_pred, labelsbest_model.classes_) disp ConfusionMatrixDisplay(confusion_matrixcm, display_labelsbest_model.classes_) disp.plot(cmapplt.cm.Blues) plt.title(混淆矩阵 - 测试集) plt.show()如果模型表现不佳考虑特征工程尝试创建交互项、多项式特征需谨慎防过拟合或使用PCA/LDA降维。处理类别不平衡如果“高钾”和“铅钡”样本数相差很大在模型中使用class_weightbalanced参数或使用过采样/欠采样技术。超参数调优使用GridSearchCV或RandomizedSearchCV搜索最优参数组合。from sklearn.model_selection import GridSearchCV param_grid { n_estimators: [50, 100, 200], max_depth: [None, 10, 20, 30], min_samples_split: [2, 5, 10], min_samples_leaf: [1, 2, 4] } grid_search GridSearchCV(RandomForestClassifier(random_state42), param_grid, cv5, scoringaccuracy, n_jobs-1) grid_search.fit(X_train, y_train) print(f最佳参数: {grid_search.best_params_}) print(f最佳交叉验证分数: {grid_search.best_score_:.4f})6. 未知类型预测与模型稳定性分析训练好最终模型后就可以对题目中那些“类型未知”的文物进行预测了。这一步很简单# 假设 X_unknown 是未知类型文物的特征数据 unknown_predictions best_model.predict(X_unknown) unknown_probabilities best_model.predict_proba(X_unknown) # 获取预测概率 # 将结果保存 result_df pd.DataFrame({ 文物编号: unknown_ids, # 对应的文物编号 预测类型: unknown_predictions, 预测概率: [max(probs) for probs in unknown_probabilities] # 取最大概率作为置信度 })但竞赛论文不能只给出一个干巴巴的预测结果。你需要进行敏感性分析或不确定性分析来增强结论的可靠性。常用方法交叉验证预测分布对未知样本用交叉验证的每个模型都预测一次看预测结果是否稳定。如果所有折的模型都给出相同预测则信心足如果分歧大则结果不确定。Bootstrap重采样从训练集中有放回地多次抽样训练多个模型用这些模型对未知样本预测统计预测类型的分布。概率解释predict_proba给出的概率本身是一种不确定性度量。你可以设定一个置信度阈值如0.8低于此阈值的预测标记为“不确定”并在论文中讨论。from sklearn.model_selection import cross_val_predict from sklearn.base import clone # 使用交叉验证获取对训练集每个样本的“样本外”预测模拟模型不确定性 cv StratifiedKFold(n_splits5, shuffleTrue, random_state42) all_predictions [] for train_idx, val_idx in cv.split(X, y): model_clone clone(best_model) model_clone.fit(X.iloc[train_idx], y.iloc[train_idx]) fold_preds model_clone.predict(X.iloc[val_idx]) all_predictions.extend(fold_preds) # 可以分析这些预测与最终模型预测的一致性在论文中你需要展示对关键未知样品的预测结果并附上置信度或稳定性分析说明你的预测是稳健的。7. 关联分析与深层次探索题目可能还要求分析“不同类别玻璃文物化学成分之间的关联关系”。这可以从多个角度入手7.1 统计检验分析组间差异类似于风化分析但现在是按“玻璃类型”分组比较各成分的差异。使用方差分析ANOVA针对多组或Kruskal-Wallis H检验非参数版ANOVA。这可以定量验证我们之前从箱线图看到的直观差异是否统计显著。from scipy.stats import f_oneway, kruskal def anova_by_type(component): groups [df[df[类型]t][component].dropna() for t in df[类型].unique()] # 先进行正态性和方差齐性检验决定用ANOVA还是Kruskal-Wallis # ... (检验代码省略) # 假设满足条件使用ANOVA f_stat, p_val f_oneway(*groups) return p_val # 对每个成分计算p值 anova_results {col: anova_by_type(col) for col in component_cols}7.2 主成分分析与可视化主成分分析PCA是一种强大的探索工具。它可以将十几种相关的化学成分降维到2-3个主成分这两个主成分包含了原始数据的大部分方差。然后在二维或三维散点图上可视化所有样品按“类型”和“风化”着色可以非常直观地看到不同类别的样品是否在空间中形成聚类。from sklearn.decomposition import PCA from sklearn.preprocessing import StandardScaler # 标准化数据PCA对尺度敏感 scaler StandardScaler() X_scaled scaler.fit_transform(df[component_cols]) # 应用PCA降至2维 pca PCA(n_components2) X_pca pca.fit_transform(X_scaled) # 创建可视化 plt.figure(figsize(10,8)) scatter plt.scatter(X_pca[:,0], X_pca[:,1], cpd.factorize(df[类型])[0], cmapviridis, alpha0.7, edgecolorsk, s80) plt.colorbar(scatter, label玻璃类型) plt.xlabel(f主成分1 ({pca.explained_variance_ratio_[0]:.2%}方差)) plt.ylabel(f主成分2 ({pca.explained_variance_ratio_[1]:.2%}方差)) plt.title(古代玻璃样品PCA降维可视化按类型着色) # 可以再画一个按风化状态着色的图进行对比 plt.show() # 查看主成分的载荷每个原始成分对主成分的贡献 loadings pd.DataFrame(pca.components_.T, columns[PC1, PC2], indexcomponent_cols) print(loadings.sort_values(PC1, ascendingFalse))通过PCA图你可能发现“高钾”和“铅钡”玻璃在主成分1上被清晰分开而“风化”样品可能沿着主成分2方向有偏移。载荷矩阵则告诉你是哪些原始成分如PbO, BaO, K2O, Na2O主导了这种分离这为你的化学解释提供了直接依据。7.3 聚类分析验证除了有监督的分类也可以使用无监督的聚类方法如K-Means、层次聚类来看看数据本身是否能自然聚成几类并与已知的标签对比这可以验证分类的合理性。from sklearn.cluster import KMeans # 使用PCA降维后的数据或标准化后的原始数据 kmeans KMeans(n_clusters2, random_state42) # 假设我们预期两类 cluster_labels kmeans.fit_predict(X_scaled) # 比较聚类结果与真实类型 comparison pd.DataFrame({真实类型: df[类型], 聚类标签: cluster_labels}) ct pd.crosstab(comparison[真实类型], comparison[聚类标签]) print(ct)如果聚类结果与真实类型高度一致说明化学成分数据本身确实能很好地区分这两种玻璃这从另一个角度支撑了你的分类模型。8. 完整流程整合与论文写作要点将以上所有步骤串联起来就形成了一个完整的分析流程。在竞赛论文写作中你需要清晰地呈现这个逻辑链问题重述与分析用自己语言简述问题并拆解为几个子问题。模型假设与符号说明列出合理的假设如“各化学成分含量测量独立”、“风化过程对同类型玻璃影响模式一致”并定义文中用到的主要符号。数据预处理展示清洗、缺失值处理、探索性分析的过程和发现。模型建立与求解子问题一风化影响描述统计检验方法与结果附图表。子问题二分类鉴别阐述特征选择、模型比较、参数调优、最终模型评估的全过程。流程图和结果表如准确率对比、混淆矩阵很重要。子问题三预测与敏感性分析给出预测结果并展示稳定性分析如Bootstrap置信区间、交叉验证一致性。子问题四关联分析展示PCA可视化、聚类分析或相关性分析结果并给出化学或考古学解释。模型评价与推广讨论模型的优点如准确率高、可解释性强、局限性如样本量小、可能过拟合以及改进方向。可以简要说明模型方法同样适用于其他类似的文物成分分析或材料鉴别问题。避坑技巧与心得代码与论文分离但关键结果对应论文里不要贴大段代码只展示核心流程图、结果图和表。但务必确保你论文中的每一个数字、每一张图都能从你的Python代码中复现出来。养成使用random_state参数固定随机种子的习惯确保结果可重现。图表美观专业Matplotlib/Seaborn的默认样式比较基础。花点时间调整颜色、字体、标签让图表看起来更专业。使用一致的配色方案如Set2, tab10等分类色板。解释重于结果国赛评委很多是数学或相关领域的老师他们看重你对问题的理解和用数学工具解决问题的能力。对于模型输出的结果一定要尝试给出合理的、基于题目背景的解释。例如为什么随机森林中PbO的特征重要性最高因为铅钡玻璃以PbO和BaO为主要特征成分这与历史知识相符。团队协作与版本控制使用Git管理代码和论文版本用Jupyter Notebook或VS Code with Python插件进行开发方便分工和整合。将数据处理、模型训练、绘图等功能封装成函数或类提高代码复用性和可读性。时间管理三天时间非常紧张。建议第一天完成数据理解和预处理、基础EDA和风化分析第二天集中攻坚分类模型完成训练和评估第三天进行深入分析PCA、聚类等、预测未知数据、撰写论文初稿和修改润色。留出足够时间用于论文写作和检查。最后把所有的分析代码、生成的图表和最终数据结果整理好。一个清晰的代码结构和详尽的注释不仅方便队友检查万一需要回溯或修改时也能节省大量时间。这道题通过Python的全面应用完美诠释了如何将数学建模、数据科学和领域知识相结合解决一个具体的交叉学科问题。 
长河君 · 2026/8/27 11:55:47 · 阅读 1.2k Python数据分析实战：数学建模竞赛中的古代玻璃成分鉴别与分类
1. 项目概述当数学建模遇见文化遗产去年带队参加国赛拿到C题《古代玻璃制品的成分分析与鉴别》时第一感觉是题目出得很有意思。它把数学建模这个看似“硬核”的竞赛和考古学、材料科学这类人文与自然科学交叉的领域结合了起来。题目给了一堆古代玻璃文物的化学成分检测数据比如二氧化硅、氧化钠、氧化钾这些氧化物的含量百分比然后让你去分析这些玻璃的风化情况、判断它们的类型比如是钾玻璃还是铅钡玻璃甚至还要推测文物之间的关联性。这本质上是一个典型的数据分析问题但背景赋予了它独特的挑战和魅力。为什么说它适合用Python来做因为整个解题流程从数据清洗、探索性分析到特征工程、建立分类与预测模型再到结果可视化Python的生态库Pandas, NumPy, Scikit-learn, Matplotlib/Seaborn几乎提供了“一站式”的解决方案。你不需要在多个软件间来回切换一个Jupyter Notebook或者PyCharm项目就能搞定所有。对于参赛队伍来说这意味着更高的效率和更流畅的协作。这篇分享我就结合去年的解题实战把用Python处理这类成分分析鉴别问题的完整思路、核心代码以及那些容易踩的“坑”系统地梳理一遍目标是让你读完就能上手复现并掌握处理类似多维数据分析问题的通用方法。2. 解题核心思路与整体设计面对这道题首要任务是理解数据并拆解问题。题目数据通常是Excel或CSV格式每一行代表一个玻璃样品列包括文物编号、表面风化情况“风化”或“未风化”、玻璃类型“高钾”、“铅钡”等以及十几种化学成分的含量。我们的目标可以分解为几个子问题1) 分析风化对化学成分的影响2) 根据成分对玻璃类型进行鉴别分类3) 对未知类型的文物进行预测4) 分析不同类别文物的化学成分关联。这形成了一个标准的数据挖掘流水线。我的整体设计遵循“数据驱动”和“模型可解释性”并重的原则。对于数学建模竞赛尤其是国赛评委非常看重你分析问题的逻辑性和模型的物理/化学意义而不仅仅是最终的准确率。因此方案设计如下数据预处理与探索性数据分析EDA这是基石。处理缺失值古玻璃数据常有成分缺失进行描述性统计可视化成分分布初步观察风化与类型对成分的影响。风化影响分析将数据按风化状态分组运用统计检验如T检验、Mann-Whitney U检验定量分析各成分在风化前后是否有显著差异并结合箱线图等可视化手段进行阐述。玻璃类型鉴别分类建模这是核心。将已知类型的样品作为训练集构建分类模型。这里的关键是特征化学成分的筛选与处理以及模型的选择与评估。未知类型预测与敏感性分析用训练好的模型对未知类型的样品进行预测并分析模型决策的稳定性例如通过交叉验证或扰动数据观察结果变化。相关性分析与深层次探索利用统计方法如卡方检验、方差分析或机器学习方法如决策树规则提取探究不同因素之间的关联完成题目要求的其他分析。Python在这个流程中的角色是“全能工具箱”。Pandas负责数据操纵SciPy/Statsmodels负责统计检验Scikit-learn负责机器学习建模Matplotlib/Seaborn负责让一切结果变得直观。注意竞赛中一定要先吃透题目背景。例如要知道铅钡玻璃和高钾玻璃是古代中国两种主要体系其成分差异有历史渊源。这能帮你理解数据并在解释模型结果时提供更有说服力的依据避免陷入纯数字游戏。3. 数据预处理清洗、转换与探索拿到数据后的第一步不是急着跑模型而是花足够的时间“认识”你的数据。这一步往往决定了后续分析的成败。3.1 数据读取与初步审视import pandas as pd import numpy as np import matplotlib.pyplot as plt import seaborn as sns # 设置中文显示和图形样式 plt.rcParams[font.sans-serif] [SimHei, DejaVu Sans] plt.rcParams[axes.unicode_minus] False sns.set_style(whitegrid) # 读取数据 df pd.read_excel(古代玻璃制品成分数据.xlsx) print(f数据形状: {df.shape}) print(df.info()) print(df.head())运行df.info()会立刻告诉你是否有缺失值。古玻璃数据中某些成分可能因检测限或风化流失而缺失显示为NaN。df.head()则让你对数据样貌有个直观感受。3.2 缺失值处理与特征工程处理缺失值需要谨慎。对于化学成分直接删除缺失样本可能损失信息简单用均值填充可能引入偏差。我采用的策略是分析缺失模式使用missingno库的矩阵图看缺失是随机分布还是集中在某些特定类型或风化状态的样品上。分组建模填充如果缺失与类别明显相关例如所有风化玻璃的某个成分都缺失则按组如“高钾-风化”、“铅钡-未风化”计算中位数进行填充这比全局均值更合理。标记缺失对于某些关键成分也可以考虑将“是否缺失”作为一个新的布尔特征加入模型有时缺失本身就有信息量。# 示例按玻璃类型和风化状态分组填充某成分的缺失值 df[二氧化硅(SiO2)] df.groupby([类型, 风化])[二氧化硅(SiO2)].transform( lambda x: x.fillna(x.median()) ) # 对于分组后中位数仍为NaN的情况可能该分组数据全缺失用全局中位数兜底 if df[二氧化硅(SiO2)].isnull().any(): df[二氧化硅(SiO2)].fillna(df[二氧化硅(SiO2)].median(), inplaceTrue)此外特征工程可能包括创建总和特征所有氧化物百分比之和应接近100%考虑误差。计算总和可以检查数据一致性但通常不作为特征加入模型因为会导致多重共线性。创建比率特征例如K2O/Na2O钾钠比对于区分高钾和钠钙玻璃可能有化学意义。类型编码将“类型”、“风化”等分类变量转换为数值标签LabelEncoder或独热编码OneHotEncoder。3.3 探索性数据分析EDAEDA的目标是形成假设指导后续建模。我会做以下几类可视化成分分布直方图看每个化学成分的分布是否正态有无异常值。component_cols [col for col in df.columns if col not in [文物编号, 类型, 风化]] df[component_cols].hist(bins30, figsize(20, 15)) plt.suptitle(化学成分含量分布直方图) plt.show()风化与类型对成分影响的箱线图这是最直观的对比。fig, axes plt.subplots(4, 4, figsize(20, 20)) axes axes.ravel() for idx, col in enumerate(component_cols[:16]): # 假设有16种成分 sns.boxplot(x类型, ycol, hue风化, datadf, axaxes[idx]) axes[idx].set_title(f{col} by 类型 and 风化) axes[idx].tick_params(axisx, rotation45) plt.tight_layout() plt.show()从箱线图可以初步判断哪些成分在“高钾”和“铅钡”玻璃间差异明显风化是否使得某种成分如Na2O、K2O含量系统性降低。相关性热力图查看化学成分之间的相关性。高相关性共线性可能影响某些线性模型的稳定性。plt.figure(figsize(16, 12)) corr_matrix df[component_cols].corr() sns.heatmap(corr_matrix, annotTrue, fmt.2f, cmapcoolwarm, center0) plt.title(化学成分相关性热力图) plt.show()实操心得EDA阶段要多看图多思考化学解释。比如你可能会发现PbO和BaO高度正相关这很合理因为它们常共同出现在铅钡玻璃中。这个发现可以支持你在后续特征选择时考虑使用主成分分析PCA来降维或者直接选取代表性成分。4. 风化影响的定量统计分析题目要求分析风化对化学成分的影响这需要用统计检验来提供量化证据。我们不能仅凭“看起来不一样”就下结论。4.1 统计检验方法选择正态性与方差齐性检验首先对每个成分按“风化”和“未风化”分组用Shapiro-Wilk检验正态性用Levene检验方差齐性。参数检验 vs. 非参数检验如果数据满足正态分布且方差齐使用独立样本T检验否则使用Mann-Whitney U检验秩和检验。对于小样本数据非参数检验更稳健。4.2 Python实现与结果解读from scipy import stats def compare_weathering(component): 比较某成分在风化与未风化组间的差异 weathered df[df[风化]风化][component].dropna() unweathered df[df[风化]未风化][component].dropna() # 1. 正态性检验 (Shapiro-Wilk) _, p_normal_w stats.shapiro(weathered) _, p_normal_uw stats.shapiro(unweathered) normal (p_normal_w 0.05) and (p_normal_uw 0.05) # 原假设数据来自正态分布 # 2. 方差齐性检验 (Levene) _, p_var stats.levene(weathered, unweathered) equal_var p_var 0.05 # 原假设方差相等 # 3. 选择检验方法 if normal and equal_var: # T检验 t_stat, p_val stats.ttest_ind(weathered, unweathered, equal_varTrue) test_used 独立样本T检验 else: # Mann-Whitney U检验 u_stat, p_val stats.mannwhitneyu(weathered, unweathered, alternativetwo-sided) test_used Mann-Whitney U检验 # 计算效应量 (Cohens d for T, Rank-biserial correlation for U) if test_used 独立样本T检验: pooled_std np.sqrt(((len(weathered)-1)*weathered.std()**2 (len(unweathered)-1)*unweathered.std()**2) / (len(weathered)len(unweathered)-2)) d (weathered.mean() - unweathered.mean()) / pooled_std effect_size fCohens d {d:.3f} else: # 秩双列相关 r 1 - (2 * u_stat) / (len(weathered) * len(unweathered)) effect_size fRank-biserial r {r:.3f} return { 成分: component, 检验方法: test_used, p值: p_val, 效应量: effect_size, 风化组均值: weathered.mean(), 未风化组均值: unweathered.mean(), 差异显著(p0.05): p_val 0.05 } # 对每个成分进行检验 results [] for col in component_cols: results.append(compare_weathering(col)) results_df pd.DataFrame(results) print(results_df.sort_values(p值))4.3 结果分析与表述得到检验结果后如何写在论文里列出表格将成分、检验方法、p值、效应量、均值等制成清晰表格。重点解读挑选p值显著如0.05且效应量较大的成分进行重点说明。例如结果可能显示Na2O、K2O在风化后含量显著降低p0.001效应量大这与化学常识一致因为碱金属氧化物易溶于水流失。而SiO2、Al2O3等可能变化不显著因为它们化学性质稳定。结合可视化在论文中附上这些关键成分的箱线图或小提琴图让统计结论有直观支撑。注意事项不要滥用p值。p0.05只意味着“在统计上显著”但不一定代表“在实际化学意义上影响巨大”。一定要结合效应量如Cohen‘s d和专业知识共同判断。一个成分变化了0.1%即使p值显著其实际意义也可能有限。5. 玻璃类型鉴别分类模型构建与评估这是题目的核心即建立一个能根据化学成分准确预测玻璃类型的分类器。我们以已知类型的样品为训练集。5.1 特征选择与数据准备首先分离特征X和标签y。特征就是化学成分含量标签是玻璃类型。# 假设我们已经处理好了缺失值并创建了特征DataFrame X和标签Series y from sklearn.model_selection import train_test_split # 划分训练集和测试集用于评估模型泛化能力 X_train, X_test, y_train, y_test train_test_split(X, y, test_size0.2, random_state42, stratifyy) # stratifyy 确保训练集和测试集中各类别比例与原数据集一致直接使用所有成分作为特征可能存在问题1) 特征间可能存在多重共线性2) 某些成分对分类贡献很小是噪声。因此特征选择很重要过滤法计算每个特征与标签之间的相关性对于数值标签用相关系数对于分类标签可以用ANOVA F值。SelectKBest可以帮我们选择排名靠前的K个特征。包裹法如递归特征消除RFE它使用一个基模型如逻辑回归、SVM来迭代选择特征。效果更好但计算量更大。嵌入法使用带正则化的模型如L1正则化的逻辑回归、基于特征重要度的树模型训练然后根据模型系数或重要性选择特征。对于此题我推荐结合树模型的特征重要性和领域知识进行选择。先跑一个随机森林看看哪些成分最重要。from sklearn.ensemble import RandomForestClassifier rf RandomForestClassifier(n_estimators100, random_state42) rf.fit(X_train, y_train) # 获取特征重要性 importances rf.feature_importances_ feature_names X_train.columns feat_imp_df pd.DataFrame({特征: feature_names, 重要性: importances}).sort_values(重要性, ascendingFalse) plt.figure(figsize(10,6)) sns.barplot(x重要性, y特征, datafeat_imp_df) plt.title(随机森林特征重要性排序) plt.tight_layout() plt.show()5.2 模型选择与训练没有哪个模型永远最好需要尝试和比较。对于成分数据通常特征数~15不算多样本量~100也有限以下模型值得尝试逻辑回归简单、可解释性强。但前提是特征间线性可分且需要处理共线性。支持向量机在高维空间表现好尤其是使用RBF核时。但对参数和尺度敏感。随机森林对非线性关系、特征交互捕捉能力强不易过拟合还能给出特征重要性。是此类问题的“常胜将军”。XGBoost/LightGBM更强大的梯度提升树精度可能更高但需要更多调参。建议构建一个简单的模型流水线进行比较from sklearn.linear_model import LogisticRegression from sklearn.svm import SVC from sklearn.ensemble import RandomForestClassifier from sklearn.preprocessing import StandardScaler from sklearn.pipeline import Pipeline from sklearn.model_selection import cross_val_score # 定义模型和管道 models { LR: Pipeline([(scaler, StandardScaler()), (clf, LogisticRegression(max_iter1000, random_state42))]), SVM: Pipeline([(scaler, StandardScaler()), (clf, SVC(kernelrbf, random_state42))]), RF: RandomForestClassifier(n_estimators100, random_state42) } # 使用交叉验证评估 cv_scores {} for name, model in models.items(): scores cross_val_score(model, X_train, y_train, cv5, scoringaccuracy) cv_scores[name] scores.mean() print(f{name} 平均交叉验证准确率: {scores.mean():.4f} (/- {scores.std():.4f}))5.3 模型评估与优化准确率不是唯一指标尤其是当类别不平衡时。要查看混淆矩阵、精确率、召回率和F1-score。from sklearn.metrics import classification_report, confusion_matrix, ConfusionMatrixDisplay # 选择表现最好的模型例如随机森林在测试集上最终评估 best_model RandomForestClassifier(n_estimators100, random_state42) best_model.fit(X_train, y_train) y_pred best_model.predict(X_test) print(测试集分类报告:) print(classification_report(y_test, y_pred)) # 绘制混淆矩阵 cm confusion_matrix(y_test, y_pred, labelsbest_model.classes_) disp ConfusionMatrixDisplay(confusion_matrixcm, display_labelsbest_model.classes_) disp.plot(cmapplt.cm.Blues) plt.title(混淆矩阵 - 测试集) plt.show()如果模型表现不佳考虑特征工程尝试创建交互项、多项式特征需谨慎防过拟合或使用PCA/LDA降维。处理类别不平衡如果“高钾”和“铅钡”样本数相差很大在模型中使用class_weightbalanced参数或使用过采样/欠采样技术。超参数调优使用GridSearchCV或RandomizedSearchCV搜索最优参数组合。from sklearn.model_selection import GridSearchCV param_grid { n_estimators: [50, 100, 200], max_depth: [None, 10, 20, 30], min_samples_split: [2, 5, 10], min_samples_leaf: [1, 2, 4] } grid_search GridSearchCV(RandomForestClassifier(random_state42), param_grid, cv5, scoringaccuracy, n_jobs-1) grid_search.fit(X_train, y_train) print(f最佳参数: {grid_search.best_params_}) print(f最佳交叉验证分数: {grid_search.best_score_:.4f})6. 未知类型预测与模型稳定性分析训练好最终模型后就可以对题目中那些“类型未知”的文物进行预测了。这一步很简单# 假设 X_unknown 是未知类型文物的特征数据 unknown_predictions best_model.predict(X_unknown) unknown_probabilities best_model.predict_proba(X_unknown) # 获取预测概率 # 将结果保存 result_df pd.DataFrame({ 文物编号: unknown_ids, # 对应的文物编号 预测类型: unknown_predictions, 预测概率: [max(probs) for probs in unknown_probabilities] # 取最大概率作为置信度 })但竞赛论文不能只给出一个干巴巴的预测结果。你需要进行敏感性分析或不确定性分析来增强结论的可靠性。常用方法交叉验证预测分布对未知样本用交叉验证的每个模型都预测一次看预测结果是否稳定。如果所有折的模型都给出相同预测则信心足如果分歧大则结果不确定。Bootstrap重采样从训练集中有放回地多次抽样训练多个模型用这些模型对未知样本预测统计预测类型的分布。概率解释predict_proba给出的概率本身是一种不确定性度量。你可以设定一个置信度阈值如0.8低于此阈值的预测标记为“不确定”并在论文中讨论。from sklearn.model_selection import cross_val_predict from sklearn.base import clone # 使用交叉验证获取对训练集每个样本的“样本外”预测模拟模型不确定性 cv StratifiedKFold(n_splits5, shuffleTrue, random_state42) all_predictions [] for train_idx, val_idx in cv.split(X, y): model_clone clone(best_model) model_clone.fit(X.iloc[train_idx], y.iloc[train_idx]) fold_preds model_clone.predict(X.iloc[val_idx]) all_predictions.extend(fold_preds) # 可以分析这些预测与最终模型预测的一致性在论文中你需要展示对关键未知样品的预测结果并附上置信度或稳定性分析说明你的预测是稳健的。7. 关联分析与深层次探索题目可能还要求分析“不同类别玻璃文物化学成分之间的关联关系”。这可以从多个角度入手7.1 统计检验分析组间差异类似于风化分析但现在是按“玻璃类型”分组比较各成分的差异。使用方差分析ANOVA针对多组或Kruskal-Wallis H检验非参数版ANOVA。这可以定量验证我们之前从箱线图看到的直观差异是否统计显著。from scipy.stats import f_oneway, kruskal def anova_by_type(component): groups [df[df[类型]t][component].dropna() for t in df[类型].unique()] # 先进行正态性和方差齐性检验决定用ANOVA还是Kruskal-Wallis # ... (检验代码省略) # 假设满足条件使用ANOVA f_stat, p_val f_oneway(*groups) return p_val # 对每个成分计算p值 anova_results {col: anova_by_type(col) for col in component_cols}7.2 主成分分析与可视化主成分分析PCA是一种强大的探索工具。它可以将十几种相关的化学成分降维到2-3个主成分这两个主成分包含了原始数据的大部分方差。然后在二维或三维散点图上可视化所有样品按“类型”和“风化”着色可以非常直观地看到不同类别的样品是否在空间中形成聚类。from sklearn.decomposition import PCA from sklearn.preprocessing import StandardScaler # 标准化数据PCA对尺度敏感 scaler StandardScaler() X_scaled scaler.fit_transform(df[component_cols]) # 应用PCA降至2维 pca PCA(n_components2) X_pca pca.fit_transform(X_scaled) # 创建可视化 plt.figure(figsize(10,8)) scatter plt.scatter(X_pca[:,0], X_pca[:,1], cpd.factorize(df[类型])[0], cmapviridis, alpha0.7, edgecolorsk, s80) plt.colorbar(scatter, label玻璃类型) plt.xlabel(f主成分1 ({pca.explained_variance_ratio_[0]:.2%}方差)) plt.ylabel(f主成分2 ({pca.explained_variance_ratio_[1]:.2%}方差)) plt.title(古代玻璃样品PCA降维可视化按类型着色) # 可以再画一个按风化状态着色的图进行对比 plt.show() # 查看主成分的载荷每个原始成分对主成分的贡献 loadings pd.DataFrame(pca.components_.T, columns[PC1, PC2], indexcomponent_cols) print(loadings.sort_values(PC1, ascendingFalse))通过PCA图你可能发现“高钾”和“铅钡”玻璃在主成分1上被清晰分开而“风化”样品可能沿着主成分2方向有偏移。载荷矩阵则告诉你是哪些原始成分如PbO, BaO, K2O, Na2O主导了这种分离这为你的化学解释提供了直接依据。7.3 聚类分析验证除了有监督的分类也可以使用无监督的聚类方法如K-Means、层次聚类来看看数据本身是否能自然聚成几类并与已知的标签对比这可以验证分类的合理性。from sklearn.cluster import KMeans # 使用PCA降维后的数据或标准化后的原始数据 kmeans KMeans(n_clusters2, random_state42) # 假设我们预期两类 cluster_labels kmeans.fit_predict(X_scaled) # 比较聚类结果与真实类型 comparison pd.DataFrame({真实类型: df[类型], 聚类标签: cluster_labels}) ct pd.crosstab(comparison[真实类型], comparison[聚类标签]) print(ct)如果聚类结果与真实类型高度一致说明化学成分数据本身确实能很好地区分这两种玻璃这从另一个角度支撑了你的分类模型。8. 完整流程整合与论文写作要点将以上所有步骤串联起来就形成了一个完整的分析流程。在竞赛论文写作中你需要清晰地呈现这个逻辑链问题重述与分析用自己语言简述问题并拆解为几个子问题。模型假设与符号说明列出合理的假设如“各化学成分含量测量独立”、“风化过程对同类型玻璃影响模式一致”并定义文中用到的主要符号。数据预处理展示清洗、缺失值处理、探索性分析的过程和发现。模型建立与求解子问题一风化影响描述统计检验方法与结果附图表。子问题二分类鉴别阐述特征选择、模型比较、参数调优、最终模型评估的全过程。流程图和结果表如准确率对比、混淆矩阵很重要。子问题三预测与敏感性分析给出预测结果并展示稳定性分析如Bootstrap置信区间、交叉验证一致性。子问题四关联分析展示PCA可视化、聚类分析或相关性分析结果并给出化学或考古学解释。模型评价与推广讨论模型的优点如准确率高、可解释性强、局限性如样本量小、可能过拟合以及改进方向。可以简要说明模型方法同样适用于其他类似的文物成分分析或材料鉴别问题。避坑技巧与心得代码与论文分离但关键结果对应论文里不要贴大段代码只展示核心流程图、结果图和表。但务必确保你论文中的每一个数字、每一张图都能从你的Python代码中复现出来。养成使用random_state参数固定随机种子的习惯确保结果可重现。图表美观专业Matplotlib/Seaborn的默认样式比较基础。花点时间调整颜色、字体、标签让图表看起来更专业。使用一致的配色方案如Set2, tab10等分类色板。解释重于结果国赛评委很多是数学或相关领域的老师他们看重你对问题的理解和用数学工具解决问题的能力。对于模型输出的结果一定要尝试给出合理的、基于题目背景的解释。例如为什么随机森林中PbO的特征重要性最高因为铅钡玻璃以PbO和BaO为主要特征成分这与历史知识相符。团队协作与版本控制使用Git管理代码和论文版本用Jupyter Notebook或VS Code with Python插件进行开发方便分工和整合。将数据处理、模型训练、绘图等功能封装成函数或类提高代码复用性和可读性。时间管理三天时间非常紧张。建议第一天完成数据理解和预处理、基础EDA和风化分析第二天集中攻坚分类模型完成训练和评估第三天进行深入分析PCA、聚类等、预测未知数据、撰写论文初稿和修改润色。留出足够时间用于论文写作和检查。最后把所有的分析代码、生成的图表和最终数据结果整理好。一个清晰的代码结构和详尽的注释不仅方便队友检查万一需要回溯或修改时也能节省大量时间。这道题通过Python的全面应用完美诠释了如何将数学建模、数据科学和领域知识相结合解决一个具体的交叉学科问题。