最小二乘法拟合的MATLAB和Excel实现

最小二乘法拟合的MATLAB和Excel实现

摘要：生活生产中我们会遇到各种各样的数据处理，然而这些数据并不像理想实验中得到的数据，有的是一元或多元函数的分布，有的是一次或多次函数的分布，这就需要我们首先观察数据的散点图，进而选择合理的选择函数进行拟合，同时分析计算该拟合得到的误差，找出最优的拟合方式。本文从数学上对最小二乘法原理进行了阐述，并通过MATLAB和Excel 完成数据的拟合，在进行数据拟合中使用的一次函数拟合和多项式拟合，并对不同的拟合方式进行了比较，到了不同拟合方式下的拟合函数和拟合误差。同时对MATLAB和Excel数据拟合方式进行了对比。

关键字：最小二乘法 MATLAB Excel 数据拟合

Abstract:we will encounter a variety of data processing in production life .However these data is not the data as we expect in ideal experiment;some distribution is a univariate or multivariate functions, some is one or more times function.So we should observe the scatter data chart,and then choose the reasonable selection function fitting, make an error analysis and find out the best way of fitting. This paper expound the principle of least square mathematically,complete data fitting by MATLAB and Excel,and use a function fitting and polynomial fitting.we also compare the different fitting methods,the fitting function and fitting error by the way of MATLAB and Excel.

Keywords: Least squares MATLAB Excel Data fitting

引言

工程试验中我们常常遇到这样的问题，试验中我们会得到各种各样的数据，不同的数据之间存在着这样那样的关系，如何把得到的试验数据用函数关系式来得到不同组数据之间的关系，并且在经过数据处理后得到的函数能够客观准确的描述数据与数据数据之间的关系。如何选择数据的拟合方式，是线性拟合还是非线性拟合？是一次还是多次拟合？如何使相关系数R2接近1？通常我们用最小二乘法来确定拟合曲线和和该数据的经验公式。本文从数学角度给出最小二乘法的推导过程，从现实生活中给出曲线拟合的实际应用，同时用MATLAB和Excel两种不同的方法对数据进行曲线拟合的实现。

1 最小二乘法原理

