数学建模竞赛中常用软件的操作(3)

hold on plot(time,f) hold off

输出7多项次的系数向量为：

p =

0.0002 -0.0113 0.2817 -3.6307 25.5345 -94.6710 153.8757 -1.2592

如果要求血液中酒精浓度的最大值，可用以下指令： roots(polyder(p))

重庆文理学院的一个参赛队于2004年关于“饮酒驾车”问题获得一等奖后，研制了“酒后安全驾车时刻表”和“人体内酒精浓度反推软件”两个产品，并在重庆永川交警大队得到使用。

Matlab提供了两个非线性最小二乘拟合函数：lsqcurvefit和lsqnonlin。两个命令都要先建立M-文件fun.m，在其中定义函数f(x)，但两者定义f(x)的方式是不同的。

（2）lsqcurvefit函数

该函数用来进行非线性拟合，其调用格式为：

x = lsqcurvefit (?fun‘,x0,xdata,ydata,options);

其中，fun为事先建立的拟合函数F(x,xdata)，其中自变量x表示拟合函数中的待定参数，xdata为已知拟合节点的x坐标，x0为待定参数x的迭代初始值，xdata,ydata为已知数据点的x和y坐标，options是一些控制参数。

lsqcurvefit函数用来求含参数x（向量）的向量值函数

F(x,xdata)={f(x,data1), f(x,data2),…, f(x,datan)}

中的参数x（向量），使得

?[f(x,xdata)?ydata]iii?1n2

最小。

例如，根据表5.22中的数据，利用lsqcurvefit函数拟合y(x)?a?be表5已知数据点 x y×10-3 100 4.54 200 4.99 300 5.35 400 5.65 500 5.90 600 6.10 700 6.26 800 6.39 900 6.50 1000 6.59 ?0.02kx。

首先建立拟合函数的M文件myfit1.m，其内容如下： function f = myfit1(x,xdata)

f = x(1)+x(2)*exp(-0.02*x(3)*xdata); 其中x(1),x(2),x(3)分别表示拟合曲线中的参数a,b,k。

然后在命令窗口输入： >> xdata = 100:100:1000;

ydata = 1e-03*[4.54,4.99,5.35,5.65,5.90,6.10,6.26,6.39,6.50,6.59]; x0=[0.2,0.05,0.05];

x = lsqcurvefit ('myfit1',x0,xdata,ydata) f= myfit1(x,xdata) （3） lsqnonlin函数

该函数用来进行非线性拟合，其调用格式为：

x = lsqnonlin (?fun‘,x0, options);

其中fun为事先建立的拟合函数f(x)，其中自变量x表示拟合函数中的待定参数，x0为待定参数x的迭代初始值，options是一些控制参数。

由于lsqnonlin中定义的拟合函数的自变量是x，所以已知参数xdata,ydata应写在该函数中。

lsqnonlin函数用来求含参量x（向量）的向量值函数

f(x)={f1(x), f2(x),…，fn(x)}

中的参量x，使得

?fi?1n2i(x)最小。其中fi(x)=f(x,xdatai,ydaytai)=F(x, xdatai)- ydatai。

例如，根据表5中的数据，利用lsqnonlin函数拟合y(x)?a?be?0.02kx。 首先建立拟合函数的M文件myfit2.m，其内容如下： function f=myfit2(x) xdata=100:100:1000;

ydata=1e-03*[4.54,4.99,5.35,5.65,5.90,6.10,6.26,6.39,6.50,6.59]; f=x(1)+x(2)*exp(-0.02*x(3)*xdata)- ydata; 然后在命令窗口输入：

>> x0=[0.2,0.05,0.05]; x = lsqnonlin('myfit2',x0) f = myfit2(x)

该结果与lsqcurvefit函数拟合的结果相同。

3. 回归分析

MATLAB中提供了一些线性和非线性回归分析函数。 （1） regress函数 一般地，称

?Y?X??? ?2?E(?)?0,COV(?,?)??In为高斯—马尔柯夫线性模型（k元线性回归模型）。

?y1??1x11?...??1x21Y???，X???...??......???y?n??1xn1x12x22...xn2...x1k???0???1???????...x2k?1?，????，???2? ?...??...?......??????...xnk???k???n?Matlab提供了多元线性回归函数regress，采用的是最小二乘估计，其调用格式有：

b = regress (y,x)

返回值为线性模型y = x*b的回归系数向量。其中x为n×(k+1)矩阵，行对应于观测值，列对应于预测变量，y为n×1 向量，为因变量，一元线性回归可取k=1。

[b,bint,r,rint,stats] = regress(y,x,alpha)

其中bint是回归系数的区间估计，r是残差，rint是置信区间，stats是用于检验回归模型的统计量，有四个数值：相关系数r^2，F值，与F对应的概率P和误差方差估计，alpha是显著性水平（缺省的时候为0.05）。相关系数r^2越接近1，说明回归方程越显著；与F对应的概率P

