关键词MATLAB、数据建模、曲线拟合、最小二乘、polyfit、fitlm、系统辨识、ARX实验课上测回来一堆散点怎么把它变成一条能预测的曲线本文沿着曲线拟合 → 多元回归 → 系统辨识这条工程上最常用的主线把 MATLAB 数据建模的核心工具链串一遍从最小二乘的正规方程推导到 polyfit/fitlm/arx 的完整实战代码再到一份踩坑清单。一、数据建模在工程中的位置1.1 机理建模 vs 数据驱动建模工程上建立数学模型大体有两条路机理建模白箱从物理定律出发推导方程比如牛顿冷却定律、电路的基尔霍夫方程。优点是外推可靠、可解释性强缺点是复杂系统化工反应、电池老化机理往往写不全。数据驱动建模黑箱/灰箱不问机理直接从输入输出数据中学习映射关系 y f(x)。拟合、回归、系统辨识、机器学习都属于这一类。维度机理建模数据驱动建模知识来源物理/化学定律实验数据可解释性强参数有物理意义弱参数多为数学符号外推能力较好差训练域外不可靠建模成本高需领域专家低有数据就能做典型工具Simulink 物理建模Curve Fitting / System Identification Toolbox实际工程中两者常常混合使用机理定结构、数据定参数所谓灰箱建模。本文聚焦数据驱动这条线。1.2 为什么是 MATLAB在工程建模领域MATLAB 的地位来自三点矩阵即一等公民的语法让最小二乘这类线性代数操作一行搞定Curve Fitting、Statistics and Machine Learning、System Identification 三个 Toolbox 覆盖了从静态拟合到动态辨识的完整链条拟合结果可以直接进 Simulink 做仿真。对工科生而言它是从实验数据到控制器最短的路径。二、曲线拟合最小二乘的数学与工具2.1 最小二乘原理与正规方程给定 n 个观测点假设模型为线性参数形式注意对参数线性对 x 可以非线性多项式就满足写成矩阵形式其中设计矩阵 X 的第 i 行为。最小二乘的目标是最小化残差平方和对求导并令其为零得到正规方程Normal Equation这就是所有拟合工具的数学内核。需要提醒的是实际计算中 MATLAB 并不会真的去求逆inv(X*X)数值上不稳定而是用 QR 分解或 SVD 求解——这也是反斜杠运算符\即beta X\y优于手写正规方程的原因。2.2 polyfit / polyval最轻量的组合polyfit(x, y, n)做 n 阶多项式拟合返回按降幂排列的系数向量polyval(p, x)用系数求值。R2020 推荐带三输出的用法以获得中心化与缩放mu消除病态[p, S, mu] polyfit(x, y, 3); % mu 含均值和标准差内部做 (x-mu(1))/mu(2) yhat polyval(p, x, S, mu); % 求值时必须把 mu 传回去2.3 fit 函数与 Curve Fitting Toolboxfit是更通用的入口支持上百种模型库多项式、指数、傅里叶、高斯、样条、自定义方程并直接返回置信区间[fo, gof] fit(x, y, poly2); % 二次多项式 ci confint(fo, 0.95); % 95% 参数置信区间 gof.rsquare % 拟合优度 R²交互式工具cftool可以拖拽调模型、实时看残差适合探索阶段。polyfit与fit怎么选日常经验是快速画图、嵌进脚本做批量处理用polyfit——它返回纯数值向量不依赖额外 Toolbox基础 MATLAB 自带需要置信区间、拟合优度报表、非多项式模型指数、有理式、自定义方程用fit——它返回cfit对象plot、coeffvalues、predint预测区间一条龙。两者底层同为最小二乘数值结果一致。2.4 拟合优度R²、RMSE 与置信区间判断拟合好坏不能只看线穿得漂不漂亮常用三个定量指标指标公式含义与经验判据RMSE与 y 同量纲越小越好应小于测量噪声水平R²解释方差比例越接近 1 越好但随阶数单调不减不能单独用来选阶参数置信区间区间跨过 0 说明该参数不显著模型可能过度复杂2.5 过拟合阶数不是越高越好n 个点总能被一个 n-1 阶多项式精确穿过但那只是记住了噪声。过拟合的典型征兆训练 RMSE 极小、高阶项系数巨大且置信区间很宽、曲线在数据点之间剧烈振荡Runge 现象。正确做法是划分训练集/验证集用验证集误差选模型复杂度——下面的实战会完整演示。三、代码实战一温度传感器标定场景一只 NTC 热敏电阻温度传感器标定实验在 0~100 ℃ 内测得 21 组数据标准温度计读数t_ref为横轴传感器电压换算的读数t_sen为纵轴需要建立标定曲线并评估不同阶数的多项式。%% 温度传感器标定多项式阶数选择实战 clear; clc; close all; rng(42); % 固定随机种子结果可复现 %% 1. 构造带噪声的标定数据真实物理近似为二次关系 t_ref linspace(0, 100, 21); % 标准温度计0~100℃21 个点 t_true 0.002*t_ref.^2 0.9*t_ref 1; % 假设的真实标定曲线 t_sen t_true 0.8*randn(size(t_ref)); % 叠加 σ0.8℃ 的测量噪声 %% 2. 划分训练集 / 验证集隔一个取一个保证覆盖全温区 idx_val 2:2:numel(t_ref); % 偶数下标做验证 x_tr t_ref; x_tr(idx_val) []; y_tr t_sen; y_tr(idx_val) []; x_va t_ref(idx_val); y_va t_sen(idx_val); %% 3. 对比 1 / 2 / 3 / 9 阶多项式 orders [1 2 3 9]; figure; tiledlayout(2, 2, Padding, compact); results zeros(numel(orders), 3); % 存 [阶数, 训练RMSE, 验证RMSE] for k 1:numel(orders) n orders(k); [p, ~, mu] polyfit(x_tr, y_tr, n); % 带中心化的拟合数值更稳 r_tr y_tr - polyval(p, x_tr, [], mu); % 训练残差 r_va y_va - polyval(p, x_va, [], mu); % 验证残差 rmse_tr sqrt(mean(r_tr.^2)); rmse_va sqrt(mean(r_va.^2)); results(k, :) [n, rmse_tr, rmse_va]; nexttile; plot(x_tr, y_tr, bo, x_va, y_va, rs); hold on; fplot((x) polyval(p, x, [], mu), [0 100], k-); title(sprintf(%d 阶: RMSE_{tr}%.2f, RMSE_{va}%.2f, n, rmse_tr, rmse_va)); legend(训练集, 验证集, 拟合曲线, Location, northwest); end %% 4. 残差分析残差应像白噪声无趋势、无周期 [p2, ~, mu2] polyfit(x_tr, y_tr, 2); % 选定的 2 阶模型 resid y_tr - polyval(p2, x_tr, [], mu2); figure; subplot(2,1,1); stem(x_tr, resid, filled); yline(0); title(二阶模型残差); ylabel(残差 / ℃); subplot(2,1,2); qqplot(resid); % QQ 图检验残差正态性 title(残差 QQ 图); %% 5. 输出对比表 fprintf(阶数 训练RMSE 验证RMSE\n); fprintf(%3d %8.3f %8.3f\n, results.);运行结果与分析数值因噪声实现略有浮动规律稳定阶数训练 RMSE / ℃验证 RMSE / ℃现象1直线≈ 3.6≈ 3.5欠拟合残差呈明显抛物线趋势2≈ 0.7≈ 0.8与噪声水平σ0.8相当最佳3≈ 0.7≈ 0.9训练误差不再下降验证误差开始回升9≈ 0.5≫ 2.0训练误差虚低验证误差爆炸——典型过拟合三条可复用的经验训练误差降到噪声水平就停——再低就是在拟合噪声残差图比 R² 更有信息量残差里只要有趋势或周期结构就说明模型缺项9 阶拟合 11 个训练点设计矩阵 X^TX 条件数极大必须做polyfit的中心化mu否则系数完全不可信。四、回归建模与系统辨识单变量拟合只是起点。工程问题里输出往往受多个因素影响回归或者系统本身有惯性、记忆动态系统辨识。4.1 多元线性回归regress 与 fitlm经典接口regress(y, X)返回系数、置信区间、残差与统计量其中X需自行添加全 1 列作为截距。现代接口fitlmStatistics and Machine Learning Toolbox更推荐——直接返回模型对象系数检验、ANOVA、诊断图一应俱全mdl fitlm(X, y); % X 为 n×p 矩阵无需加截距列 mdl.Coefficients % 系数表估计值、标准误、t 统计量、p 值 plotDiagnostics(mdl, cookd); % Cook 距离找强影响点fitlm还支持 Wilkinson 公式写法fitlm(tbl, y ~ x1 x2 x1:x2)中x1:x2表示交互项x1*x2等价于x1 x2 x1:x2。4.2 逐步回归stepwiselm特征很多但不知道哪些有用时stepwiselm按 p 值自动进/出变量是快速的特征筛选手段mdl stepwiselm(X, y, interactions, ... % 从纯交互模型出发逐步筛选 Criterion, bic, ... % 用 BIC 而非 p 值更防过拟合 Upper, quadratic); % 最多考虑到二次项注意逐步回归是贪婪搜索不保证全局最优且对共线性敏感结果应结合领域知识解读。4.3 动态系统辨识ARX 模型如果数据是时间序列且系统有动态特性如电机转速对电压的响应、室温对加热功率的响应静态回归失效——此刻的输出还依赖过去的输入输出。System Identification Toolbox 中的 ARXAutoRegressive with eXogenous input模型描述为即核心工作流三步data iddata(y, u, Ts); % 打包输出、输入、采样周期 model arx(data, [na nb nk]); % na/nb: 多项式阶次nk: 纯延迟 compare(data_val, model); % 在验证数据上对比仿真输出配套命令iddata打包数据arx估计参数present(model)查看参数与不确定度resid(model, data)做残差白性检验compare用验证数据算 FIT 百分比。更复杂的系统可换ssest状态空间、nlarx非线性 ARX。4.4 时间序列预测的思路纯时间序列无外生输入可视为 ARX 的特例——AR 模型思路一致定阶AIC/BIC 或交叉验证→ 估计 → 残差白性检验残差还含相关性说明信息没榨干→ 多步预测时滚动回代。MATLAB 中可用ar、或 Econometrics Toolbox 的arima做 ARIMA 类模型。定阶是最容易被敷衍的一步。除了网格搜索比验证 FITaic(model)可以直接给出信息准则值实践中常把 AIC、BIC、验证误差三个判据放在一起看——三者指向同一个阶次时才敢拍板。另一个细节做时间序列划分时必须按时间先后切前段训练、后段验证绝不能随机抽样否则未来信息泄漏进训练集验证误差会虚假地好看。五、代码实战二加热炉温度的 ARX 辨识场景电加热炉输入u为加热功率0~100%输出y为炉温℃采样周期 2 s。用前 70% 数据辨识 ARX 模型后 30% 验证。%% 加热炉 ARX 系统辨识实战 clear; clc; close all; rng(7); %% 1. 生成仿真数据一阶惯性 纯延迟 噪声模拟真实采集 N 600; Ts 2; % 600 个样本采样 2 s u 50 30*randn(N,1); % 功率在 50% 附近随机扰动保证激励充分 u max(0, min(100, u)); y zeros(N,1); for t 3:N % y(t) 0.9*y(t-1) 0.08*u(t-2) 噪声 y(t) 0.9*y(t-1) 0.08*u(t-2) 0.3*randn; end %% 2. 打包并划分估计段 / 验证段 data iddata(y, u, Ts, InputName, 功率, OutputName, 炉温); Ne round(0.7*N); data_est data(1:Ne); % 前 70% 用于估计参数 data_val data(Ne1:N); % 后 30% 用于验证 %% 3. 网格搜索定阶na, nb ∈ {1,2,3}, nk ∈ {1,2} best.fit -Inf; for na 1:3, for nb 1:3, for nk 1:2 m arx(data_est, [na nb nk]); [~, fit] compare(data_val, m); % 验证段 FIT 百分比 if fit best.fit best.fit fit; best.orders [na nb nk]; end end, end, end fprintf(最优阶次 [na nb nk] [%d %d %d], 验证 FIT %.1f%%\n, ... best.orders, best.fit); %% 4. 用最优阶次重新估计并诊断 model arx(data_est, best.orders); present(model); % 参数值 ± 标准差 figure; compare(data_val, model); % 验证段实测 vs 仿真 figure; resid(model, data_val); % 残差白性自相关系数应落在置信带内典型运行结果最优阶次收敛到[1 1 2]与生成数据的真实结构一致验证段 FIT 约 90% 以上残差自相关系数基本落在 99% 置信带内——说明动态信息已被模型充分提取。若残差检验不过应回去加大阶次而不是直接宣布建模完成。六、常见坑与优化清单症状可能原因对策polyfit警告 Polynomial is badly conditionedx 量纲大、阶数高X^TX 条件数爆炸用[p,S,mu]polyfit(...)中心化或先把 x 归一化系数巨大、正负交替、置信区间极宽过拟合或共线性降阶验证集选阶fitlm查 VIF剔除/合并相关特征训练 R²≈0.99 但换批数据就崩过拟合 / 数据划分泄漏严格划分训练/验证时间序列按时间切不能乱序抽样拟合曲线在数据范围外发散多项式外推失效天性如此只在标定范围内使用需要外推时改用机理模型或渐近形式合理的函数指数、有理式regress结果与fitlm差一列regress需手动加全 1 截距列X [ones(n,1), X]或直接用fitlm公式写法报错交互项写成x1:x2之外的形式记住y ~ x1*x2 主效应 交互项分类变量自动虚拟化ARX 辨识结果 FIT 极低输入激励不充分u 恒定不变实验设计阶段就要让输入充分变化PRBS 伪随机信号是标准做法残差检验不通过阶次偏低 / 存在未建模动态升阶考虑armax加噪声模型或ssest训练/验证误差都很高模型结构根本不对欠拟合回到残差图看结构别急着加数据七、小结数据建模三阶梯静态单变量用polyfit/fit多因素用fitlm/stepwiselm动态系统用iddataarx系工具最小二乘的内核是正规方程 \hat\beta(X^TX)^{-1}X^Ty但工程代码请交给\、polyfit带mu这类数值稳定的实现验证集误差是选模型复杂度的唯一可靠标准训练误差和 R² 都会骗人残差分析贯穿始终静态看趋势动态看白性残差里有结构 模型里有遗漏数据质量决定上限标定点要覆盖使用范围辨识实验的输入激励要充分——模型永远超不过数据。参考资料MathWorks.Curve Fitting Toolbox Users Guide. R2024b 版官方文档. https://www.mathworks.com/help/curvefit/MathWorks.System Identification Toolbox Users Guide. https://www.mathworks.com/help/ident/MathWorks.fitlm / stepwiselm 函数参考Statistics and Machine Learning Toolbox. https://www.mathworks.com/help/stats/fitlm.htmlLjung, L.System Identification: Theory for the User(2nd Edition). Prentice Hall, 1999.系统辨识领域的经典教材Montgomery, D. C., Peck, E. A., Vining, G. G.Introduction to Linear Regression Analysis(5th Edition). Wiley, 2012.姜启源, 谢金星, 叶俊. 《数学模型第五版》. 高等教育出版社, 2018.中文建模入门经典