机器学习在热电材料发现中的应用：数据分割与特征选择策略

1. 项目概述：当机器学习遇见热电材料发现

在材料科学领域，寻找下一个“明星”材料——比如能将废热高效转化为电能的热电材料——传统上依赖于“试错法”或基于物理原理的复杂计算，过程既昂贵又缓慢。近年来，机器学习（ML）如同一把锋利的“数据手术刀”，为我们开辟了一条从海量数据中直接挖掘材料性能规律的捷径。然而，这条路并非坦途。许多研究者都遇到过这样的困境：模型在训练集上表现优异，各项指标（如R²、RMSE）近乎完美，但一旦面对实验室里合成出的全新材料，预测结果却差强人意。这背后，往往不是算法不够先进，而是数据处理的根基没有打牢。

问题的核心在于两个常被忽视但至关重要的环节：数据如何分割，以及特征如何选择。想象一下，如果你只用亚洲人的面部数据训练一个人脸识别模型，然后让它去识别欧洲人，效果必然不佳。同样，如果训练集只覆盖了材料化学空间中某个狭窄的区域，模型就无法学会应对更广阔、更多样的未知材料。另一方面，特征选择就像为模型配备“感官”。如果输入了上百个冗余甚至无关的原子属性（如原子半径、电负性、蒸发热等），模型不仅会“消化不良”，计算负担剧增，还可能被噪声误导，学到虚假的关联，这就是过拟合。

因此，一个稳健的机器学习工作流，其起点必须是构建一个能够真实反映问题复杂性的无偏数据集，并为其配备一套精炼而有效的“特征语言”。本文将以half-Heusler（半赫斯勒）热电材料的发现为例，深入探讨我们如何构建一个融合了主成分分析（PCA）数据分割与集成化特征选择策略的机器学习框架。这个框架的核心目标不是追求训练集上的最低误差，而是确保模型在面对完全陌生的化学组成时，依然能做出可靠预测。我们最终将这个框架应用于对约6.6亿种潜在成分的高通量虚拟筛选，并成功锁定了一批具有高热电优值（zT）潜力的全新候选材料。接下来，我将拆解这个过程中的每一个关键决策、实操细节以及我们踩过的坑。

2. 核心挑战与设计思路：为何传统方法会“翻车”？

在深入技术细节之前，我们必须先理解在热电材料机器学习中，为什么常规操作容易失效。这有助于我们后续的每一个设计选择都有的放矢。

2.1 数据分割的陷阱：化学空间的代表性

热电实验数据集通常具有层次结构：一个材料（如TiNiSn）可以有多种掺杂成分（如Ti₀.₉Zr₀.₁NiSn），而每种成分又在多个温度下被测量。最常见的错误分割方式是按数据点随机划分。这会导致什么后果？假设一个材料的所有测量数据点（不同温度、不同掺杂浓度）被随机分到了训练集和测试集，那么模型在测试时，本质上还是在“回忆”它从同一个材料的其他数据点中学到的模式。这测试的是模型的插值能力，而非面对全新化学组成的外推能力。我们的目标是发现新材料，因此测试集必须包含模型从未“见过”的完整材料。

更隐蔽的陷阱是化学空间的偏见。如果训练集里的材料都集中在元素周期表的某个区域（例如，全是过渡金属基材料），而测试集是另一类材料（例如，稀土金属基材料），那么模型在测试集上的糟糕表现，反映的不是模型不好，而是训练数据没有覆盖到测试集所在的化学空间。因此，一个“公平”的划分，必须确保训练集和测试集在由元素特征所张成的高维化学空间中，具有相似且全面的分布。

2.2 特征选择的迷宫：从“维数灾难”到物理可解释性

对于half-Heusler材料（化学式通常为ABC），我们为每个A、B、C位点及其可能的掺杂元素（A‘， B’， C‘）选取了19个原子属性（如原子序数、电负性、蒸发热等）。仅此一项，就构成了一个114维的特征向量。再加上3个掺杂浓度和1个温度，总维度高达118维。然而，我们的数据集仅有1521个数据点（来自108种材料）。这陷入了典型的“维数灾难”：特征数量远多于样本数。

