利用最小二乘法进行数据拟合

______________________________________________________________________________________________________________
精品资料
利用最小二乘法进行数据拟合:

例1.
在某个低温过程中，函数y依赖于温度()C的试验数据如下表：
i

1 2 3 4

i
y
0.8 1.5 1.8 2.0

已知经验公式的形式为2yab，根据最小二乘法原理编制MATLAB程序求
出,ab，并做相应的理论分析。
解：
在两个观测量中，往往总有一个量精度比另一个高得多，为简单起见把精

度较高的观测量看作没有误差，并把这个观测量选作x，而把所有的误差只认为是y的误差。
设x和y的函数关系由理论公式
y＝f（x；c1，c2，……cm）
（0-0-1）

给出，其中c1，c2，……cm是m个要通过实验确定的参数。对于每组观测数据（
xi，
yi）i＝1，2，……，N。都对应于xy
平面上一个点。若不存在测量误差，则这些数据点都

准确落在理论曲线上。只要选取m组测量值代入式（0-0-1），便得到方程组
yi＝f（x；c1，c2，……cm）
（0-0-2）

式中i＝1，2，……，m.求m个方程的联立解即得m个参数的数值。显然N数不能确定。
在N>m的情况下，式（0-0-2）成为矛盾方程组，不能直接用解方程的方法求得m个
参数值，只能用曲线拟合的方法来处理。设测量中不存在着系统误差，或者说已经修正，则
y的观测值yi围绕着期望值 摆动，其分布为正态分布，则y
i
的

概率密度为
______________________________________________________________________________________________________________
精品资料










22212,......,,;exp21imiiiicccxfy

yp



,

式中i是分布的标准误差。为简便起见，下面用C代表（c1，c2，……cm）。考虑各次
测量是相互独立的，故观测值（y1，y2，……cN）的似然函数











NiiiNNCxfyL12221;

2
1exp...21




.

取似然函数L最大来估计参数C，应使


min;1122NiiiiCxfy


（0-0-3）

取最小值：对于y的分布不限于正态分布来说，式（0-0-3）称为最小二乘法准则。若

为正态分布的情况，则最大似然法与最小二乘法是一致的。因权重因子2/1ii，故式
（0-0-3）表明，用最小二乘法来估计参数，要求各测量值yi的偏差的加权平方和为最小。
根据式（0-0-3）的要求，应有

mkCxfycccNiiiik,...,2,10;1ˆ122











从而得到方程组

mkCCxfCxfyccNikiii,...,2,10;;1ˆ12






（0-0-4）

解方程组（0-0-4），即得m个参数的估计值mcccˆ,...,ˆ,ˆ21，从而得到拟合的曲线方程

mcccxfˆ,...,ˆ,ˆ;21
。

然而，对拟合的结果还应给予合理的评价。若yi服从正态分布，可引入拟合的x2量，




NiiiiCxfyx1
2

2
2

;

1


（0-0-5）

把参数估计mccccˆ,...,ˆ,ˆˆ21代入上式并比较式（0-0-3），便得到最小的x2值




Niiiicxfyx1
2

2
2

min
ˆ

;

1


（0-0-6）

可以证明，2minx服从自由度v＝N-m的x2分布，由此可对拟合结果作x2检验。
由x2分布得知，随机变量2minx的期望值为N-m。如果由式（0-0-6）计算出2minx接近N-m（例
______________________________________________________________________________________________________________
精品资料
如mNx2min），则认为拟合结果是可接受的；如果22minmNx，则认为拟合
结果与观测值有显著的矛盾。
把2yab化为线性二乘拟合，即：yab 于是就可以写编程序了。先
写一个主程序zxecf.m，然后写一个被调用程序phi_k.m
Matlab程序如下：
function S = zxecf(x,y,n,w)
global i;global j;
if nargin < 4
w = 1;
end
if nargin < 3
n = 1;
end
phi2 = zeros(n+1);
for i = 0:n
for j = 0:n
phi2(i+1,j+1) = sum((w.*phi_k(x,i)).*phi_k(x,j));
end
end
phif = zeros(n+1,1);
for i = 0:n
phif(i+1) = sum((w.*phi_k(x,i)).*y);
______________________________________________________________________________________________________________
精品资料
end
S = phi2\phif;

function y = phi_k(x,k)
if k == 0
y=ones(size(x));
else
y=x.^k;
end

现在就在Matlab命令窗口输入如下命令：
x = [1 2 3 4];
y = [0.8 1.5 1.8 2.0];
plot(x,y,'*');
S = zxecf(x,y./x);
hold on
a = S(1);
b=S(2);
x1 = 1:0.001:4;
y1 = a*x1+b*x1.^2;
plot(x1,y1,'r-');
title('最小二乘法拟合图像');
______________________________________________________________________________________________________________
精品资料
xlabel('x轴');ylabel('y轴');
gtext('y = ax +bx ^2的图像');
a,b

运行后得到结果如下：
a =

0.9250

b =
-0.1050
______________________________________________________________________________________________________________

精品资料
11.522.533.54
0.8
1
1.2
1.4
1.6
1.8
2
2.2
最小二乘法拟合图像

x轴
y
轴
y = ax +bx 2的图像

相关的分析略。
______________________________________________________________________________________________________________

精品资料
Welcome To
Download !!!

欢迎您的下载，资料仅供参考！