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

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

例1.
在某个低温过程中，函数依赖于温度的试验数据如下表：

1
2
3
4

0.8
1.5
1.8
2.0
已知经验公式的形式为，根据最小二乘法原理编制MATLAB程序求出，并做相应的理论分析。
解：在两个观测量中，往往总有一个量精度比另一个高得多，为简单起见把精度较高的观测量看作没有误差，并把这个观测量选作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<m时，参数不能确定。
在N>m的情况下，式（0-0-2）成为矛盾方程组，不能直接用解方程的方法求得m个参数值，只能用曲线拟合的方法来处理。设测量中不存在着系统误差，或者说已经修正，则y的观测值yi围绕着期望值 <f（x；c1，c2，……cm）> 摆动，其分布为正态分布，则yi的概率密度为
,
式中是分布的标准误差。为简便起见，下面用C代表（c1，c2，……cm）。考虑各次测量是相互独立的，故观测值（y1，y2，……cN）的似然函数
.
取似然函数L最大来估计参数C，应使
（0-0-3）
取最小值：对于y的分布不限于正态分布来说，式（0-0-3）称为最小二乘法准则。若为正态分布的情况，则最大似然法与最小二乘法是一致的。因权重因子，故式（0-0-3）表明，用最小二乘法来估计参数，要求各测量值yi的偏差的加权平方和为最小。
根据式（0-0-3）的要求，应有

从而得到方程组
（0-0-4）
解方程组（0-0-4），即得m个参数的估计值，从而得到拟合的曲线方程。
然而，对拟合的结果还应给予合理的评价。若yi服从正态分布，可引入拟合的x2量，
（0-0-5）
把参数估计代入上式并比较式（0-0-3），便得到最小的x2值
（0-0-6）
可以证明，服从自由度v＝N-m的x2分布，由此可对拟合结果作x2检验。
由x2分布得知，随机变量的期望值为N-m。如果由式（0-0-6）计算出接近N-m（例如），则认为拟合结果是可接受的；如果，则认为拟合结果与观测值有显著的矛盾。
把化为线性二乘拟合，即： 于是就可以写编程序了。先写一个主程序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

相关的分析略。

Welcome To
Download !!!

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

精品资料