做回归建模的时候最头疼的往往不是选什么模型而是怎么处理那一堆特征。手头几十上百个变量要么一股脑全丢给模型跑得慢不说还容易被无关特征带偏要么自己拍脑袋挑几个心里又没底。mRMR最大相关最小冗余就是专门解决这个问题的经典特征选择算法它用互信息来衡量特征与目标之间的相关性再把特征之间的重复信息压到最低最终输出一份清晰的特征排序。这篇文章我会用Matlab把mRMR的回归版本完整写出来从原理推导、核心代码、效果验证到常见坑位一条龙讲清楚适合正在做回归预测、苦于特征太多或者想系统入门特征选择的朋友参考。1. mRMR算法原理为什么回归场景选它1.1 最大相关用互信息代替相关系数mRMR里面第一个“MR”是Max-Relevance最大相关性。意思是说选出来的特征和目标值之间必须有足够的关联不然这个特征就是白选。大部分人第一反应是用皮尔逊相关系数来衡量相关但回归数据里特征和目标之间往往不是简单的线性关系。比如某个特征和目标之间呈现明显的U型关系皮尔逊相关系数算出来可能接近0但事实上这个特征对预测非常关键。这时候用互信息就稳得多。互信息从信息论的角度衡量两个变量之间的依赖程度它的离散形式是I(X;Y) Σ p(x,y) * log( p(x,y) / (p(x)p(y)) )简单理解就是知道了X之后对Y的不确定性减少了多少。如果减少得越多说明X和Y越“默契”。互信息不要求任何函数形式的假设线性、非线性、周期波动都能捕捉到这是它比相关系数更适合做特征筛选的根本原因。在回归场景下互信息的计算有个绕不开的问题目标y是连续值没法直接套用离散形式的公式。所以工程上最常规的做法是先对连续值做分箱离散化再用离散公式估计互信息。这也是整篇代码里最关键的预处理步骤后面会详细说。1.2 最小冗余选出来的特征要“各干各的”mRMR的另一个“MR”是Min-Redundancy最小冗余度。它要解决的是另一个问题特征A和特征B都和目标高度相关但A和B之间也高度相关那么两个都选就是浪费。用一个生活化的例子来类比一个团队如果全是同一类型的技术专家遇到突发问题就没人能解决但如果成员技能互补整体战斗力会强很多。特征选择也是同样的道理一组特征如果彼此之间只知道重复信息模型的稳定性反而会变差一旦某个特征受噪声污染整个预测都会被带偏。所以在迭代选择特征时不仅要看新特征和目标的互信息有多大还要看它和已经选出来的特征之间的平均互信息有多高。mRMR的目标函数就是把这两项放在一起权衡最终选出“相关性强、冗余度低”的特征组合。1.3 回归和分类在实现上的关键差异很多人第一次写mRMR代码时会发现网上的开源实现大多是针对分类任务的。分类时目标变量本身就是离散的类别标签直接拿来算互信息就行。但回归任务的目标是连续值比如房价、温度、销售额这种情况下必须先把目标y也离散化。代码实现上分类版和回归版的核心迭代逻辑一致差别就在数据预处理阶段。分类版输入的是一个类别向量回归版需要额外增加一步把连续目标y用等宽分箱或等频分箱变成整数标签。这一步看似简单实际影响很大分箱数少了丢信息多了互信息估计又容易受噪声干扰后面会专门讲参数怎么选。另外mRMR还有一个需要确认的选择相关性和冗余度怎么组合。常用的有两种形式结合方式公式特点MID相关性 - 冗余度计算快差值形式适合数值尺度接近的场景MIQ相关性 / 冗余度比值形式更公平但分母接近0时不稳定我自己在回归场景下更推荐MID。原因很简单连续特征离散化后互信息绝对值普遍偏小几个特征之间的差异也比较小用差值能保留更多细节而MIQ在做除法时如果某个候选特征和已选特征的冗余度恰好接近0得分会被放大得离谱容易选出边缘特征。代码里默认用的是MID想切换MIQ的话就把最后得分改成除法再加一个很小的eps防止除零。2. Matlab实现准备与整体架构设计2.1 环境要求其实不需要额外工具箱先交代一下环境。整套代码基于Matlab基础环境实现只要不是年代过于久远的版本R2015a之前都能正常跑因为用到了Matlab自带的discretize函数。如果要做后面的效果验证需要Statistics and Machine Learning Toolbox里面提供了fitlm和fitrensemble等回归模型函数这些都是常用工具箱一般装了Matlab都会带。不需要额外安装任何第三方工具箱或mRMR专用包。市面上有一些开源的mRMR工具箱可以下载但说实话自己写一遍核心逻辑也就几十行代码而且能完全按照自己的需求改比如增加并行计算、改成MIQ形式、扩展成条件互信息等。自己掌握原理后面遇到问题才不至于两眼一抹黑。2.2 模块划分四个函数各干一件事整套实现我拆成了四个模块边界非常清晰数据离散化函数把连续特征和目标值转换成离散整数互信息计算函数输入两个离散列向量输出互信息值mRMR主函数负责迭代选择输出特征排序调用与验证脚本加载数据、调用主函数、用回归模型评估效果这样的拆分好处是容易测试。比如你可以单独拿互信息函数跑几个已知的分布检查结果是否符合直觉再交给上层使用。我也强烈建议你按照这个思路组织代码不要把所有逻辑塞进一个超大函数里后面调试和复用都会省很多事。2.3 离散化方案等宽分箱还是等频分箱离散化是回归版mRMR里最容易出错、也最影响结果的环节。常用的两种方案各有优劣等宽分箱是把取值范围等分成若干个区间用linspace生成边界然后映射到1到bins的整数。实现简单但遇到偏态分布的数据时会很吃亏比如大多数样本集中在一个小范围分出来很多箱里没有样本。等频分箱是把数据按从小到大的顺序排列把样本量均匀分配到每个箱子里。这样每个箱子的样本数大致相同联合概率估计更稳定尤其适合长尾数据。但缺点是箱子宽度不一致极端值可能被单独分到一个箱子里。我做了个对比表分箱方式优点缺点适用场景等宽分箱实现简单、语义直观偏态数据容易导致空箱均匀分布、取值区间稳定等频分箱每箱样本量均衡、概率估计稳箱宽不一致、实现稍复杂长尾分布、存在明显离群点代码里默认实现的是等宽分箱因为代码简洁、好理解。如果你的数据偏态明显建议改成等频分箱实现方式是把sort后的序号按比例映射到bins即可后面第5章会给出具体改法。3. 完整Matlab代码与逐段解析3.1 主函数mRMR回归特征选择先给主函数完整代码。这个版本做了两个关键优化一是预计算所有特征两两互信息矩阵避免每一轮迭代都重复计算二是把相关性、冗余度和最终得分都返回出来方便你观察整个选择过程。function [selectedIdx, miTarget, scoresHist] mrmr_regression_fast(X, y, numFeat, bins) % mRMR回归特征选择预计算互信息矩阵版 % 输入 % X n行×d列的连续特征矩阵 % y n×1的连续目标向量 % numFeat 需要选择的特征个数 % bins 离散化分箱数默认10 % 输出 % selectedIdx 按选择先后排序的特征索引1×numFeat % miTarget 每个特征与目标变量的互信息1×d % scoresHist 每一轮的mRMR得分numFeat×d % 用法示例 % idx mrmr_regression_fast(X, y, 20, 10); if nargin 4 || isempty(bins) bins 10; end % 1. 数据离散化 Xd discretizeData(X, bins); yd discretizeData(y, bins); % 2. 预计算所有特征两两互信息矩阵 d size(X, 2); MI zeros(d, d); for i 1:d-1 for j i1:d mij mutualInfo(Xd(:, i), Xd(:, j)); MI(i, j) mij; MI(j, i) mij; end end % 3. 计算每个特征与目标的互信息 miTarget zeros(1, d); for i 1:d miTarget(i) mutualInfo(Xd(:, i), yd); end % 4. 迭代选择 selectedIdx zeros(1, numFeat); scoresHist zeros(numFeat, d); selected []; for t 1:numFeat scores zeros(1, d); for i 1:d % 已选特征直接跳过 if ismember(i, selected) scores(i) -inf; continue; end % 平均冗余度与所有已选特征的互信息均值 if isempty(selected) redundancy 0; else redundancy mean(MI(i, selected)); end % MID形式相关性 - 冗余度 scores(i) miTarget(i) - redundancy; end [bestScore, chosen] max(scores); selected [selected, chosen]; selectedIdx(t) chosen; scoresHist(t, :) scores; fprintf(第%02d轮选入特征 %d得分 %.4f\n, t, chosen, bestScore); end end3.2 辅助函数离散化与互信息计算接下来是两个辅助函数。一个处理连续数据离散化一个计算两个离散向量的互信息。function xd discretizeData(X, bins) % 等宽分箱将连续值映射为1~bins的整数 % X可以是列向量也可以是多列特征矩阵 [n, m] size(X); xd zeros(n, m); for j 1:m col X(:, j); minv min(col); maxv max(col); % 常数特征直接归为第一个箱子 if maxv - minv 1e-12 xd(:, j) 1; continue; end edges linspace(minv, maxv, bins 1); % 稍微扩展右端点保证最大值也被分进最后一个箱子 edges(end) edges(end) 1e-12; xd(:, j) discretize(col, edges); end end function mi mutualInfo(x, y) % 计算两个离散列向量的互信息 % x、y取值均为1~bins的整数 n numel(x); ux unique(x); uy unique(y); % 构建联合分布矩阵 joint zeros(numel(ux), numel(uy)); for k 1:n ix find(ux x(k), 1); iy find(uy y(k), 1); joint(ix, iy) joint(ix, iy) 1; end % 转成概率 pJoint joint / n; pX sum(pJoint, 2); pY sum(pJoint, 1); % 计算互信息只在pJoint0时累加避免log(0) mi 0; for i 1:numel(ux) for j 1:numel(uy) if pJoint(i, j) 0 mi mi pJoint(i, j) * log(pJoint(i, j) / (pX(i) * pY(j))); end end end end3.3 踩过的坑为什么这里要预计算互信息矩阵我第一版写mRMR时是把互信息计算直接嵌在迭代内部的每一轮对每个候选特征都和当前已选特征重新算一次互信息。特征少的时候感觉不出来一旦特征数量超过几百、轮次超过几十运行时间会指数级增长。实测在500维特征上跑50轮慢得让人怀疑人生。后来改成预计算互信息矩阵思路就变成所有特征两两之间的互信息只需要算一次存进一个d×d的对称矩阵里迭代过程中全部走查表。这样计算的复杂度从O(轮次 × 候选特征数 × 已选特征数 × 样本数)降到了O(d² × 样本数)实测500维数据跑50轮时间从原来的几个小时缩到了几十秒。互信息矩阵本身是对称的因为互信息不分方向I(X;Y)和I(Y;X)是同一个值所以我只算了上三角然后对称填进去省掉了一半计算量。这个优化思路也适用于其他特征选择算法比如CFS、FCBF代码结构基本通用。3.4 调用示例把你的数据放进来跑一遍假设你已经准备好了特征矩阵X和连续目标向量y调用方式如下% 载入数据 % 以波士顿房价数据为例最后一列是目标 load(housing.mat); X data(:, 1:end-1); y data(:, end); % 选出15个最重要的特征分箱数用10 selIdx mrmr_regression_fast(X, y, 15, 10); disp(mRMR选出的特征索引按重要性先后排序); disp(selIdx);注意这里有个小细节bins不宜设得太大。如果你的样本量只有几百bins设成20以上每个箱子里的平均样本量会很少联合概率估计的方差会非常大选出来的特征稳定性就差。我用过的经验范围是5~15默认10在大多数场景下都能接受。运行中间会看到每一轮的选入打印比如“第01轮选入特征 7得分 0.2134”。第一轮选出的特征一定是和目标互信息最高的那个因为此时还没有已选特征冗余度为0。后面的轮次会开始出现相关性和冗余度之间的权衡有些单看互信息不是最高的特征反而会被选中因为它们和已选特征重叠少。4. 用回归模型验证特征子集的效果4.1 线性回归快速验证R²和RMSE对比mRMR选出来的特征到底行不行不能凭感觉得放进回归模型里横向对比。最省事的验证方式是用线性回归很快就能看到特征子集和全特征之间的差距。rng(42); cv cvpartition(size(X, 1), HoldOut, 0.2); idxTrain training(cv); idxTest test(cv); % 全特征线性回归 mdlAll fitlm(X(idxTrain, :), y(idxTrain)); yhatAll predict(mdlAll, X(idxTest, :)); r2All 1 - sum((y(idxTest) - yhatAll).^2) / sum((y(idxTest) - mean(y(idxTest))).^2); rmseAll sqrt(mean((y(idxTest) - yhatAll).^2)); % mRMR特征子集线性回归 X_sel X(:, selIdx); mdlSel fitlm(X_sel(idxTrain, :), y(idxTrain)); yhatSel predict(mdlSel, X_sel(idxTest, :)); r2Sel 1 - sum((y(idxTest) - yhatSel).^2) / sum((y(idxTest) - mean(y(idxTest))).^2); rmseSel sqrt(mean((y(idxTest) - yhatSel).^2)); fprintf(全特征R² %.4fRMSE %.4f\n, r2All, rmseAll); fprintf(mRMR子集R² %.4fRMSE %.4f\n, r2Sel, rmseSel);R²代表模型解释了多少比例的方差越接近1越好RMSE是预测误差的均方根越小越好。我实测过很多次mRMR子集在线性回归里的R²往往能追平甚至超过全特征。原因不复杂无关特征和不稳定特征进入线性模型后等效于给模型引入了额外噪声反而把系数估计搞坏了。mRMR这一步等于做了个“降噪”。4.2 不只是线性和随机森林、XGBoost怎么配合线性回归只是验证模版实际项目中mRMR最常用的位置是在树模型和集成模型之前做预筛。比如随机森林回归Matlab自带fitrensemble就可以跑mdlBag fitrensemble(X_sel(idxTrain, :), y(idxTrain), ... Method, Bag, NumLearningCycles, 100); yhatBag predict(mdlBag, X_sel(idxTest, :)); r2Bag 1 - sum((y(idxTest) - yhatBag).^2) / sum((y(idxTest) - mean(y(idxTest))).^2);树模型本身有特征重要性机制听起来好像不需要mRMR但实际项目中两者是互补关系。树模型的特征重要性是训练完成后才知道的而mRMR在训练之前就能帮你把维度降下来。当你有上千个特征、样本量又不大的时候直接训练随机森林很容易过拟合先过一遍mRMR把特征压到几十个树模型的性能通常会提升一大截。至于XGBoost、LightGBM这类梯度提升模型Matlab里要么用第三方接口要么在Python环境里跑。做法是同一个先用Matlab的mRMR选出特征子集导出成CSV再在建模环境里用这份子集训练XGBoost。这种方式特别适合当你的数据预处理在Matlab完成、建模在Python端完成的混合工作流。4.3 该选多少个特征两个实用方法mRMR输出的是特征的优先顺序但具体选多少个合适算法本身不告诉你。我常用的方法有两个第一个是看mRMR得分下降曲线。每一轮得分其实对应了边际收益跑完一轮把scoresHist的主对角线提取出来画个折线图。正常情况下前面几个特征得分很高随后快速下降然后进入一个平缓的长尾。拐点附近的位置就是特征规模的合理上限。第二个更稳是交叉验证网格搜索。用mRMR分别选出1、5、10、15、20个特征在验证集上训练同一个回归模型记录每个规模的R²或RMSE选效果最好的那个规模。这样虽然麻烦一点但结果最可靠。实际操作中我一般先做法一确定一个大致范围再在这个范围内做小规模搜索能省不少时间。5. 常见问题与排查技巧实录5.1 高频报错缺失值、常数特征、类型混乱新手跑mRMR最容易遇到的第一类问题是数据本身有NaN。discretize函数遇到NaN会直接报错或者把NaN单独分到一个奇怪的箱子里。所以加载数据后第一件事就是检查是否有缺失值if any(isnan(X(:))) || any(isnan(y)) error(数据中存在NaN请先处理缺失值); end处理方式要看场景。缺失比例很低的直接用行删除或均值填充比例高的特征建议直接删掉这个特征别等mRMR帮你去筛因为它对缺失值并不友好。第二类问题是常数特征。一个特征如果所有样本取值都相同或者取值变化极小方差接近0互信息天然为0分箱后会全部落到同一个箱子里。代码里我已经做了保护把这种特征自动归到第1箱但它在迭代中会一直得0分只会浪费你的轮次。最好的做法是在跑mRMR之前直接把这些特征过滤掉。第三类问题是数据类型。如果你的数据里有类别型变量比如性别、城市名直接用字符串或者整数编码塞进来含义会错得离谱。整数编码本身就是一种偏序但类别可能是无序的。遇到这种情况要么先做独热编码再跑mRMR要么单独处理这些类别特征不参与互信息计算。5.2 参数调优bins选多少、结果不稳定怎么办bins的选择是回归版mRMR里最核心的参数。太小比如3、4会把很多不同取值揉在一起互信息严重低估太大比如30、50联合分布里大量格子是空的概率估计噪声主导。我自己的建议是样本量在1000以下用8~10个箱1000到10000用10~15个箱超过10000可以用15~20个箱。没有绝对标准但你可以做一个快速实验bins分别取5、10、15看选出的特征排序是否稳定如果差异很大就要考虑是不是样本量太少或者bins不合适。结果不稳定是另一个高频问题。连续特征离散化后信息本身有损失再加上互信息估计的随机性不同bins参数下选出来的特征顺序可能不一样。如果你发现特征排序波动大可以尝试等频分箱替代等宽分箱偏态数据下稳定性会好很多。等频分箱的改法也很简单% 等频分箱实现 function xd discretizeEqualFreq(X, bins) n size(X, 1); xd zeros(n, size(X, 2)); for j 1:size(X, 2) [~, ~, xd(:, j)] histcounts(X(:, j), bins, BinMethod, quantile); end end另外一个不稳定来源是样本划分。如果你想在多次运行中得到比较稳定的结果固定随机种子是基本操作。这不是造假而是让实验可复现。实际项目中也可以多次运行mRMR比如10次每次用不同的bins或分箱方式统计哪些特征被选中的次数最多取交集作为最终特征集合。5.3 运行效率让mRMR在海量特征上跑起来特征很多的时候比如几千上万个特征即便预计算了互信息矩阵两两互信息的计算量也是O(d²)级别非常恐怖。这时我一般用三个手段一是对样本进行子采样。互信息估计在大样本下更准但很多时候几千样本和一万样本估计出来的互信息差距不大所以可以随机抽2000~5000个样本跑mRMR速度能快好几倍结果基本不受影响。二是用并行计算。主函数里的两两互信息循环是天然并行的把for改成parfor前提是已经启用Parallel Computing Toolbox。注意parfor循环里的变量写入方式要调整可以先分别算上三角并存储临时变量跑完再合并。三是最激进的方式先用简单过滤法粗筛一遍。比如先计算每个特征与目标的皮尔逊相关系数取绝对值排名前500的特征再用mRMR在这500个里面精挑。这和第4章说的“两级筛选”思路一样牺牲一点点精确度换来效率的大幅提升。在特征数量过万的项目里我都是这么干的。5.4 中文注释乱码环境问题还是编码问题很多人在自己的电脑上打开代码发现中文注释全部变成乱码这不是代码逻辑的问题是文件编码和Matlab默认编码不匹配。新版本Matlab默认偏好设置里的文件编码通常是UTF-8但老版本或某些Windows中文环境下默认是GBK。解决办法有两个。一个是在Matlab的“预设项 编辑器/调试器 语言”里把文件编码改为UTF-8然后重新打开文件一般能恢复正常。另一个更省事的办法是代码注释统一用英文。如果你经常在多个系统、多台机器之间拷贝代码英文注释能省去很多编码相关的折腾。代码里的打印信息倒是可以随便用中文因为那是运行时输出的字符串不受编码影响。5.5 常见问题速查表我把实际使用中遇到过的问题整理成一个速查表方便你按图索骥问题现象可能原因处理方法discretize函数报错存在NaN或数据包含Inf先清洗数据填充或删除异常值选出的特征都是常数值上的没有预过滤低方差特征提前用方差阈值过滤或增大bins观察两种bins下结果差异太大样本量过少或bins选取不当用等频分箱或者并行多次运行取稳定核心运行时间过长循环里重复计算互信息使用预计算互信息矩阵版本或子采样前几轮选的特征全来自同一类变量编码方式导致伪相关做好独热编码必要时删掉冗余编码列中文注释乱码文件编码不匹配统一UTF-8保存或改用英文注释目标值离散化后信息丢太多分箱数太少适当增加bins但不要超过样本量的1/206. 我的一些实际体会和扩展思路做特征选择这几年mRMR一直是我工具箱里的常备算法。它有各种更花哨的替代品比如基于模型的特征重要性、SHAP值、Lasso系数等但mRMR的好处在于它完全不依赖模型假设在任何建模流程前面都能当第一道筛子。尤其在数据量小、特征维度高的场景下先用mRMR把维度压下来再上线性回归或者树模型效果往往比直接全特征训练好一个档次。代码层面这个实现还能往几个方向扩展。一个是把互信息换成条件互信息也就是CMI版本在序列特征选择时能捕捉更多非线性交互另一个是把分箱步骤替换成核密度估计连续化处理互信息精度会更高但计算量也会上去。如果你的数据里有非常强的周期特征还可以考虑做时间滞后版的特征选择。最后再分享一个小技巧跑mRMR之前先把特征做一次标准化不是为了算法本身而是为了后面接入回归模型时让系数可解释性更好。mRMR的排序不会因为标准化改变太多但后续建模时你会少踩很多数值稳定性的坑。希望这套Matlab实现能帮你在特征选择的路上省下不少折腾的时间。
