兵工自动化 软件技术 O. I. Automation 2005 年第 24 卷第 3 期 Software Technique 2005, Vol. 24, No. 3 文章编号：1006-1576（2005）03-0107-02 最小二乘曲线拟合及 Matlab 实现 陈光，任志良，孙海柱 （海军工程大学 兵器工程系，湖北 武汉 430033） 摘要：采用最小二乘曲线拟合，可寻求有限测量数据及其伴随误差的变化规律。曲线拟合先确定拟合模型，再确 定函数的所属类。多项式拟合先将其化为双曲线、S 型曲线、倒指数曲线、对数曲线等拟合曲线，再求解拟合多项式 系数。并用 Matlab 编制程序，对测量数据进行拟合与仿真。 关键词：曲线拟合；最小二乘；Matlab；多项式拟合 中图分类号：TP391.9; TP301.6 文献标识码：A Curve Fitting in Least-Square Method and Its Realization with Matlab (Dept. of Weapon Engineering, Naval University of Engineering, Wuhan 430033, China) CHEN Guang, REN Zhi-liang, SUN Hai-zhu Abstract: Least-square curve fitting method is used to seek the changing rule of the definite measured data and its concomitant error. Fitting model is decided firstly, and the belonged class of fitting function is decided. For polynomial fitting method, polynomial was changed into fitting curves of hyperbola, S-figure curve, converse exponential curve and logarithm curve, etc. And polynomial coefficient was solved, the polynomial modulus was programmed by Matlab. Measured data was fitted and simulated with this program. Keywords: Curve fitting; Least-square; Matlab; Polynomial fitting 1 引言 2.1 多项式拟合 由于磁偶极子在导电介质中的传播衰减大，并 具有强非线性，测量数据与被测物理量的真值不完 全一致，需对数据进行必要的数学加工和处理，寻 求被测物理量的变化规律。解决此类问题的常用方 法是对测量数据进行曲线拟合[1]。故基于最小二乘 曲线拟合，用 Matlab 实现对导电介质中磁偶极子辐 射场的传播特性的曲线拟合与仿真。 2 最小二乘曲线拟合 给定一组测量数据{(xi,yi),i＝0,1,2,…,m}，基于 最小二乘原理，求得变量 x 和 y 之间的函数关系 f(x,A)，使它最佳地逼近或拟合已知数据[2]。f(x,A)) 称为拟合模型，A＝(a0,a1,…,an)是一些待定参数。 做法是选择参数 A 使得拟合模型与实际观测值在各 点的残差 ek＝yk-f(xk,A)的加权平方和最小，即求 f*(x)使 ， )y)x(*f)(x( i )y)x(f)(x( i min ω ω − = − 2 2 i i i i m ∑ 0i = m ∑ 0i = ω(xi)≥0，称为权，它反映数据(xi,yi)在实验中所占 数据的比重。应用此法拟合的曲线称为最小二乘拟 合曲线。 用 最 小 二 乘 法 求 拟 合 曲 线 首 先 要 确 定 拟 合 模 型 f(x)，一般来说，根据各门科的知识可以大致确 定函数的所属类，若不具备这些知识，则通常从问 题的运动规律及给定数据的散点图来确定拟合曲线 的形式。 若拟合模型 f(x,A)＝a0+a1x+…+anxn，则称其为 多项式拟合。由最小二乘法确定系数 a0,a1,…,an， 假 设 各 数 据 点 的 权 为 1 ， 令 φ(a0,a1, … ,an) ＝ ∑ ＝最小，则有： = ∑ xa xa )y a( L − + + + e n i 2 i m m 2 n 1 0 i i 0i = 0i = ϕ∂ a ∂ j = 2 m ∑ 0i = a(x j i 0 + xa 1 i + L + xa n n i − )y i = 0 j＝0,1,2,…,n 即： j1 + + xa i1 j xa( i0 m ∑ 0i = 得方程组： + L + xa n jn + i ) = m ∑ 0i = j yx i i i 2 i L L m ∑ 0i = m x ∑ 0i = x n i 1n + i ⎡ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ i x + 1m m ∑ 0i = L m x ∑ 0i = n i x x m ∑ 0i = m ∑ 0i = L x m ∑ 0i = ⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ × = 0 1 a a ⎡ ⎢ ⎢ ⎢ M ⎢ a ⎢ ⎣ ⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ m ∑ 0i = m ∑ 0i = ⎡ ⎢ ⎢ ⎢ ⎢ ⎢ M ⎢ ⎢ ⎢ ⎣ L L L x m ∑ 0i = n2 i 1n + i i n x x L m ∑ 0i = 此方程称为多项式拟合的法方程。令 ⎡ ⎡ + 1m ⎢ ⎢ ⎢ ⎢ m x ⎢ ∑ ⎢ ⎢ 0i = ⎢ ⎢ ⎢ L M ⎢ m ⎢ x ∑ m ⎢ ∑ ⎣ 0i = ⎢ ⎣ 0i = 则得：XA＝Y，从而 A＝X-1Y。 m ∑ 0i = m x ∑ 0i = LLL x x m ∑ 0i = m x ∑ 0i = m ∑ 0i = m ∑ 0i = ⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ ⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ m ∑ 0i = m ∑ 0i = n yx i yx i 1n + i 1n + i n2 i L L Y = y n i 2 i n i i i i i A X = y i yx i i yx n i i ⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ = a a ⎡ ⎢ ⎢ ⎢ M ⎢ a ⎢ ⎣ 0 1 n ⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ 收稿日期：2004-12-09；修回日期：2004-12-31 作者简介：陈光（1980-），男，吉林人，在读硕士，1999 年毕业于海军工程大学，从事水下目标探测与制导研究。 ·107· 

兵工自动化 软件技术 O. I. Automation 2005 年第 24 卷第 3 期 Software Technique 2005, Vol. 24, No. 3 由此矩阵方程解出系数向量 A，即得拟合多项 式 f(x,A)＝a0+a1x+…+anxn。值得注意的是，当 n 较 大时（n≥7），法方程的系数矩阵是高度病态的，也 就是说解的舍入误差很大，以至毫无意义。因此， 用多项式作曲线拟合时，n 不宜取得很大。 频率和电流强度激励的条件下，在不同距离处分别 测得一组电磁场强度数据，如表 1 所示。 表 1 测量数据 x y 10 0.92 0.85 0.78 0.70 0.64 0.60 0.54 0.49 0.45 0.41 5 6 7 1 8 9 2 3 4 当 n=1 时得一种简单的拟合模型 y＝a0+a1x， 称为线性拟合。这种模型计算简单、应用广泛。 2.2 可化为线性拟合的其他曲线拟合 若拟合模型为非多项式形式，则可通过变数 变 换 将 其 化 为 线 性 问 题 。 利 用 最 小 二 乘 线 性 拟 合 确 定 系 数 ， 再 利 用 逆 变 换 给 出 原 问 题 的 曲 线 拟 合 函数。如双曲线 1/y＝a+(b/x)(a>0)（图 1）。可令 y΄＝1/y，x΄＝1/x 则得 y΄＝a+bx΄。 从而由线性拟 合求出 a,b。可采用此法拟合的曲线还有 S 型曲线 y＝1/(a+be-x)，倒指数曲线 y＝aeb/x(a>0)（图 2）， 及对数曲线 y＝a+blnx 等 。 这 些 曲 线 模 型 都 可 通 过线性拟合来求得。 y y 1/a 0 1/a 0 -b/a(b<0) x 图 1 双曲线 图 2 倒指数曲线 3 Matlab 程序实现 x Matlab 是一个高级的数值分析、处理与计算软 件[3]，求解矩阵方程方便。采用 Matlab 编写的求解 多项式拟合系数的程序清单如下： function A=nihe(x,y,n) %定义多项式拟合系数求解函数, % x、y 为输入数据量，n 为拟合次数 m=length(x); %测量数据长度 X1=zeros(1,2*n); %生成 X 矩阵 for i=1:2*n X1(i)=sum(x.^i); end X2=[m,X1(1:n)]; X3=zeros(n,n+1); for j=1:n X3(j,:)=X1(j:j+n); end X=[X2;X3]; Y=zeros(1,n); %生成 Y 向量 for k=1:n Y(k)=sum(x.^k.*y); end Y=[sum(y),Y]; Y=Y'; A=X\Y; %求得拟合系数向量 A 上面的程序是采用 Matlab 的函数文件编写的， 它使用简单、灵活。求解拟合系数向量时只须在调 用的 nihe 文件中输入测量数据和拟合次数即可。 4 算例 将给定的磁偶极子天线置于导电介质中，一定 ·108· 表中，x 为传播距离(m)，y 为归一化的场强衰 减系数。根据电磁场理论知，y 和 x 的关系服从 y ＝aebx 的变化规律[4]，则应用上述的 Matlab 程序对 数据进行拟合如下： x=1:10 y=[0.92 0.85 0.78 0.70 0.64 0.60 0.54 0.49 0.45 0.41] y=log(y); A=nihe(x,y,1) 得：a＝e0.0137＝1.0，b＝-0.0903。拟合关系曲 线 y＝1.0e-0.0903x 测量数据散点与拟合曲线如图 3 所示。由图可 见，所拟合的曲线与测量数据点十分接近，拟合的 效果较为理想，这也进一步说明了该拟合方法是准 确与可行的。 散点图 拟合曲线图 1.1 1 0.9 0.8 0.7 0.6 0.5 衰 减 系 数 y 0.4 0 2 4 6 8 10 传播距离 x 图 3 测量数据散点图与拟合曲线图 5 结束语 通过 Matlab 实现对磁偶极子辐射场测量数据 的曲线拟合，可在有限的测量数据条件下精确描述 导电介质中电磁波的传播特性，为实验研究与工程 应用提供依据。基于最小二乘曲线拟合及 Matlab 实现方法简明、适用，可应用于类似的测量数据处 理和实验研究。 参考文献： [1] 周陪 森, 刘震 涛, 吴淑荣. 自动 检测与 仪表[M]. 北京: 清华大学出版社, 1987. [2] 何汉林, 魏汝祥, 李卫军. 数值分析[M]. 武汉: 湖北科 学技术出版社, 1999. [3] 何仁斌. MATLAB6 工程计算及应用[M]. 重庆: 重庆大 学出版社, 2001. [4] 牛中奇, 朱满座, 卢志远, 等. 电磁场理论基础[M]. 北 京: 电子工业出版社, 2001. [5] 易芳. 采用 MATLAB 的线性回归分析[J]. 兵工自动化, 2004, 23 (2): 68-69.