简介iVISSA是一套面向光谱特征波段筛选的算法实现与示例资源适合遥感、环境科学、农业等领域的科研人员和从业者使用。光谱数据往往包含大量冗余波段该资源聚焦于从高维光谱中挑选信息价值最高的特征波段以提升定量预测模型的精度与效率。压缩包共12个文件以9个.m脚本为主涵盖主程序、预处理、交叉验证、预测等核心函数同时附有2个.mat示例数据如大豆水分光谱数据和1份license.txt许可说明。其中ivissa.m为核心筛选算法plsnipals.m、plscvfold.m提供偏最小二乘回归与交叉验证支撑pretreat.m实现去噪、平滑等预处理便于直接运行Example1.m和Example2.m快速复现完整流程。资源包大小仅157KB结构紧凑。已有1682人学习下载适合具备一定光谱数据处理基础的读者通过源码理解波段选择的每个环节并迁移到自身的建模任务中减少数据冗余、提升预测准确性。1. iVISSA 是什么光谱特征波段筛选从“肉眼挑峰”进入迭代收缩在近红外和拉曼光谱的实际项目里一条成品样本谱通常包含几百甚至上千个波长点。很多人第一步就是把这些点全部塞进 PLS 回归训练集 R² 很好看换到新批次样品却一路跑偏。这个问题的常见根源不是仪器漂移而是光谱特征波段筛选没有做模型里混进了太多与目标性质无关的波长回归系数被噪声稀释换样品后自然站不住。iVISSA 正是专门解决光谱特征波段筛选这一类问题的算法它用重复采样和概率收缩从全谱中逐步筛出真正有效的小段波段兼顾预测精度与物理解释。适合被模型迁移性折磨的应用工程师、想系统做变量选择的研究者以及要把筛选流程集成进建模套件的开发人员。2. 做一个能落地的 iVISSA概率收缩逻辑、目标函数与选型理由2.1 VISSA 的概率收缩循环从 0.5 到 0.05/0.95 的变量投票要理解 iVISSA先看它的母体 VISSA 是怎么工作的。把光谱矩阵 X 看成 n 个样本、p 个波长的二维数组每个波长就是一个候选变量。VISSA 的核心思想不是穷举所有变量组合而是维护一个长度为 p 的概率向量初始时每个波长的入选概率都是 0.5意思是它们有均等机会进入下一轮随机子模型。每一代迭代可以拆成四步按当前概率向量做伯努利采样生成 K 个变量子集对每个子集建立 PLS 模型并做交叉验证得到 RMSECV把所有子模型按 RMSECV 排序取前 10% 作为优秀子模型统计这些优秀子模型里每个波长被选中的频率频率高的波长概率上调频率低的概率下调。下一轮迭代就基于更新后的概率向量继续采样。这个循环的精妙之处在于它把“哪个变量子集预测误差低”这个问题转化成了“哪些变量在低误差子模型里高频出现”。因为每一轮只有部分波长参与建模相当于在反复做带权重的组合实验。一个波长如果总是出现在优秀子模型里说明它对预测目标的贡献是稳定且独立的反之一个波长虽然单看可能与 y 相关但一旦和其他变量组合预测误差就上升它的概率就会逐步被压低。经过几十轮迭代真正重要的变量概率向 0.95 靠拢无关变量向 0.05 靠拢最终用阈值切一刀就得到了特征波段集。2.2 interval 还是 improvediVISSA 在光谱落地的两种主流形态iVISSA 在不同论文和仓库里可能指两种改进方向。一种是 improved VISSA改进点放在概率更新策略、目标函数或初始采样分布上比如用贝叶斯信息准则替换纯 RMSECV或者引入加权自举采样来加速收敛。另一种是 interval VISSA把搜索单位从单个变量提升到连续区间。我在光谱项目里几乎都优先按 interval VISSA 来理解理由很朴素光谱是连续信号一个特征谱区本来就对应某种化学键振动或分子结构如果算法筛出来的结果是 15 个互不相邻的孤立波长点化学上很难交代换一台同型号仪器也未必复现。interval 版本的标准做法是两级搜索。先把连续波长按固定长度切成窗口比如每 20 个波长为一段在窗口层面做一次 VISSA筛出重要的谱区段然后只把这些重要区段内的波长展开在变量尺度上再做一次 VISSA得到最终的精筛波段。这样既保留了 VISSA 用模型误差驱动筛选的优点又强制加入了光谱连续性约束结果更接近工程师能理解和向领导汇报的形式。后面第 3 章的代码就是按这个两级结构写的。2.3 为什么要选 iVISSA和 SPA、CARS、UVE 的对比光谱特征筛选的算法不少各自适用场景差别很大。SPA连续投影算法通过寻找共线性最小的投影方向来选变量速度快、确定性高但输出是离散点物理意义弱适合作为粗筛工具而不是最终方案。CARS竞争自适应重加权采样用回归系数绝对值做蒙特卡洛淘汰在化学计量学里非常流行但它依赖的阈值参数多对光谱噪声敏感结果稳定性常被吐槽。UVE无信息变量消除把人工噪声变量混入矩阵做对比擅长快速剔除完全无用的波长却不擅长从剩余变量中挑出最优子集。iVISSA 的优势在于目标函数直接就是交叉验证后的 RMSECV每一步都在为最终预测误差服务而概率收缩机制又能跳出局部最优不会像贪婪式搜索那样一条路走到黑。下面是常用方法的选择参考。方法搜索单位主要机制结果形态常见短板SPA单个变量投影寻优离散波长点忽略相位邻域信息CARS单个变量蒙特卡洛淘汰稀疏变量集参数多阈值敏感UVE单个变量噪声对比无信息变量剔除不能直接给最优子集VISSA单个变量概率收缩采样稀疏变量集单点搜索连续性弱iVISSA窗口变量两级概率收缩连续谱区计算量较大2.4 目标函数与 EDF把“筛得好”从感觉变成可计算VISSA 这类算法最大的隐患是过拟合。变量选得越多模型自由度越大交叉验证误差常常显得越低但这个低误差是模型复杂度堆出来的不是真实预测能力的提升。VISSA 论文中引入了 EDF等效自由度来约束这个过程EDF 同时考虑变量数和 PLS 潜变量数对模型复杂度的影响让变量筛选不只优化 RMSECV而是优化“RMSECV 与复杂度之间的平衡”。实际代码实现里我一般用近似方案代替严格 EDF当两个随机子模型的 RMSECV 相差不到 2% 时强制选择变量数量更少的那个在保存全局最优子集时同时记录变量数和 PLS 成分数把二者之和作为复杂度代理指标。这个近似和论文里精确 EDF 的决策差异通常很小但理解和修改成本低得多。如果后续有严格对标论文的需求再按 EDF 定义重算一遍损失函数即可。3. 用 Python 跑通 iVISSA预处理、核心循环到区间版扩展3.1 数据准备与预处理顺序为什么筛选前必须先校正光谱数据直接进筛选经常翻车。近红外光谱受固体颗粒大小、装样松紧、环境湿度影响基线漂移和光散射差异常常比目标化学成分的变化还大。如果不先做预处理iVISSA 极有可能把基线漂移区域筛成“重要特征”因为那些区域在区分训练集样本时确实有效但换一批物理形态不同的样品就完全失效。我的固定套路是先做 SNV标准正态变量变换或 MSC多元散射校正视数据散射程度决定然后视峰重叠情况决定是否做一阶导或二阶导最后再进变量筛选。SNV 实现很简单每个样本减去自身光谱均值再除以自身光谱标准差。import numpy as np def snv(x): # x: (n_samples, n_wavelengths)逐样本做标准正态变换 mean x.mean(axis1, keepdimsTrue) std x.std(axis1, keepdimsTrue, ddof1) std[std 1e-10] 1.0 # 防止常值光谱除零 return (x - mean) / std这段代码里keepdimsTrue是为了让均值向量的形状能和原始矩阵对齐ddof1用的是样本标准差和多数化学计量学软件保持一致。注意预处理参数必须在训练集上计算测试集直接套用训练集的均值和标准差否则会产生数据泄漏。后面第 6 章会再强调这一点。3.2 VISSA 最小实现概率采样、PLS 评价与收缩主循环下面是一个教学级但完整可运行的 VISSA 最小实现。它包含三个部分一个计算 RMSECV 的 PLS 交叉验证函数一个执行概率收缩的主循环和一个固定随机种子的调用示例。生产环境里我会把子模型数调到 1000 以上但演示代码先用较小数值保证能跑通。import numpy as np from sklearn.cross_decomposition import PLSRegression from sklearn.model_selection import KFold def rmsecv_pls(X, y, mask, n_components8, cv5, seed42): # 只用 mask 选中的波长列建立 PLS返回 5 折交叉验证 RMSECV Xs X[:, mask] kf KFold(n_splitscv, shuffleTrue, random_stateseed) errs [] n_comp min(n_components, Xs.shape[1]) # 潜变量数不能超过变量数 for train_idx, val_idx in kf.split(Xs): pls PLSRegression(n_componentsn_comp, scaleTrue, max_iter2000) pls.fit(Xs[train_idx], y[train_idx]) pred pls.predict(Xs[val_idx]).ravel() errs.append((y[val_idx] - pred) ** 2) return float(np.sqrt(np.mean(np.concatenate(errs)))) def vissa(X, y, max_iter30, n_submodels400, cv5, pls_comp8, thres0.7, step0.1, elite_ratio0.1): n_vars X.shape[1] prob np.full(n_vars, 0.5) # 所有波长初始入选概率相同 best_rmse np.inf best_mask None for it in range(max_iter): models [] for _ in range(n_submodels): mask np.random.binomial(1, prob, sizen_vars).astype(bool) while mask.sum() 2: # 保证每个子模型至少选 2 个变量 mask np.random.binomial(1, prob, sizen_vars).astype(bool) score rmsecv_pls(X, y, mask, n_componentspls_comp, cvcv) models.append((mask, score)) models.sort(keylambda x: x[1]) # 按 RMSECV 升序 keep max(1, int(n_submodels * elite_ratio)) elite [m for m, _ in models[:keep]] freq np.vstack(elite).astype(int).mean(axis0) # 优秀子模型中的入选频率 prob step * ((freq 0.6) - (freq 0.4)) # 概率上调或下调 prob np.clip(prob, 0.05, 0.95) # 防止概率收敛到 0 或 1 再不出来 if models[0][1] best_rmse: best_rmse models[0][1] best_mask models[0][0].copy() n_active int((prob thres).sum()) print(fiter {it1}: best RMSECV{best_rmse:.4f}, active{n_active}) final_mask prob thres return final_mask, prob, best_rmse这个实现的逻辑和论文版本是一致的。每次迭代生成 400 个随机子模型每个子模型都做一次 5 折 PLS 交叉验证用 400 个 RMSECV 的排序结果来决定波长概率的涨跌。freq 0.6和freq 0.4是两个经验阈值对应概率更新方向step0.1控制每轮变化的幅度。阈值thres0.7是最终决定入选波段的概率门限设低一点得到更宽的波段集设高一点得到更稀疏的波段集。调用方式如下np.random.seed(42) X_snv snv(X) sel, prob, rmse vissa(X_snv, y, max_iter20, n_submodels200) print(selected wavelength indices:, np.where(sel)[0])注意这里只固定了全局随机种子实际项目还需要在交叉验证、数据划分等环节统一固定种子否则每次跑出的波段集都会有差异。运行时间方面700 个波长、200 个子模型、20 轮迭代在普通笔记本上大约几分钟内可以完成如果换成 1000 个子模型和 50 轮时间会线性上升建议先用降采样或窗口化把变量数压下来。3.3 iVISSA 的区间版实现先把波长装进窗口再收缩现在把 VISSA 升级成 interval 形态。第一级以窗口为单位做变量收缩第二级在选中窗口内部做精细选择。窗口特征的构造常见做法是取窗口内所有波长的平均值这样每个窗口得到一个综合信号强度相当于把光谱数据降维后再做变量选择。def window_means(X, win20): n, p X.shape n_win p // win return np.column_stack([X[:, i*win:(i1)*win].mean(axis1) for i in range(n_win)]) # 第一级窗口级 VISSA Xw window_means(X_snv, win20) sel_win, _, _ vissa(Xw, y, max_iter10, n_submodels200) # 将选中窗口展开为原始波长候选集 candidate_idx np.concatenate([ np.arange(i*20, (i1)*20) for i in np.where(sel_win)[0] ]) # 第二级在候选波长内部做变量级 VISSA Xc X_snv[:, candidate_idx] sel_var, _, _ vissa(Xc, y, max_iter15, n_submodels300) final_wave candidate_idx[sel_var] print(final wave indices:, final_wave)第一级 VISSA 里max_iter10就够用了因为窗口数量远小于波长数量收敛得快第二级在候选区间内做变量级筛选迭代次数可以稍多。窗口大小 20 是经验值具体要根据光谱分辨率来定如果化学峰半宽只有 5 个波长点窗口取 10 更合适如果是宽缓的蛋白质峰窗口可以放到 30。选择窗口的另一个好处是减少最终结果的离散程度筛出来的波段天然是连续区间方便后续做谱峰归属。3.4 可视化输出把筛选结果落回原始谱线筛选完一定要把波段画回原始光谱上用眼睛确认结果是否合理。这个步骤能发现很多数字指标看不出的问题比如筛出了强烈的噪声区、选在基线陡坡上、或者漏掉了肉眼可见的强吸收峰。import matplotlib.pyplot as plt wavelength np.arange(X.shape[1]) plt.plot(wavelength, X_snv.mean(axis0), lw1, colorblack) for idx in final_wave: plt.axvspan(idx, idx1, colorred, alpha0.3) plt.xlabel(wavelength index) plt.ylabel(mean SNV intensity) plt.show()图中红色竖带就是最终保留的谱区。如果红带恰好覆盖了已知的目标成分吸收峰说明筛选结果有物理依据如果红带集中在噪声区或光谱边缘即使 RMSECV 很好看这个结果也值得怀疑需要回头检查预处理和参数设置。4. iVISSA 的 4 个关键参数迭代次数、子模型数、概率步长与窗口尺度4.1 max_iter用活跃变量数判断收敛而不是只看迭代到上限max_iter 不等于“跑满就是收敛”。我见过不少人设了 50 次迭代跑完直接拿结果但中间过程完全不看。更稳妥的做法是在主循环里打印active数量也就是当前概率超过最终阈值 thres 的变量个数。如果 active 数量连续 5 轮到 10 轮不再变化说明概率向量已经基本稳定可以提前终止如果到最后一轮还在上下跳动说明概率步长太大或者子模型数太少加迭代次数没有意义要先降 step 或者加 n_submodels。经验上变量数 700 左右时30 到 50 轮足够收敛如果做了窗口化第一级通常 10 轮就稳定了。4.2 n_submodels子模型少是噪声多则算力爆炸每一轮的随机子模型数直接决定统计可信度。n_submodels200 时优秀子模型只有 20 个频率统计的方差很大概率更新容易受偶然组合影响升到 1000 时优秀子模型有 100 个频率估计更稳但耗时线性增长。我的一般习惯是样本数低于 100 时用 500 个子模型样本数几百到上千时用 1000。如果项目对运行时间敏感先做 UVE 或 SPA 粗筛把变量数压到 100 到 200 个再跑 VISSA这样 400 个子模型就够用。计算时间与子模型数和迭代次数成正比优化方向永远是先降变量数而不是先降子模型数。4.3 prob_step 与 threshold概率更新步长决定收敛质量step 控制每轮概率调整的幅度。step0.1 是常见默认值意思是频率高的变量每轮概率上调 10%大约 4 到 5 轮就能从 0.5 走到 0.9 附近。step 太大容易造成振荡比如某轮由于随机抽样偶然选中了一批变量概率被猛抬下一轮这批变量又以高概率进入大量子模型结果优秀子模型全被它们占据真实重要的变量反而被压下去。如果观察到 active 数量反复上下跳把 step 降到 0.05 再跑一轮。threshold 决定最终入选变量数量一般取 0.7。想要更稀疏的特征集就用 0.85想保守一点就用 0.6。注意 threshold 设低时入选变量变多模型复杂度上升需要同步增加 EDF 惩罚或减少 PLS 成分数。4.4 window_size 与 n_components两个隐含参数决定筛选的物理尺度window_size 不是 VISSA 主循环里的参数而是 interval 版的第一级切窗宽度。它决定了筛选结果的“分辨率”窗宽太细比如 5 个点结果还是接近单点筛选连续性优势消失窗宽太粗比如 50 个点一个重要谱区里混进了不重要的波长第二级筛选压力大。我一般先观察原始光谱中最窄吸收峰的半峰宽用半峰宽的两倍左右作为窗宽。比如某体系的峰半高宽是 8 个波长点窗口就取 15 到 20。PLS 成分数 n_components 同样重要无论哪个子模型都带着同一套成分数参与交叉验证如果成分数设得过高每个子模型几乎都能拟合得很好排序就失去区分度设得过低则所有模型都欠拟合。常见做法是先用全谱 PLS 找 RMSECV 曲线的最低点再用对应成分数作为 VISSA 的固定值。参数推荐范围主要影响调参要点max_iter30 - 100概率是否收敛看 active 数连续 5 轮不变即可停n_submodels500 - 1000统计稳定性小样本 500大样本 1000step0.05 - 0.1收敛速度与振荡active 数跳动时降 stepthres0.6 - 0.9最终入选变量数结合 EDF 惩罚调整window_size10 - 30 点谱区连续粒度按半峰宽估算n_components5 - 10PLS 潜变量数用全谱 RMSCV 曲线初选5. iVISSA 实战避坑5 个高频翻车场景与排查路径5.1 同一条数据每次跑出不同波段集现象代码没改数据没换连续运行两次 iVISSA筛选出的波段索引对不上有时差几个点有时差一整段。原因VISSA 本身依赖随机采样np.random.binomial每次用不同的随机数流交叉验证里的数据划分也是随机的两处随机性叠加结果自然不同。解决在整个流程入口固定随机种子常见做法是np.random.seed(42)同时给 KFold 的random_state设固定值。更严格的项目还要把训练集、验证集划分顺序存成文件保证每次运行完全可复现。如果固定种子后结果仍然不稳定说明 n_submodels 太少加大到 1000 再跑。5.2 筛选结果全是一条条孤立单峰现象选出来的波长点彼此不相邻在光谱图上呈散点状完全看不出谱区形状化学归属无从谈起。原因直接用变量级 VISSA忽略光谱连续性。算法只关心预测误差不关心波长之间是否有物理关联所以在信息冗余的光谱里选出互相独立的单点很容易。解决改用 interval 版 iVISSA先窗口后变量两级搜索从结构上保证结果连续。如果坚持用变量级也可以在每轮概率更新后对概率向量做一次滑动平均平滑让相邻波长的概率互相影响但这会引入额外参数不如切窗口直接。5.3 训练组 RMSECV 很低外部验证组 R² 却崩了现象VISSA 跑完训练集交叉验证 RMSE 漂亮得惊人模型拿到外部测试集上预测误差翻倍R² 掉到 0.5 以下。原因筛选过程只优化了 RMSECV没有对变量数量做惩罚。变量越多模型自由度越大交叉验证误差被复杂性压低。尤其在样本量小于变量数的高维场景这种过拟合几乎是必然的。解决引入 EDF 惩罚或者按第 2.4 节的近似方案两个子模型 RMSE 相差 2% 以内时选变量更少的。同时在外部验证时使用重复交叉验证比如 5 次 5 折报告 RMSE 的标准差而不是单次数值。5.4 换一种预处理筛选出的波段区完全改变现象同一批数据用 SNV 预处理筛出 A、B 两个区域换成 MSC 后筛出的区域变成 C、D两部分几乎不重叠。原因预处理改变的是光谱的形状和方差分布而 VISSA 是模型驱动方法哪里的信号与 y 相关性最强哪里就是优势变量。不同预处理放大的信息不同结果自然不同这不一定是算法的 bug。解决先花时间确定一条固定的预处理链不轻易切换。比较不同预处理方案时不要只比较最终筛选波段的重合度而要比较筛选后 PLS 模型在外部测试集上的 RMSEP。只要 RMSEP 在同一水平波段不同是可以接受的。5.5 小样本高维数据跑不动或不敢信现象样本只有 50 个波长却有 2000 个iVISSA 跑得非常慢而且每一轮生成的子模型训练集只有几十个样本PLS 很容易在随机组合上拟合出极低的 CV 误差。原因高维小样本下随机子模型数量虽多但每次可用的信息极其有限统计噪声被样本量放大筛选结果缺乏可信度。解决先用 SPA 或 UVE 把变量数粗筛到 100 以内再跑 iVISSA 细筛。如果仍然担心稳定性把整个筛选流程重复 20 次只保留每次都被选中的变量。这样虽然会丢失一些边缘变量但留下的是最可靠的波段集。6. iVISSA 筛完之后重复投票、SPA 串联与外部验证6.1 二十次重复投票让随机收敛变成稳定集筛选结果稳定比筛选结果完美更重要。我常用的收尾操作是把 iVISSA 重复跑 20 次每次用不同的随机种子记录每个波长被选中的次数只保留出现频率不低于 15 次的变量作为最终特征集。这样做的本质是在随机搜索算法上做了一次集成抵消单次运行时的抽样波动。代码实现不复杂只需要被调主循环外面再套一层循环。votes np.zeros(X.shape[1]) for i in range(20): np.random.seed(100 i) sel, _, _ vissa(X_snv, y, max_iter30, n_submodels300) votes[sel] 1 stable votes 156.2 与 SPA 或 UVE 串联先粗筛再精筛的常规组合变量数特别多或样本量很小的时候我会先用 SPA 选 30 到 50 个候选波长然后把 VISSA 的搜索空间限制在这 50 个候选点上。这个组合既给了 VISSA 一个合理起点又保留了概率收缩对组合优化的能力。类似地UVE 可以先剔除大约 70% 的无信息变量再用 iVISSA 精筛。串联后总耗时往往比直接跑全谱 iVISSA 短得多结果稳定度也更高。6.3 用外部测试集验证筛选结果RMSEP 与 R² 的对比标准筛选和建模都只能在训练集上做测试集必须留到最后一刻才接触。这里最容易犯的细节错误是预处理泄漏SNV 的均值和标准差必须在训练集上算然后套用到测试集而不是测试集单独做一次 SNV否则等于测试信息提前进入了模型。代码上要把筛选、预处理、建模三项都封装成“训练专用函数”测试集只走一遍预测流程。from sklearn.model_selection import train_test_split X_tr, X_te, y_tr, y_te train_test_split(X, y, test_size0.2, random_state42) X_tr_s snv(X_tr) # 训练阶段筛选 建模 sel, _, _ vissa(X_tr_s, y_tr, max_iter20, n_submodels300) pls PLSRegression(n_components5).fit(X_tr_s[:, sel], y_tr) # 测试阶段套用训练集 SNV 均值/标准差 X_te_s (X_te - X_tr.mean(axis0)) / X_tr.std(axis0) pred_te pls.predict(X_te_s[:, sel]).ravel()这个三段式流程我已经形成习惯先预处理再筛选最后建模每一段都只在训练侧拟合参数。换过太多数据集后发现凡是外部验证崩掉的案例往回查八成都是测试侧悄悄做了“独立预处理”。筛选稳定性和数据泄漏问题处理干净iVISSA 这类算法才能真正成为光谱建模里的可复现工具这也是我最后想提醒你的一个原则特征筛选不是一次性的魔法而是一条需要固定种子、固定流程、反复验证的工程链路。希望帮到你。本文还有配套的精品资源点击获取