例如，线性回归模型y= a+b*k+c*l+d*t输入

>> y=[8.8818 8.9487 9.0541 9.1545 9.2693 9.4289 9.6160 9.8150 9.9825 10.1558 10.3193]';

k= [7.8381 7.9167 8.0048 8.1026 8.2556 8.5822 8.8287 9.0756 9.2175 9.4148 9.6198]'; l= [8.3871 8.3872 8.3935 8.3971 8.4025 8.4048 8.4079 8.4141 8.4261 8.4377 8.4444]'; t= [9.9551 9.9057 10.0972 9.9537 9.9370 9.9449 9.9636 10.1291 10.1573 10.2944 10.2093]';

x= [ones(size(k)) k l t]; [b,bint,r,rint,stats]=regress(y,x) 输出结果为

b = -55.4988 0.5644 7.1254 0.0222 stats = 1.0e+003 *

0.0010 2.1022 0.0000 0.0000

该结果说明y=-55.4988+0.5644k+7.1254l+0.0222t，stats中的数据说明r2=1，F=2102.2，p=0，由于p<0.05可知回归模型成立。 （2） rstool函数

该函数是多元二项式回归函数，其调用格式为

rstool(x,y, 'model',alpha)

其中x为n×m为矩阵，y为n维列向量（n为数据点数，m为元数），'model'为以下4种模型：

'linear'（线性，缺省）：y??0??1x1????mxm； 'interaction'（交叉）：y??0??1x1????mxm?1?j?k?m??jkxjxk；

'quadratic'（完全二次）：y??0??1x1????mxm?'purequadratic'（纯二次）：y??0??1x1????mxm?alpha为显著性水平，默认值为0.05。

1?j,k?m?n?jkxjxk；

x2j。

??j?1jj例如，设某商品的需求量与消费者的平均收入、商品价格的统计数据如表5.3，建立多元二项式纯二次回归模型，并预测平均收入为1000、价格为6时的商品需求量。

表6需求量、平均收入和价格统计表 需求量 收入 价格

100 1000 5

75 600 7

80 1200 6

70 500 6

50 300 8

65 400 7

90 1300 5

100 1100 4

110 1300 3

60 300 9

可以直接使用多元二项式回归，在命令窗口输入： >> x1=[1000 600 1200 500 300 400 1300 1100 1300 300]; x2=[5 7 6 6 8 7 5 4 3 9];

y=[100 75 80 70 50 65 90 100 110 60]'; x=[x1' x2'];

rstool(x,y,'purequadratic') 其输出如图2。

图2多元二项式回归

在图2中x1上面的方框中输入1000，x2上面的方框中输入6，在图形框左侧的―Predicted

Y1‖下方的数据变为88.47981，即预测出平均收入为1000，价格为6时的商品需求量为88.4791。

单击图形框左边Export，则出现图3对话框，可以将回归参数beta、剩余标准差rmse和残差residuals传送到Matlab的工作区中。

图3输出对话框

在命令窗口输入： >> beta rmse

22可知故回归模型为：y?110.5313?0.1464x1?26.5709x2?0.0001x1，?1.8475x2剩余标准差为4.5362, 说明此回归模型的显著性较好。

还可将该模型转化为多元线性回归模型，利用regress函数进求解。 可输入：

>> x1=[1000 600 1200 500 300 400 1300 1100 1300 300]; x2=[5 7 6 6 8 7 5 4 3 9];

X=[ones(10,1) x1' x2' (x1.^2)' (x2.^2)']; [b,bint,r,rint,stats]=regress(y,X); b,stats

输出的回归参数和rstool函数的结果相同。

（3） nlinfit函数

nlinfit函数用来确定非线性回归系数，调用格式为：

[beta, r, J] = nlinfit (x, y, 'modelfun', beta0)

其中输入数据x,y分别为n×p维矩阵和n维列向量，对于一元非线性回归，取p=1即可；'modelfun'为事先定义的非线性回归函数的M文件，是回归系数beta和x的函数；beta0是回归系数的初值，输出参数beta是估计出的回归系数，r为残差，J为Jacobain矩阵。

例如，根据表5中的数据，利用nlinfit函数进行非线性回归，回归函数为

y(x)?a?be?0.02kx。

首先建立回归函数的M文件myfit3.m，内容如下： function f = myfit3(beta,xdata)

f = beta(1)+beta(2)*exp(-0.02*beta(3)*xdata);

百度搜索“77cn”或“免费范文网”即可找到本站免费阅读全部范文。收藏本站方便下次阅读，免费范文网，提供经典小说综合文库数学建模竞赛中常用软件的操作(3)在线全文阅读。