最小二乘法解决人口预测问题假设人口满足指数增长类型.x(t)=xo×e^rt,r为人口增长率,x.为初始人口,x(t)为t时刻的人口总数.年:1790 1800 1810 1820 1830 1840 1850 1860 1870 1880 1890 1900人口:3.9 5.3 7.2 9.

来源:学生作业帮助网 编辑:六六作业网 时间:2024/12/19 02:22:20
最小二乘法解决人口预测问题假设人口满足指数增长类型.x(t)=xo×e^rt,r为人口增长率,x.为初始人口,x(t)为t时刻的人口总数.年:17901800181018201830184018501

最小二乘法解决人口预测问题假设人口满足指数增长类型.x(t)=xo×e^rt,r为人口增长率,x.为初始人口,x(t)为t时刻的人口总数.年:1790 1800 1810 1820 1830 1840 1850 1860 1870 1880 1890 1900人口:3.9 5.3 7.2 9.
最小二乘法解决人口预测问题
假设人口满足指数增长类型.x(t)=xo×e^rt,r为人口增长率,x.为初始人口,x(t)为t时刻的人口总数.
年:1790 1800 1810 1820 1830 1840 1850 1860 1870 1880 1890 1900
人口:3.9 5.3 7.2 9.6 12.9 17.1 23.2 31.4 38.6 50.2 62.9 76.0
要求用线性最小二乘法求x.和r.
要求:
1、解释最小二乘法原理(越详细越好)
2、写出用最小二乘法求解的详细步骤(不要只摆公式,要给出如何利用以上数据的过程)
3、给出MATLAB求解的代码

最小二乘法解决人口预测问题假设人口满足指数增长类型.x(t)=xo×e^rt,r为人口增长率,x.为初始人口,x(t)为t时刻的人口总数.年:1790 1800 1810 1820 1830 1840 1850 1860 1870 1880 1890 1900人口:3.9 5.3 7.2 9.
最小二乘法与数据拟合
一、问题
某公交公司1路车过去20个季度内的客流量(单位:百万)如下表:
季度序号 1 2 3 4 5 6 7 8 9 10
客流量 1.85 2.18 1.6 2.31 1.93 2.35 1.638 2.51 1.92 2.49
季度序号 11 12 13 14 15 16 17 18 19 20
客流量 2.01 2.68 2.05 2.77 2.23 2.95 2.46 2.87 2.35 3.07

试确定客流量与季度序号之间的函数关系,并预测未来八个季度的客流量.
二、实验目的
掌握最小二乘法的原理并会用于解决实际问题;学会用mathematica中的算符进行曲线拟合.
三、预备知识
1 最小二乘法
在许多实际问题中,往往需要根据实验测得两个变量x与y的若干组实验数据(x1,y1),…(xn,yn)来建立这两个变量的函数关系的近似式,这样得到的函数近似式称为经验公式.
通过对实验数据的处理,能够判断x、y大体上满足某种类型的函数关系y=f(x,a1,a2,…,as),但是其中s个参数a1,a2,…,as的值需要通过n组实验数据来确定,通常可以这样来确定参数:选择参数a1,a2,…,as,使得f(x,a1,a2,…,as)在x1,x2 …xn处的函数值与实验数据 y1,y2 …yn 的偏差的平方和为最小,就是使
d= (1)
为最小,这种方法称为最小二乘法.当f(x,a1,a2,…,as)是s个参数的线性函数时,利用求极值与解线性方程组的方法可以解决.
例如,若x、y大体上满足线性关系即f(x,a,b)=ax+b ,则
d(a,b)= (2)
由多元极值的求法有
(3)
解上述关于a、b的二元一次方程组得

a= , b=
从而求得经验公式y=ax+b .
d= 的大小是衡量经验公式精度的一种尺度.
线性函数是最简单最常用的经验公式,有一些实际问题,它们的经验公式可能不是线性函数,我们可以把它化为线性函数来讨论,例如y=kemx,两边取对数得lny=mx+lnk ,令z=lny,b=lnk,即可化为z=mx+b.
2 mathematica中数据拟合算符的用法
在数据处理中常常设法用一个函数按照某种法则去描述一组数据,这就是数据拟合.上面介绍的最小二乘法就是一种最常用的数据拟合方法.mathematica中最基本的数据拟合算符是 fit[ ] ,语法为
fit[数据,拟合函数的基函数列表,变量]
线性函数拟合的基函数为1,x ,n阶多项式拟合的基函数是1,x,x2,…xn .
例 一册书的成本费y与印刷的册数x有关,统计数据如下:
xi(千册) 1 2 3 4 5 6 7 8 9 10
yi(元) 10.15 5.52 4.08 2.85 2.11 1.62 1.41 1.30 1.21 1.15
试用y=a+ 去拟合上述数据.
mathematica程序及运行结果如下:
data={10.15,5.52,4.08,2.85,2.11,1.62,1.41,1.30,1.21,1.15};
fit[data,{1,1/x},x]

四、实验内容与要求
1 画出实验问题的数据图,并粗略估计这些数据与什么类型的函数比较吻合?
2 取经验公式为线性函数y=ax+b 按照最小二乘法的原理用mathematica编程解实验问题.
3 取经验公式为y=ax+b +c sin[ x]+d cos[ x] ,用mathematica中算符fit[]来求解实验问题,并与内容2的精度比较,对比实际情况,你能得出什么?

五、操作提示
1

2 拟合程序及运行结果如下:

预测程序及运行结果如下:


3 程序及运行结果如下:





