简介MATLAB实现的正交信号校正OSC算法脚本面向化学计量学与光谱分析领域的科研人员及定量校正模型开发者。OSC通过滤除红外及近红外光谱中与目标变量无关的信号成分可降低噪声干扰与基线漂移影响从而提升定量校正模型的精度与稳健性对于近红外光谱中常见的背景漂移和多重共线性问题该预处理方法能从源头削弱非目标扰动为后续建立稳定可靠的定量校正模型提供更洁净的输入信号。压缩包内共1个M文件代码体积仅2KB结构紧凑便于直接调用、移植或作为教学示例。CSDN平台已有549人学习浏览该资源适合需要快速上手OSC预处理方法并应用于光谱数据建模的MATLAB用户借助该脚本可完成光谱数据的正交信号校正处理并在此基础上改进偏最小二乘等定量校正模型的预测性能。1. 正交信号校正先分清光谱里的“有用方差”和“与浓度无关的方差”拿到一批近红外光谱建模时校正集误差很低可验证集误差高得离谱或者你把基线拉平、做了一阶导数精度还是不达标。这类问题的根源往往不是随机噪声而是光谱里存在大量的、有结构的系统变异——比如样品颗粒度变化、温度漂移、仪器响应缓慢、光程差异。这些变异与目标浓度无关却会挤占偏最小二乘PLS模型的潜在变量自由度导致模型在训练集上拟合得很漂亮却学错了方向。正交信号校正OSC就是针对这个场景设计的有监督预处理方法它利用浓度向量 y 的信息主动找到光谱 X 中与 y 正交也就是与浓度无关的方向并把它们从原始光谱里滤掉。滤除之后再做 PLS 回归主成分的载荷更干净RMSEC 和 RMSEP 的差距也不再夸张。本文会用 MATLAB 从算法原理讲到可运行的函数再结合红外和近红外光谱数据给出完整的定量校正建模流程以及其中几个最容易踩的坑。2. 正交信号校正的原理为什么“有监督滤除”比直接平滑有效2.1 光谱里的方差不是铁板一块近红外光谱记录的是样品对光的吸收和散射综合结果除了目标化学组分的含量信息还叠加了物理效应和仪器状态。常见的干扰包括光散射引起的基线漂移和整体斜率变化环境温度不同造成的氢键峰位移样品装填密度差异带来的不重复光谱路径仪器光源老化导致的谱段强度缓慢变化。这些干扰的共同点是数值上很有规律甚至占光谱总方差的很大比例但它们的本质是“系统性”的不是随机白噪声。对这类干扰单纯用小波去噪或 Savitzky-Golay 平滑是无效的因为平滑只会让低频率的漂移变得更“干净”地存在。标准正态变量变换SNV和多元散射校正MSC能解决一部分散射问题但是它们只从 X 自身出发完全不知道浓度 y 长什么样所以滤除的方向可能会连真实有用的化学信息一起损失掉。2.2 OSC 与其他预处理方法的定位差异下表把 OSC 和常见近红外预处理方法在“是否使用 y 信息”和“解决什么问题”上做一个对比。方法是否用到 y核心假设典型用途主要局限SNV / MSC否散射变异可以用均值或参考谱估计消除颗粒度、光程引起的基线偏移对化学峰与散射峰重叠严重时效果有限一阶/二阶导数否基线漂移是低频信号分离交叠峰去除慢变基线放大高频噪声需要搭配平滑Savitzky-Golay否信号在局部窗口内可用多项式近似平滑、求导不能处理系统性物理变异OSC是与 y 正交的方差是干扰在建模前定向滤除 X 中与浓度无关的系统变异正交成分数过大会滤掉有效信息造成过拟合从这个表能看出OSC 是所有常见预处理里少见的“有监督”方法。这个属性既是优点也是风险优点是有目的性它可以把与 y 无关的那些“大而强”的系统峰直接消掉而不是绕道走风险是你再拿同一组 y 去做交叉验证选主成分数误差估计容易偏乐观因此模型验证时的样本划分一定不能把测试集信息混进来。2.3 OSC 的数学流程每一步都在做正交投影经典 OSC 算法由 Wold 等人提出后来有不少变种核心步骤是一致的。给定校正集光谱矩阵 Xm 行样本、n 列波长和浓度向量 ym 行假设要滤除 ncomp 个正交成分。第一步中心化。分别对 X 和 y 做均值中心化得到 Xc 和 yc这一步是 PLS 建模的标准预处理避免模型被均值带偏。第二步寻找与 y 相关性最大的方向。对当前剩余光谱矩阵 X_work 和 yc 做一次单主成分的 PLS 回归取出权重向量 w。这个 w 表示的是当前光谱中最能给浓度 y 带来解释的方向。第三步正交化。计算得分向量 t X_work * w然后把 t 向 yc 的正交子空间投影t_orth t - yc * (yc * t) / (yc * yc)这一步是 OCS 的灵魂。因为 t 来自 PLS 的第一权重它天然与 y 强相关但经过上式投影后t_orth 与 yc 的内积为零即不包含任何有关浓度的线性信息。第四步计算载荷并扣除。用 p X_work * t_orth / (t_orth * t_orth) 计算该方向在光谱空间中的载荷然后从 X_work 中减掉这个成分X_work X_work - t_orth * p减完之后循环回到第二步继续寻找下一个正交方向。通常滤除 1 到 3 个正交成分就够用滤得越多剩余光谱的有效信息损失风险也越大。2.4 为什么正交是“金标准”你可能想到一个问题既然要滤除与 y 无关的干扰为什么不直接用 PCA 去掉方差最大的几个主成分PCA 确实能抓到大的变异方向但它不区分这些方向的方差和 y 有没有关系。比如某条样品装填密度的变化方差很大但浓度可能也小幅影响了密度PCA 会把这条主成分直接删掉损失部分浓度信息。OSC 加上了“与 y 正交”这个约束它的意思是只要某个方向和浓度没有任何线性关系不管它方差多大都删掉反过来哪怕方差很小只要与 y 相关就保留下来。这也是 OSC 和后来 OPLS正交偏最小二乘思想一脉相承的原因OPLS 内部实际上是把回归模型分解成预测部分和正交部分而 OSC 是在回归前做正交投影两者对谱图的解释方式很接近。3. 用 MATLAB 实现 OSC最小可用函数与正确投影方法3.1 基于 plsregress 的 osc 函数在正式写函数之前先说一个关键决策点。MATLAB 自带统计和机器学习工具箱里的plsregress函数可以输出 PLS 权重矩阵。OSC 的每一步要提取“当前 X 与 y 相关性最强”的方向稳妥的做法是把这步交给plsregress做而不是自己算X*y。后者在某些极端尺度下不稳定而且 PLS 的权重本身考虑了 X 和 y 的协方差结构收敛更快。下面给出一个可直接保存成osc.m的函数实现。function [X_osc, W, T, P, x_mean, M] osc(X, y, ncomp) % OSC 正交信号校正基于plsregress提取PLS权重 % 输入: % X - 校正集光谱矩阵m行样本n列波长 % y - 校正集浓度列向量m行1列 % ncomp - 需要滤除的正交成分数建议值 1~3 % 输出: % X_osc - 校正后的光谱矩阵与X同尺寸 % W - 每个正交方向使用的权重向量n行ncomp列 % T - 每个正交方向的得分向量m行ncomp列供诊断图使用 % P - 载荷矩阵n行ncomp列 % x_mean- X的列均值用于新样本中心化 % M - 线性投影矩阵n行n列X_osc (X-x_mean)*M x_mean m size(X, 1); x_mean mean(X); y_mean mean(y); Xc X - x_mean; yc y - y_mean; X_work Xc; W zeros(size(X, 2), ncomp); T zeros(m, ncomp); P zeros(size(X, 2), ncomp); M eye(size(X, 2)); for i 1:ncomp % 在当前剩余光谱上做一次单PLS取权重向量的第一列 [~, ~, ~, ~, ~, ~, stats] plsregress(X_work, yc, 1); w stats.W(:, 1); w w / norm(w); % 计算得分并投影到与yc正交的子空间 t X_work * w; t_orth t - yc * ((yc * t) / (yc * yc)); % 最小二乘载荷使X_work近似等于t_orth * p p X_work * t_orth / (t_orth * t_orth); % 扣除该正交成分 X_work X_work - t_orth * p; % 累积线性投影矩阵每步等效右乘 (I - w*p) M M * (eye(size(X, 2)) - w * p); W(:, i) w; T(:, i) t_orth; P(:, i) p; end % 输出时把均值加回保持光谱的物理尺度 X_osc X_work x_mean; end这个函数有几点值得说明。plsregress的第 8 个输出stats.W才是权重矩阵前面几个输出分别是 X 的得分、Y 的得分、载荷等别拿错列。权重向量要做单位化否则后续正交化公式的量纲会乱。正交化公式里的(yc*t) / (yc*yc)是个标量表示 t 在 yc 方向上的分量把它减掉后t_orth和 yc 的内积严格为零。p的计算采用的是普通最小二乘回归系数这样t_orth * p才是 X_work 在该方向上的最佳一维逼近。参数设计上ncomp1在绝大多数近红外定量任务中就有效果因为最大的一块系统干扰占一个方向如果做了第一个正交成分后验证集误差没有下降再尝试ncomp2。超过 3 个通常不是更好的选择而是过拟合的先兆。3.2 新样本怎么处理不能重新算 OSC很多人第一次用 OSC 会犯一个错误把训练集和验证集拼在一起对整批光谱做预处理再划分训练集和验证集。这是最典型的信息泄漏因为验证集样本已经被用于确定权重向量 w 和正交化方向得到的验证误差是假的模型真正上线后表现会大幅缩水。正确做法是把 osc 函数当作一个“在训练集上训练、在验证集上应用”的过程。训练集做完 OSC 后会得到x_mean和投影矩阵 M验证集的光谱不能再用plsregress和 y 去计算权重只能直接用同一个 M 做线性变换。实现方式如下。function X_new_osc osc_apply(X_new, x_mean, M) % OSC新样本变换函数 % 输入: % X_new - 新样本光谱矩阵行数不限 % x_mean - 训练集光谱均值由osc函数输出 % M - 线性投影矩阵由osc函数输出 % 输出: % X_new_osc - 预处理后的光谱 X_new_osc (X_new - x_mean) * M x_mean; end这里 M 的意义是OSC 的每一轮迭代在数学上都是对当前光谱做一次线性变换 X_work X_work - (X_work * w) * p即右乘一个矩阵(I - w*p)。多轮迭代等于连续右乘若干个这样的矩阵所以合并成一个 M。这样新样本不需要知道训练集的 y也不需要重复迭代一次矩阵乘法完成任务。相比需要在预测时重新计算得分并扣减的写法这个方案在验证集上更稳定也不容易引入尺度错误。3.3 用模拟数据验证函数行为为了确认函数写对了可以用一组带强正交干扰的模拟光谱做最小测试。下面的脚本生成 60 个样本400 个波长点浓度向量 y 服从均匀分布。光谱由三个洛伦兹峰组成其中前两个峰与浓度成正比第三个峰与浓度完全无关并且这个无关峰叠加了较大幅度。这样一个理想的测试集里OSC 理论上应该滤除第三个峰。% 测试osc函数的最小复现脚本 rng(42); m 60; n 400; wl linspace(1000, 2400, n); % 浓度向量 y 2 rand(m, 1) * 5; % 两个与浓度相关的化学峰 chem_peak1 8 * exp(-((wl - 1600)./60).^2); chem_peak2 5 * exp(-((wl - 1900)./50).^2); X_chem y * (0.6 * chem_peak1 0.3 * chem_peak2); % 一个与浓度无关的强干扰峰比如水汽吸收的残留 interf_peak 20 * exp(-((wl - 1800)./150).^2); X_interf randn(m, 1) * erase(interf_peak); % 避免变量被复用 X X_chem bsxfun(times, X_interf, interf_peak) 0.05 * randn(m, n);然后调用osc并比较原光谱和校正后光谱在第三个峰位置的方差变化。ncomp 1; [X_osc, W, T, P, x_mean, M] osc(X, y, ncomp); % 看正交得分与浓度的相关系数 corr_T_y abs(corr(T, y)); fprintf(正交得分与y的相关系数: %.4f\n, corr_T_y); % 比较原始光谱和校正光谱在干扰峰中心处的方差 idx_center 401; % 对应1800 nm附近 var_before var(X(:, idx_center)); var_after var(X_osc(:, idx_center)); fprintf(干扰峰方差: 原始 %.2f, 校正后 %.2f\n, var_before, var_after);如果代码正确corr_T_y应该非常小接近 1e-14 量级说明 T 和 y 严格正交。干扰峰附近光谱方差明显变小说明该方向被有效扣除。这个测试脚本的价值在于验证 OSC 数学实现有没有问题而不需要先准备真实数据。4. 实战OSC 改进近红外定量校正模型的完整流程4.1 数据划分要放在 OSC 前面无论用真实数据还是模拟数据模型评估流程都是一样的。先随机划分样本再用训练集去训练 OSC 和 PLS 模型最后把验证集的光谱用训练集得到的参数一次性变换送入模型预测。下面这段代码演示的是完整流程。% 读取数据假设变量为X_all和y_all % load(nir_data.mat); % X_all为nSample x nVary_all为nSample x 1 rng(2026); idx randperm(size(X_all, 1)); ncal 80; Xcal X_all(idx(1:ncal), :); ycal y_all(idx(1:ncal)); Xval X_all(idx(ncal1:end), :); yval y_all(idx(ncal1:end)); % 在训练集上做OSC nosc 2; [Xcal_osc, W, T, P, x_mean, M] osc(Xcal, ycal, nosc); % 验证集用同一个线性投影 Xval_osc osc_apply(Xval, x_mean, M);需要留意数据划分的随机种子。近红外数据经常是按采集时间排序的如果直接取前一部分做校正集后一部分做验证集仪器的漂移和样品温度变化会导致验证误差被高估。随机划分的另一个好处是校正集和验证集的浓度范围接近OSC 的正交化更稳定。4.2 用 PLSR 建模并对比误差指标预处理完成后用plsregress建立定量校正模型并计算常见评价指标。手动计算 RMSEC、RMSEP 和 R² 的方法如下。% 不经过OSC的对照模型 nLV 6; [~, ~, ~, ~, beta_raw] plsregress(Xcal, ycal, nLV); yhat_cal_raw [ones(ncal, 1), Xcal] * beta_raw; yhat_val_raw [ones(size(Xval, 1), 1), Xval] * beta_raw; % 经过OSC的模型 [~, ~, ~, ~, beta_osc] plsregress(Xcal_osc, ycal, nLV); yhat_cal_osc [ones(ncal, 1), Xcal_osc] * beta_osc; yhat_val_osc [ones(size(Xval_osc, 1), 1), Xval_osc] * beta_osc; % 计算均方根误差和决定系数 rmse_cal_raw sqrt(mean((ycal - yhat_cal_raw).^2)); rmse_val_raw sqrt(mean((yval - yhat_val_raw).^2)); rmse_cal_osc sqrt(mean((ycal - yhat_cal_osc).^2)); rmse_val_osc sqrt(mean((yval - yhat_val_osc).^2)); r2_val_raw 1 - sum((yval - yhat_val_raw).^2) / sum((yval - mean(yval)).^2); r2_val_osc 1 - sum((yval - yhat_val_osc).^2) / sum((yval - mean(yval)).^2);得到的四个指标可以汇总成一张对比表大概长这样。模型预处理方式RMSECRMSEPR²(验证)PLS仅均值中心化0.721.650.841PLSOSC(2)均值中心化0.581.020.931这个表反映的是典型趋势OSC 通常让 RMSEC 和 RMSEP 同时下降且 RMSEP 的下降幅度大于 RMSEC。如果看到 RMSEC 大幅下降但 RMSEP 反而上升优先检查是不是 OSC 的正交成分数太多或者验证集样本没有用训练集投影处理。4.3 两个参数互相影响nosc 和 nLV这节有一个常见误区认为做 OSC 是为了减少 PLS 主成分数 nLV。实际上不一定。OSC 可能滤掉一个大干扰但剩余光谱里可能还需要相同数量的潜在变量来表述浓度信息nLV 的减少与否取决于干扰占据的是否是前几个主成分。真实实验中更常见的现象是nLV 不变或减少 1 个预测精度明显提升如果 nLV 减少太多说明 OSC 把部分有效信息当成正交成分滤掉了。调参顺序建议先定 nLV对原始光谱做 10 折交叉验证找到 RMSEP 最低的 nLV。固定这个 nLV 后再尝试 nosc 1, 2, 3用验证集 RMSEP 做选择。注意交叉验证不能在完整数据集上做必须在训练集内部再做一次嵌套划分否则会低估验证误差。5. 把 OSC 用好三个验证技巧和一个警惕5.1 画正交得分与 y 的散点图OSC 函数的输出 T 是得分矩阵每一列代表一个正交方向。正确实现的 OSC 里T 的每一列与 y 的相关系数都应当接近零。画散点图的用途在这里如果发现某个正交得分和 y 存在明显线性趋势说明该成分没有真正正交可能是正交化公式里 yc 没有中心化或者plsregress的权重取的列不对。代码很简单。figure; plot(y, T(:, 1), o); xlabel(浓度 y); ylabel(正交得分 t_1); title(正交得分与浓度关系诊断);图上如果看到点均匀散布成水平带说明正交性合格如果看到斜向的带状趋势回到osc.m检查第三步的yc*t是否用的中心化后的 y。5.2 与 SG 平滑的组合顺序OSC 和 Savitzky-Golay 平滑组合时应先做 SG 再做 OSC。原因是 SG 平滑会降低高频噪声让 OSC 计算权重向量时更专注于低频的系统变异反过来如果先做 OSC残差中的噪声在 SG 平滑后可能会产生新的伪峰干扰后续建模。组合参数上SG 窗口一般取 11 到 15 个数据点多项式阶数 2 或 3。窗口太大会抹掉细窄的化学峰太小则噪声压制不足。5.3 警惕“滤过头”用 Q 残差和 Hotelling T² 判断随着 nosc 的增加剩余光谱的方差会不断下降PLS 的校正误差也随之下降但验证误差往往先降后升。这个拐点就是滤除过度出现的信号。另一种靠谱的判断工具是 Hotelling T² 和 Q 残差在做完 OSC 后将校正集光谱投影到 PLS 主成分空间中计算每个样本的 Q 残差如果某个样本的 Q 残差突然变得异常大说明该样本的有效化学信息可能被正交成分过度扣除。实际项目里我发现 nosc 大于 3 后遇到这种情况的频率显著增高所以我的默认配置是 nosc1只有验证误差没有改善时才加到 2极少使用 3 以上。5.4 一个容易被忽略的细节预测时均值中心化的顺序osc_apply里先减x_mean再乘 M最后加回x_mean。这个步骤常被简化成直接乘 M导致预测偏差。原因在于 M 是在中心化后的空间里构造的它不包含均值还原能力。如果你打算把 OSC 嵌入到工业化实时预测流程中建议把均值中心化、OSC 投影、PLS 回归的 beta 系数合并成一个大矩阵使新样本只需要做一次矩阵乘加运算。合并后的系数矩阵等于M * beta_p其中beta_p是 PLS 回归系数向量这样在线部署时对每一条光谱的计算成本可以忽略不计。本文还有配套的精品资源点击获取