不加选择地使用所有特征会带来多重问题：

1. 计算负担：训练时间，尤其是超参数优化时间，呈指数级增长。
2. 过拟合风险：模型容易记住噪声而非规律，泛化能力差。
3. 多重共线性：许多原子属性之间高度相关（例如，原子半径和原子量），这会使模型系数不稳定，难以解释。

因此，特征选择不是“可选项”，而是“必选项”。但如何选？简单的过滤法（如按与zT的相关性排序）速度快，但忽略了特征间的相互作用。封装法（如递归特征消除）考虑模型性能，但计算成本极高，且结果严重依赖于特定的模型和超参数。我们需要一个既能降维、又能保持甚至提升预测性能的稳健策略。

2.3 整体工作流设计

基于以上挑战，我们设计了如下图所示的端到端工作流。它的核心创新在于将数据分割的公正性与特征选择的系统性深度耦合，并最终服务于高通量筛选。

原始数据集 (108种材料， 1521个数据点) ↓ 基于化学空间（114维元素特征）的PCA分析 ↓ PCA空间内基于距离质心的排序与模块化采样 ↓ 生成5个在化学空间上均匀分布的训练/测试折 (Folds) ↓ (对每个Fold独立进行) 特征预处理 → 相关性过滤 → k-Best特征筛选 + 贝叶斯超参优化 ↓ 训练多个模型 (RF, XGBoost, NN, SISSO) ↓ 在独立测试集上评估，计算泛化误差 ↓ 对所有Fold的预测进行集成，得到最终稳健模型 ↓ 应用于约6.6亿候选成分的高通量筛选 ↓ 输出高潜力新材料候选列表

这个流程确保了从数据准备到最终预测的每一步，都以模型泛化能力为最高指导原则。

3. 实操详解一：基于PCA的无偏数据分割

这是整个流程的基石。我们的目标不是随机分割，而是在化学特征空间中进行分层抽样，确保每个“折”都能代表整体。

3.1 构建材料的化学“指纹”

首先，我们将每种half-Heusler材料（例如TiNiSn）表示为一个固定向量。我们使用其组成元素的19个原子属性，为A、A‘、B、B’、C、C‘六个位点各生成一个19维向量，拼接成一个114维的“化学指纹”。这个指纹与具体的掺杂浓度和测量温度无关，它只表征材料的元素身份。这样，每种材料在化学空间中就是一个唯一的点。

3.2 PCA降维与空间映射

直接将114维空间用于距离计算和可视化是困难且不直观的。我们采用主成分分析（PCA）进行降维。PCA会找到数据方差最大的几个正交方向（主成分）。我们保留累计贡献率达到95%的前13个主成分。这13维空间在最大程度上保留了原始化学空间的分布结构，同时极大地简化了后续操作。

注意：这里95%的阈值是一个经验选择。在实际操作中，我们需要绘制碎石图，观察特征值下降的拐点。通常选择拐点之前的主成分，在信息保留和维度精简之间取得平衡。盲目追求99%的方差可能会引入大量噪声维度。

3.3 模块化采样生成五折交叉验证集

接下来是关键步骤。在13维PCA空间中，我们计算每种材料点到所有点质心的欧氏距离，��按距离从小到大排序。这个距离反映了材料在化学空间中的“边缘”程度。

然后，我们采用一种模块化采样策略来创建5个独立的训练/测试分割（即5折）：

• Fold 1 测试集：包含排名为第1、6、11、16...的材料。
• Fold 2 测试集：包含排名为第2、7、12、17...的材料。
• 以此类推，Fold 5 测试集：包含排名为第5、10、15、20...的材料。

这种做法的精妙之处在于：

1. 全面覆盖：每个Fold的测试集都包含了距离质心近、中、远的材料，确保了该Fold的测试集能广泛代表整个化学空间。
2. 相互独立：任何一个材料都只出现在一个Fold的测试集中，确保了交叉验证的独立性。
3. 挑战性：测试集不是“好预测”的简单样本，而是包含了化学空间中的各种极端和典型情况，能真正考验模型的泛化能力。