当由实验提供了大量数据时,不能要求拟合函数)(x S 在数据点（x,y ） 处的偏差,即i i i y x S -=)(δ(i=1,2,…,m) 严格为零,但为了使近似曲线尽量反映所给数据点的变化趋势 ,需对偏差有所要求.通常要求偏差平方和

∑∑∑===-=-=m

i i i i i m

i m

i y x S y x S 0

22

*

2

])([min ])([δ （1-1）

这里

)()()()(1100x a x a x a x S n n ??? ++= (n < m) （1-2） 这就是一般的最小二乘逼近，用几何语言说，这就称为曲线拟合的最小二乘法。[1]用最小二乘法求拟合曲线的问题，就是在形如（1-2）式的S （x ）中求一函数)(*x S y =使偏差的平方和最小。它转化为求多元函数

∑∑==-=m

i i n

j i j j i n x f x a x a a a I 0

20

10)]()()[(),,(?ω （1-3）

的极小点),,(**1*0n a a a 的问题。有求多元函数极值的必要条件，有

,0)()]()()[(2020

=-=??∑∑==i k m i i n

j i j j i k x x f x a x a I

??ω k=0,1，

···n （1-4） 若记

n

,1,0,)()()(),()

()()(),(0

0 =≡==∑∑==k d x x f x f x x x m

i k i k i i k i k i j m

i i k i ?ω???ω?? （1-5）

则可改写为矩阵形式

d Ga = （1-6） 其中T n T n d d d d a a a a ),,(,),,(1010 ==，

()()()()

()

()()(

)(

)

??

?

???

???

?????=

???

??????????????

?n n n

n n

n

G ,...

,,............,...,,,...,,1

11

1010

1

00

2 最小二乘法应用

例1 已知数据 x 1.2 2.8 4.3 5.4 6.8 7.9 y

2.1 11.5 28.1 41.9 72.3 91.4 i ω

1

1

1

1

1

1

试用最小二乘法拟合这组数据的曲线。

解：根据所给的数据，在坐标纸上标出，见右图，从途中看到各点在一条直线附近，故选择线性函数作拟合曲线，即令x a a x S 10)(+=,这里m=5，n=1，,)(,1)(10x x x ==??故

51

.1595),(,

3.247),(,58.165),(,

4.28),(),(,

6),(5

15

02

50115

01105

000===========∑∑∑∑∑=====i i i i i i i i i i i i i i i f x f f f x x ω?ω?ω??ω????ω??

由法方程得线性方程组

{51

.159558.1654.283

.2474.

286110

=+=+a a a

解得a 0=-23.3498，a 1=13.6408。于是所求拟合曲线为

x x S 6408.133498.23)(1+-=

同理，可以分别求出二次和三次拟合的方程)()(32x S x S 和分别为：

8143

.3*2774.4*6690.2*0899.0)(8328

.0*4740.0*4424.1)(2

3

322+-==-+=x x x x S x x x S

3数据拟合结果的实现和误差的分析

在例1中我们运用数学公式进行了最小二乘法的拟合运算[2]，这种方法在计算过程中太过于复杂，下面介绍两种方法对上例进行数据拟合，首先用MATLAB 进行编程处理。拟合的程序代码见附录。该实验程序得到的拟合图形如下：

图3-1 数据的散点图和不同类型的拟合曲线

在MATLAB 中的得到的拟合函数及均方偏差结果如下：

001

+1.363550e =814

.3*2774.4*669.2*0899.0001

+1.645645e =8328

.0*4740.0*442.1002

+2.734957e 35

.23*64.13232332222211δδδ+-+=-+==-=x x x y x x y x y 从得到的拟合函数和均方偏差可知，用一次拟合得到的一次函数的均方偏差比较的大，也就是用多项式函数拟合得到的均方偏差更小，随着多项式的次数越多，均方偏差趋近于零。

图3-2 用Excel 表格进行数据拟合结果图

用Excel 表格拟合得到的拟合函数及相关系数结果如下:

9978.08143

.3*2774.4*6690.2*0899.09973

.08328.0*4740.0*4424.19549

.0350.23*641.132323322222

11=+-===-+==-=R x x x y R x x y R x y

同理，用Excel 得到的拟合函数和MATLAB 中得到的函数结果相同，在Excel 中进行数据拟合时用相关系数R 2表示误差，(0

可以看出运用MATLAB 或Excel 进行对数据进行曲线拟合比人工计算大大减少了计算量，运用软件处理数据已经是当前时代解决数学问题的重要方法之一。

参考文献

[1] 李庆阳，王能超，易大义.数值分析（第五版）.北京：清华大学出版社，2008 [2] 姜健飞，胡良剑，唐俭.数值分析及其MATLAB 试验.北京：科学出版社，2004

完成人签名：

学院：研究方向：

学号：导师姓名：

完成时间：2014年12月5日

附录

MATLAB编程代码如下：

x=[1.2 2.8 4.3 5.4 6.8 7.9]; y=[2.1 11.5 28.1 41.9 72.3 91.4];

subplot(2,2,1)

stem(x,y);

xlabel('x')

ylabel('y')

title('数据散点图')

p=polyfit(x,y,1); % 拟合出的一次函数的系数

yc=y-polyval(p,x); % 计算误差

s2=sum(yc.^2); % 误差的平方和

disp(sprintf('一次误差的平方和=%d',s2));

x1=linspace(min(x),max(x)); % 绘图用到的点的横坐标

y1=polyval(p,x1); % 拟合曲线的纵坐标

subplot(2,2,2);

plot(x,y,'o',x1,y1); % 绘图，原始数据+拟合曲线

legend('原始数据','一次曲线',2); % 图示

s=char(vpa(poly2sym(p,'x'),4)); % 一次函数式转换为字符串，vpa转换小数，保留4位有效数字

title(['y=' s]);

p=polyfit(x,y,2); % 拟合出的二次函数的系数

yc2=y-polyval(p,x); % 计算误差

s2=sum(yc2.^2); % 误差的平方和

disp(sprintf('二次误差的平方和=%d',s2));

x2=linspace(min(x),max(x)); % 绘图用到的点的横坐标

y2=polyval(p,x2); % 拟合曲线的纵坐标

subplot(2,2,3);

plot(x,y,'o',x2,y2); % 绘图，原始数据+拟合曲线

legend('原始数据','二次曲线',2); % 图示

s=char(vpa(poly2sym(p,'x'),4)); % 二次函数式转换为字符串，vpa转换小数，保留4位有效数字

title(['y=' s]);

p=polyfit(x,y,3); % 拟合出的三次函数的系数

yc3=y-polyval(p,x); % 计算误差

s2=sum(yc3.^2); % 误差的平方和

disp(sprintf('三次误差的平方和=%d',s2));

x3=linspace(min(x),max(x)); % 绘图用到的点的横坐标

y3=polyval(p,x3); % 拟合曲线的纵坐标

subplot(2,2,4);

plot(x,y,'o',x3,y3); % 绘图，原始数据+拟合曲线

legend('原始数据','三次曲线',2); % 图示

s=char(vpa(poly2sym(p,'x'),4)); % 三次次函数式转换为字符串，vpa转换小数，保留4位有效数字

title(['y=' s]);