递推最小二乘算法详解:原理、MATLAB实现与实验报告指南
简介:递推最小二乘(RLS)算法是动态系统参数在线估计的经典方法,本资源以Python实现为核心,配套Word版详细实验报告,面向自动化、信号处理及相关专业学生与工程师,适合课程实验、课设或工程预研。压缩包共3个文件,.py源码包含初始化、预测、更新、遗忘因子调节等关键步骤,结构清晰且可直接运行;.txt文件对程序使用要点和运行环境作简要说明;.doc报告系统讲解算法背景、递推公式推导、矩阵更新细节、代码逐段解析以及实验结果分析,同时结合遗忘因子对估计精度的影响进行讨论,形成从原理到实践的完整闭环,既适合教学演示,也适合工程验证。资源包仅316KB,轻量实用,便于按需阅读和调试。目前已有466人学习下载,尤其适合希望在短时间内理解RLS机制并动手验证的读者,借助报告中的曲线和结论可快速校验算法效果,是学习自适应滤波与在线辨识的有力工具。 拿到“递推最小二乘算法程序代码及word版详细实验报告.rar”这种压缩包时,我的第一反应不是赶紧解压,而是先问一句:递推最小二乘算法到底是用来干什么的,我现在能不能用自己的话讲清楚?作为系统辨识和自适应控制课程里的常驻角色,递推最小二乘(Recursive Least Squares, RLS)解决的问题非常实际:系统模型结构已知,但参数未知,数据还在源源不断到达,我们必须一边接收新观测值,一边在线刷新对参数的估计。这篇内容我不打算只替你盘点资源包,而是从实践角度完整梳理RLS的原理推导、程序代码模块、实验报告写法,以及跑这类程序时容易踩的坑。
1. 拿到这份rar资源后,先别急着解压跑代码
1.1 压缩包里大概率有什么
这种标题的资源包,内容构成一般比较固定:一个主程序文件(多数是MATLAB的.m文件,有时是Python脚本)、一份word版实验报告、可能还带几个数据文件或图片。rar只是一个打包容器,真正值钱的是里面代码的组织方式和报告里对公式、步骤、结果图的完整记录。
我见过太多人拿到这类资源后,第一件事就是解压、运行、截图、导出PDF,然后当成自己的实验报告交上去。这种用法其实特别浪费。代码能不能跑通是一个层面,能不能讲清楚递推式里每个变量为什么这样更新,是另一个层面。资源包里的代码和报告,应该被当成一份“参考答案”,而不是“作业本体”。
1.2 三分钟判断这个资源和你是否匹配
解压之前,先花三分钟检查几件事:
- 代码是MATLAB还是Python?如果是MATLAB,你的环境里是否有MATLAB或者Octave可以打开;
- 实验报告里的模型结构是什么?最常见的是一阶/二阶ARX模型,也有可能是差分方程形式;
- 报告和代码里的符号约定,和你课程教材是否一致?比如有人用θ表示参数向量,有人用a、b分开表示;
- 代码里是否含有数据生成部分,还是只读取现成数据文件。
这几个问题不确认清楚,直接跑代码很容易出现“代码报错找半天原因,结果是版本不兼容”或者“报告结论和代码结果对不上”的情况。宁可先花几分钟把目录结构和符号约定理清楚,也别盲跑。
1.3 一个高效的学习顺序建议
如果你是想通过这份资源真正学会递推最小二乘,我的建议顺序是这样:先通读实验报告里的原理部分,把算法流程框架建立起来;再运行一遍示例代码,观察参数估计曲线的收敛过程;然后逐行阅读主循环代码,把每一行和递推公式对应起来;最后把仿真对象换成你自己的系统,或者换一组参数重新辨识。
这个顺序的核心逻辑是“先见森林,再见树木”。先知道算法要解决什么问题、流程是什么模样,再去扣代码里的细节,才不会迷失在一堆矩阵运算里。
2. 递推最小二乘的原理:从“一次算完”到“来一组修一组”
2.1 先回顾一下离线最小二乘
递推最小二乘并不是一个全新的独立算法,它就是离线最小二乘的在线版本。离线情况下,假设我们已经收集了N组输入输出数据,模型可以写成:
Y = Φθ + e
其中Y是输出向量,Φ是回归矩阵,θ是待辨识的参数向量,e是残差。最小二乘的核心思路是让残差平方和最小:
J = (Y - Φθ)^T (Y - Φθ)
对θ求导并令导数为零,就得到经典的解:
θ_hat = (Φ^T Φ)^(-1) Φ^T Y
这个公式本身很简单,问题在于它是一次性计算。每一组新数据到来时,如果都要重新构建Φ和Y,再求一次矩阵逆,计算开销会随数据量增加越来越大。更重要的是,在实时控制、在线监测这类场景里,我们不可能等所有数据采集完再解方程,数据是不断产生的,处理器必须在每个采样周期内完成参数更新。这正是递推最小二乘存在的理由。
2.2 递推公式是怎么推出来的
递推的思路是:当新样本(φ(k+1), y(k+1))到来时,利用上一时刻已有的估计结果,加一个修正项得到新估计,而不是重新求解整个最小二乘问题。
为了说清楚这个修正项,先定义P(k) = (Φ_k^T Φ_k)^(-1),它相当于当前时刻的“协方差矩阵”,衡量的是对当前估计的置信程度。新数据到来后,P的逆矩阵可以写成:
P(k+1)^(-1) = P(k)^(-1) + φ(k+1)φ(k+1)^T
利用矩阵求逆引理,可以从这个式子推出三行核心递推公式:
增益矩阵:K(k+1) = P(k)φ(k+1) / (1 + φ(k+1)^T P(k) φ(k+1))
参数更新:θ_hat(k+1) = θ_hat(k) + K(k+1)[y(k+1) - φ(k+1)^T θ_hat(k)]
协方差更新:P(k+1) = (I - K(k+1)φ(k+1)^T)P(k)
括号里的欧几里得向量形式就是整个RLS的骨架。如果只用一句话说清楚:新估计 = 旧估计 + 增益 × 预测误差。预测误差是新观测值和旧模型预测值之间的差值,增益K决定这个差值有多少比例被用来修正参数。
我通常用一个导航类比的例子帮助理解:旧参数估计就像原本规划的路线,新观测数据就像实时路况信息。如果路况信息质量很高且当前路线可信度不足,增益K就大,路线会大幅修正;如果当前路线已经很可靠,K就小,新信息只是微调。P矩阵就是这个“可信度”的量化载体。
2.3 遗忘因子是怎么进来的
普通递推最小二乘给所有历史数据相同的权重,这适合时不变系统。但工程中很多对象会随工况变化,比如电机温升后电阻值改变,或者飞行器气动参数随马赫数变化。这种情况下,我们希望算法“忘记”过于陈旧的数据,更相信最近时刻的观测。
做法是在损失函数里引入指数遗忘因子λ。每增加一个新数据,旧数据的权重都会乘以λ,整体形成指数衰减。递推公式变成:
K = P(k)φ(k+1) / (λ + φ(k+1)^T P(k) φ(k+1))
P(k+1) = (I - Kφ(k+1)^T)P(k) / λ
当λ=1时,退化为普通递推最小二乘,所有历史数据权重相同;λ越小,旧数据被遗忘得越快,算法跟踪时变参数的能力越强,但代价是对噪声更敏感,稳态波动更大。
实际调试时,λ的取值范围通常在0.95到1之间。0.98是个不错的起点,快速时变系统可以尝试0.95,高噪声环境则建议0.995以上。具体选多少,要看你更在乎收敛速度还是稳态精度。
3. 递推代码的模块拆解:一份可以直接改着用的模板
3.1 初始化参数
下面这份代码是我自己常用的RLS模板,以二阶ARX模型为例。先看初始化部分:
% 参数初始化 n = 2; % 待辨识参数个数,这里以 y(k) = a1*y(k-1) + b1*u(k-1) 为例 theta_hat = zeros(n,1); % 参数估计初值 P = 1e6 * eye(n); % 协方差矩阵初值 lambda = 0.98; % 遗忘因子这里有几个关键点。theta_hat初值在没有先验信息时直接取零向量是合理的,但P矩阵初值一定要取大,比如1e4到1e6倍的单位矩阵。P表示对当前参数估计的置信程度,P越大代表“我越不信当前估计”,因此前几步会有较大的修正量,让参数快速从初始猜测收敛到真值附近。如果P初值设置太小,比如直接取eye(n),算法会以为初始参数已经很准了,新数据的修正作用很小,结果就是收敛极慢,甚至在一段时间内看不到明显变化。
3.2 构造回归向量与递推主循环
核心递推循环大概是这样的:
N = length(y); % 数据长度 theta_rec = zeros(n, N); % 记录每一时刻的参数估计值,方便画图 for k = 2 : N % 构造回归向量,结构与模型形式强相关 phi = [-y(k-1); u(k-1)]; % 预测误差 e = y(k) - phi' * theta_hat; % 增益矩阵 K = P * phi / (lambda + phi' * P * phi); % 参数更新 theta_hat = theta_hat + K * e; % 协方差更新 P = (eye(n) - K * phi') * P / lambda; % 记录结果 theta_rec(:, k) = theta_hat; end这段代码里,最容易被忽视的是回归向量phi的构造。它必须和你的模型结构严格对应。比如模型是y(k) = a1y(k-1) + b1u(k-1),那phi就是[-y(k-1); u(k-1)];如果是二阶ARX模型,phi就要包含y(k-1)、y(k-2)、u(k-1)、u(k-2)四项,同时n也要改成4。很多同学直接把别人模板里的phi拿过来跑,结果参数怎么都收敛不到真值,问题往往就出在模型结构对不上。
另外,在MATLAB里,phi' * P * phi得到的是一个标量,K的维度是n×1,这样参数更新时theta_hat + K*e的维度才能匹配。在Python的numpy里也一样,只不过需要更小心矩阵维度和广播机制。
3.3 数据生成与结果验证
为了验证算法,先得有一组已知真值的仿真数据。我经常用这样的方式生成:
% 生成仿真数据:真实系统取 a1=-1.5, b1=0.3 a1_true = -1.5; b1_true = 0.3; N = 500; u = randn(N,1); % 激励信号用白噪声 y = zeros(N,1); for k = 2:N y(k) = -a1_true * y(k-1) + b1_true * u(k-1) + 0.01 * randn; end激励信号这里不能偷懒用常数或者斜坡。递推最小二乘要求输入信号满足持续激励条件,简单理解就是输入要足够“丰富”,包含足够的频率成分,才能把系统的动态特性激发出来。白噪声和伪随机序列都是常用选择。如果输入一直是直流,回归矩阵会病态,参数无法唯一辨识,这是实验里最经典的翻车原因之一。
跑完递推后,画参数收敛曲线、预测输出对比曲线,再算一下估计值与真值之间的均方根误差,整个实验流程就完整了。
3.4 模板使用时的几个注意事项
- 回归向量要按模型结构改,不能只改参数个数;
- 记录估计结果的theta_rec要预先分配矩阵,避免循环里动态扩容拖慢速度;
- Python实现时,P更新完成后最好做一次对称化处理,P = (P + P.T) / 2,防止数值误差破坏对称正定性;
- 如果数据长度很长,注意lambda过小时协方差矩阵可能数值爆炸,需要定期检查P是否出现NaN或Inf。
4. 实验报告的核心结构:如何把“详细”落到实处
4.1 原理部分不要只堆公式
一份word版实验报告,原理部分最容易犯的毛病是公式连抄三页,却没有一句人话。写递推最小二乘的原理,关键是说清楚“为什么需要递推”。可以从三个角度展开:第一,实时性要求,工业现场不可能等全部数据采集完成再计算;第二,计算量控制,递推形式避免了不断重复的矩阵求逆;第三,时变系统需要遗忘因子,让估计能跟随参数变化。
实验目的也别写“学习递推最小二乘算法”这种空话,而是写得更具可检验性,比如“掌握递推最小二乘的迭代流程,能够通过仿真数据分析遗忘因子对参数估计收敛性能的影响”。这样的目的写着明确,老师看着也舒服。
4.2 仿真结果要能讲出“故事”
仿真结果部分,至少要有三样东西:参数估计收敛曲线、预测输出与真实输出对比图、误差曲线。收敛曲线要能看到参数从初始值出发逐渐逼近水平线;预测输出对比图要能体现“模型输出跟随真实输出”的效果;误差曲线最好同时画出预测误差和参数误差两种。
关键是要对每一张图都有明确的观察结论,不能只贴图不解释。比如这样写:“由图2可以看出,a1估计值从初始值0开始,经过约80个采样点后收敛到-1.5附近,稳态波动范围小于±0.02;遗忘因子λ=0.98时,参数收敛速度高于λ=1的情况,但稳态波动略大。”这种表述说明你真的观察过数据,而不是随手截图。
不同遗忘因子的对比可以用表格呈现:
| 遗忘因子λ | 收敛速度 | 稳态波动 | 时变跟踪能力 |
|---|---|---|---|
| 0.95 | 快 | 较大 | 强 |
| 0.98 | 中等 | 中等 | 中等 |
| 1.00 | 慢 | 较小 | 弱 |
4.3 容易被忽略的加分项
如果想让报告更出彩,可以在仿真中设计一个参数突变场景。比如前250个采样点真实a1=-1.5,第251个采样点开始真实a1突然变成-1.0,然后观察不同遗忘因子下参数估计的跟踪速度。这个实验能非常直观地体现遗忘因子的价值,也是面试或答辩时很容易展开聊的一个点。
另外,把递推结果和离线最小二乘结果做一个对比,说明两者在数据量足够大时应该趋于一致,也是一个很有深度的分析角度。报告中可以写出离线最小二乘得到的参数作为“终值参考”,再用递推结果的最终值去比较,差异能反映递推算法的收敛性能。
4.4 word版报告排版上的实际问题
- 公式用Word公式编辑器输入,不要贴截图,排版更清晰;
- 图用矢量图或者高分辨率PNG,插入后注意缩放比例,确保坐标轴文字能看清;
- 图和表都要有编号和标题,正文里要先有引用,比如“如图1所示”;
- 如果资源包里的报告是别人写好的,务必自己重新生成图表和数据后再用。出于学术诚信考虑,这既是对自己的保护,也是学习的必要环节。
5. 调试递推最小二乘程序时踩过的坑
5.1 参数发散或出现NaN
这是新手最容易遇到的情况。估计值曲线突然冲到离谱数值,或者直接变成NaN/Inf。我复盘过几次,原因基本逃不出这几类:第一,P矩阵在长时间递推后失去正定性,这通常是因为λ太小或者数值舍入误差累积;第二,数据里有缺失值或异常值,喂进了NaN;第三,激励信号太弱,回归矩阵接近病态。
解决办法依次是:把λ往回调到0.98以上;每次更新P之后做一次对称化处理;检查数据是否干净;如果信号确实很弱,考虑加一点随机扰动作为激励。还有一种相对稳妥的变体是带正则化的RLS,在P逆矩阵更新时加一个小正则项,能有效抑制数值发散。
5.2 参数收敛不到真值
参数持续在某个错误值附近波动,最常见的原因是模型结构与实际系统不匹配。你拿一阶模型去辨识一个二阶系统,参数当然无法收敛到“真值”,因为所谓真值在一阶模型里根本不存在。另一个常见原因是信噪比太低,测量噪声淹没了系统动态特征,这时无论什么辨识算法都很难给出准确结果。
我的调试建议是:先用离线最小二乘对整段数据做一次辨识,如果离线结果都明显偏离真值,说明模型结构或数据质量有问题,先解决这个;如果离线结果没问题但递推结果不理想,再检查初始化、遗忘因子、数据顺序等递推特有的环节。这个排查顺序能帮你快速缩小问题范围。
5.3 P初值和遗忘因子是耦合的
P初值取1e6还是1e3,对算法前几步影响会很大。P初值过小会让参数更新幅度变小,收敛变慢;P初值过大会让前几步估计值波动剧烈,但通常能快速拉回真值附近。遗忘因子和P初值会共同影响动态特性,调参时要记住:它们的本质都是在权衡“相信旧信息”和“接受新信息”。
我一般先固定λ=1、P初值取1e6,把基础收敛性跑通;然后再引入λ并逐步减小,观察时变跟踪效果。这样分步调参,能避免两个变量同时变动导致“不知道是谁的功劳”。
5.4 数据预处理的隐形坑
递推最小二乘默认数据是等间隔采样的。如果你的数据采集时间不均匀,直接做递推会在时间语义上出错,结果自然不可靠。这时候需要先做插值重采样。另外,输入输出数据如果有明显的直流分量,回归矩阵中的某些列会高度相关,影响参数辨识的数值稳定性。建议先对输入输出做去均值处理,也就是减去平均值,只保留波动成分,辨识完成后再还原到原始坐标。
6. 从“跑通代码”到“真正掌握”的最后一步
我的经验是,处理这类资源包最有效的路径永远是“先复现、再改造、后内化”。解压之后先给代码每一行写上注释,搞懂每一步对应哪个公式;然后把它改造到你的场景里,无论是换一个被辨识对象,还是把模型阶数提高,都会逼你重新审视回归向量的构造;最后把遗忘因子、P初值当成实验变量去观察影响,记录下你调参时看到的现象。
做到这一步,哪怕把别人的代码和word报告都还回去,算法也已经长在你自己身上了。递推最小二乘是很多在线辨识方法的基础,后面学卡尔曼滤波、自适应控制甚至各类梯度递推算法时,你会发现它们的思维框架惊人相似:都是“预测—修正—更新”的循环。把这个循环吃透,比单纯解压一个rar、跑通一段脚本有意义得多。
本文还有配套的精品资源,点击获取