下图是一个二维t-SNE投影的示意图（t-SNE用于可视化，实际分割在PCA空间进行），可以直观看到五个测试集的材料点（不同颜色）均匀地散布在整个化学空间中，没有形成聚集。

(注：此处为示意图，实际应根据论文中Figure 2描述)

3.4 实操心得与避坑指南

• 不要用t-SNE进行实际分割：t-SNE是一种优秀的非线性可视化工具，但它不保持全局距离，且结果受参数影响大。用它来分割数据会导致不稳定的分组。PCA是线性且确定性的，更适合作为分割的基础。
• 层次分割是关键：务必在材料级别进行分割，而不是在数据点级别。确保同一种材料的所有数据（不同温度、不同掺杂）要么全在训练集，要么全在测试集。这可以通过在分割前对数据按材料ID进行分组聚合来实现。
• 检查分布：分割后，务必检查每个Fold中训练集和测试集的zT值分布（均值、方差、直方图）。理想情况下，它们应该大致相似。如果某个Fold的测试集zT值普遍偏高或偏低，说明分割可能引入了偏差。

4. 实操详解二：两阶段特征选择策略

有了高质量的数据分割，下一步是为模型准备“营养均衡”的输入特征。我们采用了一个两阶段策略：先粗筛，再精挑。

4.1 第一阶段：基于相关性的冗余特征过滤

目标：剔除高度相关的特征，因为它们提供的是重复信息。

步骤1：计算特征-zT相关性由于数据是分层的（一个材料对应多个zT值），我们不能直接计算。我们进行了两次聚合：对每个成分取其在所有温度下的最大zT，再对每个材料取其所有成分中的最大zT。这样就得到了一个“材料-zT（最大值）”的一对一映射列表。在此基础上，计算每个114维元素特征与这个聚合zT之间的皮尔逊相关系数的绝对值。绝对值高的特征，被认为与zT的线性关联更强。

步骤2：构建特征相关性图并去冗余我们分别为A/A‘、B/B’、C/C‘三组位点构建特征相关性矩阵。因为同一组位点的候选元素池相同，其特征分布也相同。

1. 设定一个相关系数阈值r_min（例如0.8）。
2. 将每个特征视为图中的一个节点。如果两个特征之间的相关系数绝对值大于r_min，则在它们之间连一条边。
3. 使用Bron–Kerbosch算法找出图中所有的最大团。一个团就是一个完全连接的子图，意味着团内的所有特征两两高度相关。
4. 在每个团中，只保留与zT相关性绝对值最高的那个特征，剔除团内其他特征。

这个过程对每个r_min值都会产生一个精简后的特征列表。r_min的选择是个权衡：值太高（如0.95），去冗余效果弱，特征数多；值太低（如0.5），可能把有互补信息的特征也剔除，导致信息损失。下图展示了随着r_min变化，保留的特征数量k_max的变化。

(注：此处为示意图，实际应根据论文中Figure 4b描述)

4.2 第二阶段：k-Best筛选与超参数优化联合搜索

第一阶段给了我们一个去冗余后的特征排序列表（按与zT的相关性从高到低排）。但到底选前多少个特征（k）最好？这不是一个能预先回答的问题，需要与模型性能联合优化。

我们为每个Fold独立执行以下流程：

1. 网格搜索：遍历不同的r_min值（例如从-1到1取20个等间距点）和不同的k值（从1到该r_min下保留的特征总数k_max）。每一个(r_min, k)组合都定义了一个特征子集（即该r_min下去冗余列表中的前k个特征）。
2. 模型训练与验证：对于每个特征子集，我们使用该Fold的训练集（80%数据），并进一步在其内部进行10折交叉验证。在每一折内部，我们使用贝叶斯优化来搜索随机森林（RF）、XGBoost（XGB）和神经网络（NN）的最佳超参数组合。
   • 为什么用贝叶斯优化？相比网格搜索或随机搜索，贝叶斯优化利用历史评估结果来构建代理模型，预测下一个可能最优的超参数点，能以更少的迭代次数找到更优解，对于计算成本高的模型尤其合适。
