数值分析上机实验报告。
选题:曲线拟合的最小二乘法。
指导老师:
专业: 学号:
姓名: 课题八曲线拟合的最小二乘法。
一、问题提出。
从随机的数据中找出其规律性,给出其近似表达式的问题,在生产实践和科学实验中大量存在,通常利用数据的最小二乘法求得拟合曲线。
在某冶炼过程中,根据统计数据的含碳量与时间关系,试求含碳量y与时间t的拟合曲线。
二、要求。1、用最小二乘法进行曲线拟合;
2、近似解析表达式为;
3、打印出拟合函数,并打印出与的误差,;
4、另外选取一个近似表达式,尝试拟合效果的比较;
5、*绘制出曲线拟合图*。
三、目的和意义。
1、掌握曲线拟合的最小二乘法;
2、最小二乘法亦可用于解超定线代数方程组;
3、探索拟合函数的选择与拟合精度间的关系。
四、计算公式。
对于给定的测量数据(xi,fi)(i=1,2,…,n),设函数分布为。
特别的,取为多项式。
(j=0, 1,…,m)
则根据最小二乘法原理,可以构造泛函。
令。(k=0, 1,…,m)
则可以得到法方程。
求该解方程组,则可以得到解,因此可得到数据的最小二乘解。
曲线拟合:实际工作中,变量间未必都有线性关系,如服药后血药浓度与时间的关系;疾病疗效与疗程长短的关系;毒物剂量与致死率的关系等常呈曲线关系。曲线拟合是指选择适当的曲线类型来拟合观测数据,并用拟合的曲线方程分析两变量间的关系。
五、结构程序设计。
在程序结构方面主要是按照顺序结构进行设计,在进行曲线的拟合时,为了进行比较,在程序设计中,直接调用了最小二乘法的拟合函数polyfit,并且依次调用了plot、figure、hold on函数进行图象的绘制,最后调用了一个绝对值函数abs用于计算拟合函数与原有数据的误差,进行拟合效果的比较。
5.1用一元三次多项式进行拟合。
计算解析表达式系数: a1, a2, a3
t=[0 5 10 15 20 25 30 35 40 45 50 55];
y=[0 1.27 2.16 2.86 3.44 3.87 4.15 4.37 4.51 4.58 4.02 4.64];
> n=length(xi);
f=0.34364.*10.^(4)*x.^3-5.2156.*10.^(3)*x.^2+0.26340.*x+0.017839;
x=0:0.01:55;
f=0.34364.*10.^(4)*x.^3-5.2156.*10.^(3)*x.^2+0.26340.*x+0.017839;
fy=abs(f-y);
fy2=fy.^2;ew=max(fy),e1=sum(fy)/n,e2=sqrt((sum(fy2))/n)
plot(xi,y,'t*')hold on, plot(t,f,'b-')hold off
所得函数为。
运行后屏幕显示数据与拟合函数f的最大误差ew,平均误差e1和均方根误差e2及其数据点和拟合曲线y=f(x)的图形如图5.1.
ew =0.4243
e1 =0.0911
e2 =0.1467
图5.1一元三次多项式拟合曲线误差图。
5.2用一元四次多项式进行拟合:
计算多项式系数:a1, a2, a3, a4
xi=[0 5 10 15 20 25 30 35 40 45 50 55];
y=[0 1.27 2.16 2.86 3.44 3.87 4.15 4.37 4.51 4.58 4.02 4.64];
n=length(xi);
x=0:0.01:55;
f=0.6026.*10.
^(6)*x.^4-0.31918.
*10.^(4)*x.^3-0.
0029323.*x.^2+0.
23807.*x+0.060449;
x=0:0.01:55;
f=0.6026.*10.
^(6)*x.^4-0.31918.
*10.^(4)*x.^3-0.
0029323.*x.^2+0.
23807.*x+0.060449;
fy=abs(f-y);fy2=fy.^2;ew=max(fy),e1=sum(fy)/n,e2=sqrt((sum(fy2))/n)
plot(xi,y,'r*')hold on, plot(x,f,'b-')hold off
所得函数为。
运行后屏幕显示数据与拟合函数f的最大误差ew,平均误差e1和均方根误差e2及其数据点和拟合曲线y=f(x)的图形如图5.2。
ew = 0.3897
e1 = 0.1034、
e2 =0.1429
图5.2一元四次多项式拟合曲线误差图。
5.3用一元二次多项式进行拟合:
计算多项式系数:a1, a2, a3
输入程序:> syms a1 a2 a3
x=[0 5 10 15 20 25 30 35 40 45 50 55];
fi=a1.*x.^2+ a2.*x+ a3
运行后屏幕显示关于a1,a2和a3的线性方程组:
fi=[ a3, 25*a1 + 5*a2 + a3, 100*a1 + 10*a2 + a3, 225*a1 + 15*a2 + a3, 400*a1 + 20*a2 + a3, 625*a1 + 25*a2 + a3, 900*a1 + 30*a2 + a3, 1225*a1 + 35*a2 + a3, 1600*a1 + 40*a2 + a3, 2025*a1 + 45*a2 + a3, 2500*a1 + 50*a2 + a3, 3025*a1 + 55*a2 + a3]
编写构造误差平方和的matlab程序:
y=[0 1.27 2.16 2.86 3.44 3.87 4.15 4.37 4.51 4.58 4.02 4.64];
fi =[a3, 25*a1 + 5*a2 + a3, 100*a1 + 10*a2 + a3, 225*a1 + 15*a2 + a3, 400*a1 + 20*a2 + a3, 625*a1 + 25*a2 + a3, 900*a1 + 30*a2 + a3, 1225*a1 + 35*a2 + a3, 1600*a1 + 40*a2 + a3, 2025*a1 + 45*a2 + a3, 2500*a1 + 50*a2 + a3, 3025*a1 + 55*a2 + a3];
fy=fi-y;fy2=fy.^2;j=sum(fy.^2)
运行后屏幕显示误差平方和如下:
j =(100*a1 + 10*a2 + a3 - 54/25)^2 + 25*a1 + 5*a2 + a3 - 127/100)^2 + 225*a1 + 15*a2 + a3 - 143/50)^2 + 400*a1 + 20*a2 + a3 - 86/25)^2 + 900*a1 + 30*a2 + a3 - 83/20)^2 + 625*a1 + 25*a2 + a3 - 387/100)^2 + 1225*a1 + 35*a2 + a3 - 437/100)^2 + 1600*a1 + 40*a2 + a3 - 451/100)^2 + 2025*a1 + 45*a2 + a3 - 229/50)^2 + 2500*a1 + 50*a2 + a3 - 201/50)^2 + 3025*a1 + 55*a2 + a3 - 116/25)^2 + a3^2
为求使达到最小,只需利用极值的必要条件,得到关于的线性方程组,这可以由下面的matlab程序完成,即输入程序 :
> syms a1 a2 a3
j =(100*a1 + 10*a2 + a3 - 54/25)^2 + 25*a1 + 5*a2 + a3 - 127/100)^2 + 225*a1 + 15*a2 + a3 - 143/50)^2 + 400*a1 + 20*a2 + a3 - 86/25)^2 + 900*a1 + 30*a2 + a3 - 83/20)^2 + 625*a1 + 25*a2 + a3 - 387/100)^2 + 1225*a1 + 35*a2 + a3 - 437/100)^2 + 1600*a1 + 40*a2 + a3 - 451/100)^2 + 2025*a1 + 45*a2 + a3 - 229/50)^2 + 2500*a1 + 50*a2 + a3 - 201/50)^2 + 3025*a1 + 55*a2 + a3 - 116/25)^2 + a3^2;
ja1=diff(j,a1);ja2=diff(j,a2);ja3=diff(j,a3);
ja11=******(ja1),ja21=******(ja2),ja31=******(ja3),运行后屏幕显示j分别对a1, a2 ,a3的偏导数如下:
ja11 =49967500*a1 + 1089000*a2 + 25300*a3 - 217403/2
ja21 = 1089000*a1 + 25300*a2 + 660*a3 - 27131/10
ja31 = 25300*a1 + 660*a2 + 24*a3 - 3987/50
解线性方程组ja11 =0,ja21 =0,ja31 =0输入下列程序:
>a=[49967500,1089000,25300;1089000, 25300,660;25300,660,24];
b=[217403/2,27131/10,3987/50];
c=b/a, f=poly2sym(c)
运行后屏幕显示拟合函数f及其系数c如下:
c =-0.0024 0.2037 0.2305
f = 7338734818964133*x)/36028797018963968 - 5489104202452799*x^2)/2305843009213693952 + 8303449950332545/36028797018963968
故所求的拟合曲线为:
编写下面的matlab程序估计其误差,并作出拟合曲线和数据的图形。输入程序:
> xi=[0 5 10 15 20 25 30 35 40 45 50 55];y=[0 1.27 2.16 2.
86 3.44 3.87 4.
15 4.37 4.51 4.
58 4.02 4.64];
n=length(xi);
f=-0.0023805.*x.^2+0.20369.*x+0.23047;
x=0:0.01:55;f=-0.0023805.*x.^2+0.20369.*x+0.23047;
数值分析上机作业
一。上机作业任务一 用五点差分格式求解poisson方程边值问题,p 任选一题 二。上机作业任务二 用simpson求积法计算定积分。下面两种方法任选一。1 变步长复化simpson求积法2 自适应simpson求积法。三。上机作业任务三 用matlab语言编写连续函数最佳平方逼近的算法程序 函数式...
数值分析上机作业
今天的上机作业。1.lagrange 插值。给出的数值表。用lagrange 插值计算的近似值。2 newton插值。用newton插值计算x 0.41的近似值。3.插值法的全部内容。把chap 2试验。doc的全部内容作一边,粘在这个文件里 包括图形 答 插值。function f lagrang...
数值分析上机作业
第二次上机作业。一。任务 用matlab语言编写连续函数最佳平方逼近的算法程序 函数式m文件 并用此程序进行数值试验,写出计算实习报告。二。程序功能要求 在后面的附一的基础上进行修改,使其更加完善。要求算法程序可以适应不同的具体函数,具有一定的通用性。所编程序具有以下功能 1.用lengendre多...