计算两种经验公式的精度可以看出第二种较好,这与客流量呈季节被动变动的实际情况吻合.
怎样用 mathematica 拟合二元函数?
数据拟合
由一组已知数据(xk,yk)(k=1,2,…,n),求函数的近似解析式y=f(x),就是数据拟合问题,当然函数还可以是多元的.
Mathematica提供了进行数据拟合的函数:
Fit[data,funs,vars] 对数据data用最小二乘法求函数表funs中各函数的一个线性组合作为所求的近似解析式,其中vars是自变量或自变量的表.
例如:
Fit[data,{1,x},x] 求形为y=a+bx的近似函数式.
Fit[data,{1,x,x2},x] 求形为y=a+bx+cx2的近似函数式.
Fit[data,{1,x,y,x y},{x,y}] 求形为z=a+bx+cy+dxy的近似函数式.
以上出现的参数data的格式为{{x1,y1,…,f1},{x2,y2,…,f2},…}.
函数表中的函数还可以是更复杂的初等函数.
例1 由下面给出的一组数据进行线性拟合,并绘制拟合曲线.
xi 19.1 25 30.1 36 40 15.1 50
yi 76.3 77.8 79.25 80.8 82.35 83.9 85.1
In[1]:=data={{19.1,76.3},{25,77.8},{30.1,79.25},{36,80.8},
{40,82.35},{45.1,83.9},{50,85.1}};
f=Fit[data,{1,x},x]
Out[2]=70.5723+0.291456x
In[3]:= pd=ListPlot[data,DisplayFunction→Identity];
fd=Plot[f,{x,19,52},DisplayFunction→Identity];
Show[pd,fd,DisplayFunction→$DisplayFunction]
图13-49 线性拟合的示意图
Out[5]=-Graphics-
说明:上例使用一次函数得到很理想的拟合,图形如图13-49所示.
例2 由下面给出的一组数据进行二次函数拟合,并绘制拟合曲线.
xi 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9
yi 5.1234 5.3057 5.5687 5.9378 6.4337 7.0978 7.9493 9.0253 10.3627
In[1]:= data={{0.1,5.1234},{0.2,5.3057},{0.3,5.5687},
{0.4,5.9378}, {0.5,6.4337},{0.6,7.0978},
{0.7,7.9493},{0.8,9.0253},{0.9,10.3627}};
f=Fit[data,{1,x,x^2},x]
Out[2]=5.30661-1.83196x+8.17149x2
In[3]:= pd=ListPlot[data,DisplayFunction→Identity];
fd=Plot[f,{x,0,1},DisplayFunction→Identity];
Show[pd,fd,DisplayFunction→$DisplayFunction]
图13-50 使用二次函数拟合的示意图
Out[5]= -Graphics-
以上两例都是计算方法教材中的习题,利用Mathematica可以轻而易举地得到答案,并同时画出图形以便直观地了解拟合的质量.
以下是二元拟合.
例3 观察下面的二元函数拟合.
In[1]:=Flatten[Table[{x,y,1 + 5x –x y},
{x,0,1,0.2},{y,0,1,0.2}],1]
Out[1]={{0,0,1},{0,0.2,1},{0,0.4,1},
{0,0.6,1},{0,0.8,1},{0,1.,1},
{0.2,0,2.},{0.2,0.2,1.96},{0.2,0.4,1.92},
{0.2,0.6,1.88},{0.2,0.8,1.84},{0.2,1.,1.8},
{0.4,0,3.},{0.4,0.2,2.92},{0.4,0.4,2.84},
{0.4,0.6,2.76},{0.4,0.8,2.68},{0.4,1.,2.6},
{0.6,0,4.},{0.6,0.2,3.88},{0.6,0.4,3.76},
{0.6,0.6,3.64},{0.6,0.8,3.52},{0.6,1.,3.4},
{0.8,0,5.},{0.8,0.2,4.84},{0.8,0.4,4.68},
{0.8,0.6,4.52},{0.8,0.8,4.36},{0.8,1.,4.2},
{1.,0,6.},{1.,0.2,5.8},{1.,0.4,5.6},
{1.,0.6,5.4},{1.,0.8,5.2},{1.,1.,5.}}
In[2]:=Fit[%,{1,x,y,x y},{x,y}]
Out[2]=1.+5. x+7.77156×10-16 y -1. x y
In[3]:=Chop[%]
Out[3]= 1.+ 5. x -1. x y
说明:在上例的In[1]中,首先生成二元函数1+5x-xy在0≤x≤1,0≤y≤1时的一个数据表,然后In[2]由这些数据反过来求二元函数,说明Fit可以求解多元问题.In[3]使用函数Chop去掉系数很小的项,以此消除误差.
函数Chop的一般形式为:
Chop[expr,δ] 去掉表达式expr的系数中绝对值小于δ的项,δ的默认值为10-10.
最后这个例子用于说明Fit的第二个参数中可以使用复杂的函数,不限于1,x,x2等基本类型.
例4 观察下面使用初等函数组合进行的拟合.
In[1]:= ft=Table[N[1+2Exp[-x/3]],{x,10}]
Out[1]={2.43306,2.02683,1.73576,1.52719,1.37775,
1.27067,1.19394,1.13897,1.09957,1.07135}
In[2]:=Fit[ft,{1,Sin[x],Exp[-x/3],Exp[-x]},x]
Out[2]= 1. -4.44089×10-15e-x +2.e-x/3+2.22045×10-16Sin[x]
In[3]:=Chop[%]
Out[4]=1. +2. e-x/3