3. 性能评估：记录在10折交叉验证中平均验证集RMSE最低的(r_min, k)组合。
4. 确定最终特征集：选择使验证集RMSE最小的(r_min, k)组合对应的特征子集，作为该Fold、该模型的最终特征集。
5. 最终训练：用这个最优特征子集和对应的最优超参数，在整个Fold的训练集（80%）上重新训练最终模型。

这个过程将特征选择完全嵌入到模型训练流程中，找到的是针对特定模型和特定数据子集的最优特征组合，而不是一个全局通用的“最佳特征列表”。

4.3 针对SISSO的迭代变量选择

SISSO（Sure Independent Screening and Sparsifying Operator）是一种符号回归方法，旨在找到zT与特征之间简洁的数学表达式。由于其算法需要从巨大的特征组合空间中搜索，直接应用k-Best方法计算上不可行。

我们采用了一种迭代变量选择（VS）算法：

1. 从一个随机的小特征子集S_a开始。
2. 用S_a运行SISSO，得到一个初始模型和RMSE。
3. 将模型中出现的特征移入一个固定集合S_b。
4. 向S_a中随机添加新的未探索特征。
5. 用S_a ∪ S_b再次运行SISSO。如果新模型的RMSE提升，则将新出现在模型中的特征从S_a移入S_b。
6. 重复步骤4-5，直到RMSE在连续多次迭代中不再提升。

这个方法能高效地从上百个特征中筛选出对构建SISSO描述符最关键的那一小部分。

4.4 实操心得与避坑指南

• 相关性过滤的“站点特异性”：全局看原子半径和原子量可能相关性不高，但在特定晶格位点（如C位）的有限元素集合中，它们可能高度相关。必须分位点计算相关性矩阵，否则会漏掉重要的冗余信息。
• k-Best搜索的起点：从k=1开始搜索是必要的。有时，仅仅一两个强相关特征就能达到不错的基线性能，增加特征反而引入噪声。这能帮助我们理解预测任务的“信息下限”。
• 超参数优化与特征选择的耦合：这是最耗资源但也是最关键的一步。绝对不能先固定一组“经验”超参数，再做特征选择。因为不同的特征集合下，模型的最佳超参数（如树的深度、学习率、正则化强度）可能完全不同。必须将它们联合优化。
• 记录与可视化：将每个(r_min, k)组合对应的验证RMSE绘制成等高线图或热图，可以直观地看到性能平台区域。这有助于确认选择的稳健性，避免选到某个尖锐的局部最优点。

5. 模型训练、集成与高通量筛选

经过严谨的数据分割和特征选择，我们得到了针对每个Fold和每种算法（RF, XGB, NN, SISSO）优化好的模型。接下来是整合与预测。

5.1 模型评估与集成预测

对于每个Fold，我们已经在独立的20%测试集上评估了模型性能。但我们的目标是得到一个单一的、稳健的预测模型用于筛选。

我们采用了一种基于逆误差加权的集成平均方法：

1. Fold内集成（针对同一算法）：对于一个候选材料，五个Fold的模型会给出五个预测值zT_pr_t,F。我们不是简单平均，而是根据每个Fold模型在673K（我们设定的工作温度）插值后的测试集RMSE（e_F）进行加权平均。误差越小，权重越高。zT_pr_t = Σ (zT_pr_t,F / e_F) / Σ (1/e_F)这样，在未知数据上表现更稳定的Fold拥有更大话语权。
2. 算法间集成：不同算法（SISSO, NN, XGB, RF）各有优劣。我们计算每个算法在所有数据上重构后的总体RMSE（e_t），然后再次使用逆误差加权，将四个算法的预测值融合为一个最终的“技术加权预测值”zT_pr_pon。

这种两级集成策略，既考虑了数据子集划分带来的不确定性，也融合了不同算法视角的预测，通常能获得比单一模型更稳定、更准确的结果。

5.2 高通量筛选流程

我们构建了一个巨大的虚拟化学空间，包含所有可能的ABC基体及其A‘， B’， C‘位点掺杂组合，浓度从0到0.5以不同步长离散化，总计约6.6亿种成分。直接让模型预测这么多点计算量巨大，因此需要前置过滤。

