简介：2025年妈妈杯B题完整参赛资料，面向数学建模参赛者及生物信息学爱好者，围绕结肠癌基因表达图谱展开分析。论文综合运用GB指数、BP神经网络、小波变换与贝叶斯估计等方法，完整呈现从无关基因筛选、特征基因提取到数据去噪与未知基因探索的建模过程。资源包仅含1个docx文档，压缩包大小2.56MB，内容覆盖问题重述、基本假设、符号说明、模型建立与求解等论文结构，并附带处理流程与对应代码结果，便于直接对照学习。目前已有626人学习下载。文档在问题一中给出GB综合指数筛选114个信息基因的阈值选择方法；问题二利用相关性剔除与MIV方法确定含12个基因的最优组合；问题三用MATLAB小波工具箱去噪后保留61个基因并提取8个特征基因；问题四结合聚类与Bayes估计给出探索未知基因的思路，对理解多学科方法在基因研究中的应用很有帮助。
1. 2025妈妈杯B题基因信息提取：从2000个基因到12个标签的降维路线
2025年妈妈杯B题的完整论文+代码结果包里，最值钱的不是那套漂亮的排版，而是它把“基因表达谱分类”这个高维问题拆成了一条可复现的流水线：先用GB综合指数把2000个基因砍到114个，再用BP神经网络配合MIV逐步剔除锁定12个最优基因组合，最后用小波去噪反推验证。做这类题最怕的就是“基因越多信息越多”的直觉——真实情况是大部分基因在正常人和结肠癌患者两组样本里的分布几乎没有差异，留着它们只会让分类器无所适从。这份资源适合正在备赛数学建模、以及需要做高维特征筛选但不想从零开始造轮子的人。下面按论文的实际求解顺序，把每个环节的算法逻辑、参数依据和复现时容易翻车的地方拆开讲。
2. GB综合指数选基因：Gini指数与Bhattacharyya距离怎么互补、阈值怎么定
2.1 单看一个指标为什么不够
问题一的本质是剔除“无关基因”。判断一个基因是否与分类相关，最朴素的做法是看它在两类样本中的均值差异，但均值相同不代表分布相同——方差差异同样携带分类信息。论文把Gini指数和Bhattacharyya距离组合成GB综合指标，核心动机就是让两个指标互相补位。
Gini指数衡量的是表达值在不同类别间的“不均衡程度”，值越小说明该基因在不同类别中的表达分布越分离；Bhattacharyya距离则同时考虑了均值的差异和方差的差异，值越大代表两类样本的可分性越好。单独用Gini，容易被表达水平整体偏低或偏高的基因干扰；单独用Bhattacharyya距离，又容易漏掉那些均值接近但分布形状差异明显的基因。两个指标各排一次序、再取交集，能有效降低单指标排序带来的误选风险。
2.2 Gini指数：21级离散化与升序排序的细节
计算Gini指数之前，论文先把每个基因的表达值离散化到0—20共21个等级。这一步很关键：原始芯片数据是连续值，直接套Gini公式会受离群值影响，离散化相当于做了一个非参数化的稳健变换。
离散化公式可以理解为：对每个基因g，找到它在所有样本中的最大值max和最小值min，然后把当前样本值n映射到0-20的整数区间。映射时加了0.5的补偿量再取整，本质是四舍五入而不是向下截断。这样做的好处是等级边界不会把恰好处于中间值的样本错误归入低一级。
离散化之后，每个基因对类别k的Gini指标为：
[ Gini(k) = 1 - \sum_{j=0}^{20} p_{ij}^2 ]
其中p_ij是基因i在等级j上属于类别k的相对频率。整体Gini值再按两类样本数加权平均。理解这个指标的关键点在于：当某一类别的所有样本都集中在同一个表达等级时，p_ij接近1，Gini(k)接近0，说明这个基因对分类贡献最大；反之如果样本散布在多个等级，Gini(k)变大，分类价值就低。所以后面排序采用升序，Gini值越小越靠前。从论文给出的部分基因数据表可以看到，排名靠前的基因Gini值普遍在0.87-0.92区间，听起来数值接近，但排序后仍然能拉开差距。
2.3 Bhattacharyya距离：0.05阈值为什么能把无关基因划走
Bhattacharyya距离的计算公式综合了两类样本的均值和方差：
[ B = \frac{1}{4} \cdot \frac{(\mu_1 - \mu_2)^2}{\sigma_1^2 + \sigma_2^2} + \frac{1}{2} \ln\left(\frac{\sigma_1^2 + \sigma_2^2}{2\sigma_1\sigma_2}\right) ]
前一项衡量均值差异的贡献，后一项衡量方差差异的贡献。距离越大，两类分布的重叠越少，基因的可分性越好。
论文统计了2000个基因的Bhattacharyya距离分布，78.55%的基因距离落在0到0.05之间。这是一个非常有说服力的信号：这些基因在正常人和结肠癌患者样本中的均值和方差都没有明显区别，属于典型的“无关基因”。于是阈值定在0.05，小于该值的基因直接剔除。
这里有个值得留意的细节：0.05不是通过交叉验证选出来的，而是基于分布频数的一个自然截断点。实际操作中这个阈值提供了很大的容错空间，78.55%的基因都挤在0-0.05区间，即使阈值稍微浮动，对筛选结果的稳定性影响也不大。在建模比赛中，这种基于数据分布观察来定参数的做法比拍脑袋定阈值更容易在论文里写清楚。
2.4 两组300取交集：114个信息基因的计算流程
论文把m取为基因总量的15%，也就是300。这个比例的依据可以从两条排序曲线看出来：Gini指数升序排列后，曲线在前300个基因之后进入平缓区；Bhattacharyya距离降序排列后，前300个基因贡献了绝大部分的分类距离。
具体流程是：先把2000个基因按Gini值升序取前300，再按Bhattacharyya距离降序取前300，然后取两组备用基因的交集。交集内的基因在两个指标上同时表现优异，共得到114个信息基因。排序时以Bhattacharyya距离排名为主、Gini排名为辅，这样最终得到的GB综合排名实际上是以“可分性”为第一优先级、以“分布不均衡度”为第二优先级的综合排序。
这一步做完，维度从2000降到114，降幅超过94%。从信息保留角度看，后续问题二的12基因最优组合就是从这114个基因里继续筛出来的，说明这114个基因确实兜住了核心分类信息。
3. BP神经网络+MIV逐步剔除：12个最优基因组合的确定过程
3.1 强相关冗余剔除：Pearson阈值0.85时保留46个基因
进入114个基因之后，第一步是剔除冗余。基因之间存在调控关系，表达水平会呈现相关性，两个高度相关的基因同时留在特征集合里，并不会增加分类信息量，反而给后续的MIV计算增加干扰。
论文用Pearson相关系数计算114个基因两两之间的相关性，然后设置不同阈值做“两两冗余”分析：相关系数超过阈值的基因对，剔除其中GB综合排名较低的那个。
这里有一个值得注意的实验过程：论文分别测试了1、0.9、0.85、0.8、0.75、0.725六个阈值，对应的剩余基因数量分别是114、83、46、30、17、10。分类错误数对应为2、2、3、5、5、6。阈值从1降到0.85时，基因数量从114锐减到46，但分类错误数保持2不变；再往下压到0.8，基因只剩30个，错误数却跳到3。这说明0.85是一个临界点——它去掉了冗余信息但没有损伤分类能力。最终的46个基因集合作为下一阶段的输入。
这个阈值选择过程非常典型：先做敏感性分析，再选“分类能力不下降的最小特征集”。而不是机械地把阈值调到0.8或0.9。
3.2 MIV算法拆解：每次都加10%扰动再训练
MIV（Mean Impact Value，平均影响值）是这一问的技术核心。它的思想很直观：如果某个基因对分类结果真的重要，那么把它的表达值人为上下扰动10%，网络输出应该发生明显变化；反之，如果怎么扰动输出都不动，说明这个基因可有可无。
具体计算流程可以分为四步：
第一步，用当前候选基因子集训练一个BP神经网络，输入节点数等于基因个数，输出节点为样本类别。训练结束后网络权重固定。
第二步，对训练样本中的每一个特征，分别在其原值基础上加10%和减10%，构造出两组新样本P1和P2。注意P1和P2的样本数量与原始样本一致，只是某一列的表达值整体偏移。
第三步，把P1和P2分别送入已经训练好的网络做仿真预测，得到两个输出矩阵A1和A2。
第四步，计算IV = A1 - A2，再将所有样本上的IV值取平均，得到该基因的MIV值。MIV的符号代表影响方向，绝对值代表影响强度。
以下是MATLAB风格的MIV计算骨架：
% 假设net是已经训练好的BP网络，P是训练样本矩阵（行=样本，列=基因） % 需要对第col个基因计算MIV P1 = P; P2 = P; P1(:, col) = P(:, col) * 1.1; % 原始值加10% P2(:, col) = P(:, col) * 0.9; % 原始值减10% A1 = sim(net, P1'); % 仿真得到加扰动后的输出 A2 = sim(net, P2'); % 仿真得到减扰动后的输出 MIV = (A1 - A2) / size(P, 1); % 按样本数平均这段代码的核心逻辑是：网络训练好之后不再更新权重，只改变输入列的值来观察输出变化。10%的扰动幅度不是固定的，如果你的基因表达值动态范围很大，可以改成5%或15%，但注意扰动太小会被网络输出的舍入误差淹没，太大则会跨过非线性区间，得到一个不真实的梯度。常见做法是先看数据的标准差，扰动幅度取标准差的一半左右，效果比较稳。
MIV比直接看网络权重更可靠的地方在于：它是在整个网络的输入输出映射关系下计算的影响值，基因之间的交互作用会被网络结构捕捉到。多个基因协同影响分类时，单独看权重矩阵是看不出来的，但扰动一个基因后，误差会通过隐藏层传播到输出层，MIV能够反映出包括交互效应在内的综合影响。
3.3 后10%逐轮剔除与BP错判率检验
直接对所有46个基因一次性算完MIV然后砍掉后10%，这种做法虽然快，但有缺陷：基因之间存在冗余关系，单独计算MIV时排在末尾的基因，可能在删掉另一个基因后变得重要。论文采用逐步剔除法来规避这个问题：
第一轮：用46个基因训练BP网络，计算每个基因的MIV，剔除绝对值最小的10%（约5个基因），保留41个。
第二轮：用剩下的41个基因重新训练BP网络，重新计算MIV，再剔除后10%。每轮都是重新训练、重新计算，而不是沿用上一轮的MIV排序。
循环往复，直到候选基因集合为空。每一轮剔除后，都用当前基因子集训练BP网络，记录分类错判数。
最终选择标准有两个：错判率最低、基因数量最少。这两个目标是有冲突的——基因越多，分类器越容易在训练集上拟合到低错误率，但泛化能力未必好。所以论文的选法是：不追求全局最低错判数，而是在错判数可接受的范围内选基因数量最少的组合。
最后锁定的12个基因组合为：M85079、T62947、R39209、R84411、T54303、M82919、H43887、X12671、H08393、M26383、R36977、R87126。这12个基因对应的错判数在论文的BP网络框架下做到了最低，同时基因数量在同样错判水平下最少。
3.4 最终12个基因组合的验证视角
拿到12个基因后，论文用自组织竞争神经网络做了进一步的分类效果检验。和BP网络不同，自组织竞争网络不需要标签信息，属于无监督竞争学习，它能把样本按表达模式自动聚成两类。用两种不同原理的分类器交叉验证12基因组合，比单独依赖BP网络的说服力强得多。
从特征选择的角度看，这12个基因组合的价值在于：它不是一个一个独立挑出来的“最强单基因”堆在一起的组合，而是经过逐轮MIV剔除后整体表现最优的功能组合。基因之间可能有相互补偿效应，单独排名靠前的基因组合在一起反而未必是最优的。这也是为什么论文坚持“逐轮剔除+重新训练”而不是一次性排序截断的根本原因。
4. 小波去噪模型：MATLAB小波工具箱处理基因信号的三步流程
4.1 为什么去噪要选小波而不是均值/中值
问题三把每个基因的表达值序列看作一个信号：有用信号x_i叠加噪声n_i。噪声来源包括芯片制造、荧光标记、杂交过程等，论文假设噪声是零均值高斯白噪声。
对比均值去噪和中值去噪，小波变换的优势在于它在时域和频域同时具有局部化能力。基因芯片数据的特点是样本数量少（62个样本）、基因维度高（2000个基因），均值滤波会抹平真实表达峰值的毛刺，中值滤波对高斯噪声的抑制效果一般。小波去噪可以把信号分解成低频近似部分和高频细节部分，真实表达趋势集中在前者，噪声主要分布在后者，两者在频域上可以区分开。
从样本数量的角度讲，小波变换对短序列信号也能处理到较深的分解层数，这是它适合基因数据的原因之一。
4.2 分解-加阈值-重建：一条命令行下来的操作
小波去噪的三步流程在MATLAB小波工具箱里对应一组明确的函数调用，整个流程可以用如下代码结构表示：
% Step 1: 信号分解，选择小波基和分解层数 [C, L] = wavedec(expr_signal, 3, 'db4'); % C为各层小波系数，L为对应长度 % 分解层数3，小波基db4是Daubechies系列中常用选择 % Step 2: 对细节系数加阈值 % 阈值规则可选 sqtwolog / rigrsure / heursure sigma = median(abs(C)) / 0.6745; % 噪声标准差估计 thr = thselect(expr_signal, 'heursure'); % 启发式阈值 C_filt = wthresh(C, 's', thr); % 软阈值处理细节系数 % Step 3: 信号重建 recon_signal = waverec(C_filt, L, 'db4');逻辑说明：wavedec把原始表达信号分解成三层，每层有对应的近似系数和细节系数。噪声主要集中在细节系数上，阈值函数thr的选择直接影响保留多少细节。wthresh的's'表示软阈值处理——小于阈值的系数直接置零，大于阈值的系数向零收缩。软阈值处理后的信号更平滑，硬阈值（'h'）保留的细节更多但容易出现震荡。
参数说明：分解层数取3层，对62个样本的表达序列来说是合适的，分解层数太多会把近似系数也拆碎，太少则无法有效分离噪声频带。小波基选db4，它在平滑性和局部性之间的平衡较好，是生物信号处理里的常用选择。阈值规则用heursure是兼顾软阈值和风险估计的折中方案，sqtwolog更激进，rigrsure更保守，实际操作时可以三个规则各跑一遍，对比去噪后筛出的基因数量再做决定。
4.3 去噪后重新筛基因：61个和8个怎么解读
去噪完成后，论文把去噪后的基因表达数据重新走了一遍问题一的GB指数筛选流程，对比结果值得深究：去噪后的数据在做基因分类时保留了61个基因，比原始数据少53个；特征基因进一步提取后得到8个，比问题二的12个还少4个。
这个结果的正确解读方式不是“去噪后筛选更严格了所以基因更少了”，而是：原始数据里有一部分基因的分类信息其实是噪声伪影。噪声的存在让某些基因在两类样本中的分布差异被放大了，这些基因通过了GB指数筛选，但并不代表真实的生物学分类信号。去噪之后，这些虚假差异被抹平，基因筛出的数量和特征基因数量都显著下降，说明去噪让筛选结果更贴近真实表达信号。
这个发现对论文的价值很大：它证明了问题二得到的12个基因组合里可能混入了噪声驱动的冗余基因，而去噪后的8个基因才是更纯粹的分类标签。换言之，小波去噪不仅是预处理手段，更是一种对现有筛选结果有效性的反向验证。
5. 复现基因筛选论文的常见问题与避坑记录
5.1 Gini指数算出来全接近1，排序形同虚设
现象：对原始表达值直接算Gini指数，2000个基因的Gini值全部集中在0.95-0.999区间，排序后前300个和后300个之间几乎没有区分度，筛选结果和随机抽样差不多。
原因：没有先做离散化。原始芯片表达值范围很大且分布偏斜，直接套Gini公式时，每个等级上的频率都被摊薄了，p_ij²之和非常接近0，Gini值自然趋近于1。离散化到0-20的21个等级，本质上是把连续分布的尾部差异压缩进了等级边界，让Gini指数对“类别间分布差异”更敏感。
解决：先按论文公式做0-20等级离散化，再算Gini值。离散化时要注意个别离群样本会把max/min区间拉得过大，导致大多数样本挤在低等级区间。常见做法是在离散化前先做一次1%和99%分位的截尾，或者先用Z-score归一化再映射到0-20区间，离散后的分布会更均匀，Gini值的区分度明显上升。
5.2 Bhattacharyya距离出现NaN
现象：用公式计算Bhattacharyya距离时，个别基因的输出是NaN或Inf，导致排序时报错。
原因：在所有62个样本中，某些基因的表达值可能完全相同（尤其在离散化之前检查原始值，或某一类样本中所有值都相等），此时σ2为零，公式第二项的分母出现0，对数项无定义。
解决：在计算协方差之前先检查每个基因在两类样本中的方差，如果方差小于1e-10，直接把这个基因视为无分类信息并划入无关基因集合。另一种做法是给σ1²和σ2²各加一个极小量epsilon（比如1e-8）防止除零。这个处理不影响正常基因的距离值，但能避免整个批次的计算中断。
5.3 BP网络收敛到极小误差，MIV排序反而失真
现象：为12个基因训练BP网络时，训练误差能达到1e-14量级，训练集分类完全正确，但用MIV值排序筛选出的基因和直观认知不符，排名靠前的基因换到测试集上效果很差。
原因：这是典型的过拟合场景——隐层节点数设置过大、训练轮数过长时，网络记住了训练样本的个体细节而不是分类规律。MIV的扰动信号在网络处于过拟合状态时会包含大量噪声信息，导致影响值排序失去统计意义。论文中Performance逼近1.9e-14，说明这类数模论文中的超参数设置本身也踩在这条线上。
解决：控制隐层节点数，经验公式是输入节点数的1/3到1/2量级，12个基因对应5-6个隐层节点。训练时保留一部分样本做验证集，验证误差不再下降就提前停止，而不是等训练误差压到0再停。MIV不是越精细越好——扰动幅度10%时，如果网络对训练样本记忆得太精确，输出差值反映的其实是单个训练样本的局部梯度，而不是基因的普适影响。
5.4 相关系数阈值越低基因越少，但分类错误却反弹
现象：按论文的表3做阈值敏感性测试，阈值从0.85降到0.8后，基因数量从46降到30，但分类错误数反而从2升到3；继续降0.725，基因只剩10个，错误数升到6。
原因：两个基因相关系数高，不代表它们携带的信息完全重复。强行剔除掉所有高相关基因中的一方，可能把一组中某个“单独看冗余、组合后有补位作用”的基因误删。阈值压得越低，误删风险越大。
解决：不要只看阈值与基因数量、错误数的关系，要把每一轮阈值下保留的基因集合保存下来，做分类错误率的置信区间分析。或者换个思路：先保留高相关基因对，在MIV剔除阶段再统一处理。论文最终选0.85，因为它兼顾了降维比例和分类能力，这是基于多轮实验数据选出来的，而不是拍脑袋定的。
5.5 小波默认阈值把真实表达差异也抹平了
现象：用小波工具箱默认的sqtwolog阈值做去噪，去噪后重新筛基因，得到的基因数量比预期少太多，连问题二锁定的12个基因里有价值的几个也被过滤掉了。
原因：sqtwolog阈值是固定阈值，它基于噪声方差的通用估计，对以“整体分布差异”为核心的基因筛选任务来说过于激进。基因表达数据的有效信号不集中在前几个低频系数里，部分分类信息恰好分布在幅值较小的细节系数中，固定阈值把这些细节一起置零了。
解决：改用heursure或rigrsure自适应阈值，它们会根据信号的局部风险来调整阈值强度。另一个可行方案是分层设置阈值：对第一层细节系数用较大阈值，对深层细节系数用较小阈值，因为深层细节中可能包含真实的表达趋势。每次去噪后都重新跑一遍GB指标筛选对比，如果保留基因数和分类错误率同时下降，说明阈值方向不对。
6. 把妈妈杯B题的方案迁移到其他赛题：验证与改造技巧
6.1 用交叉验证确认基因组合不是偶然
MIV逐步剔除法有一个隐患：整个筛选过程是在同一份数据上训练和评估的，最后选出的12个基因可能只是在这62个样本上表现最优的过拟合组合。确认它不是偶然的方式是重采样验证：随机打乱样本顺序，以不同的训练/测试划分重新跑完整套流程，看选出的基因组合是否稳定。如果每次跑出来的12基因集合重合率超过70%，说明这套筛选结果对数据的敏感度不高；如果重合率很低，那问题二的结果需要标注“特定划分下的最优组合”而不是普适结论。
6.2 用更现代的树模型特征重要性做交叉验证
论文用BP网络是因为赛题背景限定在神经网络框架内。迁移到其他数据集时，我一般会用LightGBM的特征重要性排序和MIV排序做对比：同一个基因集合，如果在两种模型下的排名高度一致，说明这些基因的分类贡献是模型无关的；如果冲突严重，则说明基因之间交互复杂，单一模型的特征排序不可靠。这种交叉验证成本很低，但能大幅提升论文结论的可信度。
6.3 去噪前后对比图的画法
问题三的数据结果非常适合做对比图：用两条曲线分别表示去噪前后2000个基因的Bhattacharyya距离降序分布。如果去噪有效，你会看到距离大的基因数量减少、曲线整体下移，说明部分“分类优势”是噪声造成的假象。再叠一张Gini指数排序曲线，就能直观说明去噪对两类指标的差异化影响。这种图放到论文里，比单独贴去噪前后的数字更有说服力。
另外提醒一句：问题四的贝叶斯聚类部分在论文里篇幅最短，但其思想很实用——用已有信息基因作为先验，结合质心法确定初始聚类中心，再通过Bayes后验概率更新样本归属。实际复现时可以把它当作半监督聚类来处理：已知基因标签是带标记样本，未知基因是待标记样本，每次迭代更新聚类中心和先验概率直到收敛，效果远好于一次性分类。
记得我第一次复现这类基因筛选流程时，一上来就盯着MIV绝对值最大的基因看，换个数据矩阵结果就失灵了。从那以后，我拿到任何特征筛选结果都会强制走一遍重采样验证：把样本随机打乱三次，看选出的基因组合是否稳定。希望这份资料的拆解能帮你少走这段弯路，完整论文和代码结果在文章末尾可以获取，对照着跑一遍比只看文字理解深得多。
本文还有配套的精品资源，点击获取