1. 稳定性过滤：我们利用Open Quantum Materials Database (OQMD)和Materials Project (MP)数据库，计算了每个ABC三元化合物的形成能凸包距离。只保留距离低于经验阈值（OQMD: 0.15 eV/atom， MP: 0.9 eV/atom）的稳定或亚稳基体。这一步淘汰了绝大多数热力学上极不稳定的组合。
2. 电子计数过滤：half-Heusler材料通常需要满足一定的价电子数范围（如17-19）以获得良好的电学性能。我们应用此规则进一步筛选。
3. 模型预测：对过滤后的候选成分，使用上述集成模型预测其在673K下的zT值。
4. 结果排序与去重：分别按每种算法和最终的加权预测值对候选成分排序，选取Top-N。同时，确保每个独特的基体元素组合（A, B, C）只保留一个最优的掺杂配方，以避免结果偏向于某个特定体系。

5.3 实操心得与避坑指南

• 温度插值的必要性：实验数据是在不同温度下测量的，而高通量筛选通常需要一个固定温度下的性能指标。我们采用二次函数拟合每种成分的zT-T实验曲线，然后插值到673K，以此作为“真实值”来评估模型在该温度下的预测误差。这比简单使用某个温度点的原始数据更合理。
• 数据库过滤的局限性：OQMD和MP的计算基于标准条件，可能漏掉一些在非平衡条件下（如快速淬火）可合成的亚稳相。阈值的选择也带有经验性。因此，筛选结果是一个“高概率”候选列表，而非绝对保证。实验验证必不可少。
• 关注“假阳性”与“假阴性”：分析那些模型预测值很高但实际可能不稳定的材料（假阳性），以及那些模型预测值低但实际有潜力的材料（假阴性），能反过来帮助我们改进特征工程或模型。例如，我们发现在预测极高zT（>1）和极低zT（<0.1）的材料时，所有模型都倾向于向均值回归（预测值更保守）。这提示数据集在极端值区域样本不足，或当前特征无法充分描述导致性能突变的物理机制（如共振能级、能带收敛）。

6. 结果解读与物理洞察

经过完整的流程，我们得到了什么？

6.1 模型性能对比

在重组后的完整测试集（每个数据点都来自其作为“未知数据”的那个Fold）上，各模型的RMSE在0.163到0.169之间。SISSO和神经网络（NN）表现最佳，RMSE约为0.164。这个误差水平与领域内先进工作相当，但更重要的是，我们的评估方式（基于化学空间的无偏分割）更能反映真实的泛化能力。

一个普遍现象是，模型对低zT（<0.2）材料倾向于高估，对高zT（>0.2）材料倾向于低估。这并非数据错误，而是源于数据集中一种特殊模式：某些本征材料zT极低，但微量掺杂后zT飙升数个数量级（如TaFeSb）。模型难以捕捉这种剧烈的“掺杂效应跃变”，导致对纯相预测过高，对掺杂相预测偏低。在面向掺杂材料的高通量筛选中，我们更关心zT>0.2的区域，在此区间模型的平均绝对百分比误差（MAPE）约为26-31%。

6.2 关键特征揭示：什么在驱动高性能？

通过SHAP分析和SISSO描述符解析，我们得以窥见模型的“思考”过程：

1. 温度（T）：毫无悬念是最重要的特征。在树模型（RF, XGB）中，它与zT呈单调正相关。但在神经网络中，SHAP值分布复杂，表明T通过与其他特征的复杂非线性相互作用来影响zT。
2. A位掺杂浓度（xA‘）：这是仅次于温度的第二重要特征。SHAP分析显示，低xA‘通常与低zT相关，但存在例外（一些蓝点出现在高SHAP值区域），证实了少量A位掺杂即可大幅提升某些材料性能的实验观察（如TaFeSb）。SISSO描述符中多次出现xA‘/ea_A‘（掺杂浓度与A’位元素电子亲和能的比值）的函数项，表明这种影响是非线性的，且与掺杂元素的本质紧密相关。
3. A位蒸发热（HV_A）：在模型无关的特征排序和模型相关的分析中均频繁出现。高HV_A通常意味着A位原子与周围原子键合较强，有利于高温下保持高电导率，从而提升zT。
4. 元素间差异：SISSO描述符中出现了如|Z_A - Z_A‘|（A位与A’位原子序数之差）的项，这指向了质量波动散射机制。A位与掺杂位原子质量差异越大，对声子的散射越强，可能降低晶格热导率，有益于zT。

这些发现不仅验证了已知的物理知识（如A位调控的重要性），也提出了新的假设（如特定电子亲和能与掺杂浓度的耦合效应），为实验设计提供了定量化的指导。

6.3 高通量筛选发现

应用集成模型对过滤后的化学空间进行预测，我们得到了数种高潜力的新型half-Heusler成分。例如，SISSO模型预测HfLiNiZnSbBi等成分在673K下zT可达1.24。值得注意的是，预测列表中出现了一些使用轻质元素（如Li, Mg）作为A位掺杂剂的新组合，这为通过引入强点缺陷散射来降低热导率提供了新思路。当然，这些预测需要后续的密度泛函理论计算和实验合成来验证。

7. 常见问题、挑战与未来方向

在实际操作这套流程时，会遇到不少典型问题。

7.1 数据量与质量瓶颈

• 问题：高质量、一致性的实验数据集规模有限（本研究仅1521个数据点）。小数据量限制了复杂模型（如深度神经网络）的发挥，也使得统计结论的不确定性增大。
• 应对：
  • 数据增强：对于成分数据，可基于物理规则生成“虚拟”数据点（需谨慎）。
  • 迁移学习：利用大型计算数据库（如OQMD）的DFT计算结果预训练模型，再在小型实验数据上微调。
  • 主动学习：这正是我们下一步的计划。将模型预测不确定性高或与实验偏差大的成分，优先进行实验合成与测量，然后将新数据加入训练集，重新训练模型。如此循环，用最少的实验成本最大化提升模型性能。

7.2 特征工程的局限性

• 问题：当前使用的都是组成特征，缺乏结构特征（如晶格常数、键长、能带结构）和微观结构特征（如晶粒尺寸、孔隙率）。而这些对热电性能（尤其是热导率）有决定性影响。
• 应对：整合多尺度特征。与高通量计算结合，获取关键的结构描述符。在实验中，尝试将工艺参数（如烧结温度、时间）作为特征，虽然这会使问题更复杂，但更贴近实际。

3. 模型的可解释性与物理一致性

• 问题：即使SHAP和SISSO提供了一些洞见，但像神经网络这样的“黑箱”模型内部机制仍不透明。模型可能学到一些物理上不合理的虚假关联。
• 应对：
  • 使用物理约束模型：在损失函数中加入物理约束项（如zT必须为正，Seebeck系数与电导率之间的某种关联）。
  • 符号回归的优先性：像SISSO这类能产生显式数学公式的方法，在可解释性上具有天然优势，应作为理解物理机制的首选工具。
  • 交叉验证：如果一个特征在多种不同的模型和分割方式下都显示出重要性，那么它的物理意义就更可信。

7.4 从预测到合成的“死亡之谷”

• 问题：模型预测出的高分材料，可能在热力学上不稳定，或难以合成。
• 应对：将稳定性预测（如凸包距离、声子稳定性）和可合成性预测（如反应能、前驱体分解温度）作为并行筛选模型，与性能预测模型组成多目标筛选管道。我们工作中使用的数据库过滤就是向这个方向迈出的第一步。

构建一个用于材料发现的稳健机器学习框架，远不止是调包和跑算法。它要求我们深刻理解数据背后的物理、化学层次，对模型评估抱有批判性态度，并精心设计每一个环节以避免数据泄露和过拟合。基于PCA的数据分割和严谨的特征选择策略，是确保模型从“实验室宠儿”走向“未知战场侦察兵”的关键两步。这套方法论不仅适用于热电材料，对于其他功能材料（如催化剂、电池材料、光电材料）的数据驱动发现，同样具有重要的参考价值。未来的方向无疑是走向更大规模、多模态数据的融合，以及“预测-合成-表征-反馈”的闭环主动学习，真正让机器学习成为加速材料研发的引擎。