《数学实验 第4版》课件 7.2多元回归分析

1第七章回归分析实验7.1一元回归分析实验7.2多元回归分析数学实验2实验7.2多元回归分析数学实验一、多元线性回归分析二、多元线性回归分析的matlab实现四、多元非线性回归分析三、逐步回归分析3一、多元线性回归分析

（7.3）称（7.3）多元线性回归模型，假设预测对象为y与m(m≥2)个影响因素之间有以下线性关系其中y称为因变量或响应变量，称为回归变量，

称为回归系数．

是随机误差，一般假设，是未知参数.

更一般地有，是未知的待定系数，（7.4）其中

是已知的函数，

也称为（广义）多元线性回归模型.4多元线性回归分析实验的主要任务是：1.用试验值（样本观测值）对待定系数做出估计；2.对建立的回归方程和每个回归变量进行显著性检验；3.利用回归方程对y作预测，或给定y值，对回归变量作控制．51.回归系数的估计为了估计回归系数，做了n组实验得到数据代入多元线性回归模型得到矩阵形式为其中矩阵X为已知的样本数据矩阵，称为资料矩阵；6B为未知的列向量(回归系数)；服从独立同分布，即应用最小二乘法估计可得到回归系数的估计值，

设为

可得（样本）多元线性回归方程：代入一组观测值

,通过回归方程可计算出,称之为回归预测值．72.线性回归方程的显著性检验多元线性回归模型的检验包括：(1)回归方程的显著性检验；(2)检验每个回归变量与因变量之间的相关程度；(3)检验模型本身是否存在自相关等.首先看回归方程的整体显著性检验方法.多元线性回归分析选取的统计量同实验7.1，此处略.（1）回归方程的拟合优度（R2检验法）R这里被称为复相关系数或全相关系数，即决定系数的算术平方根。

复相关系数R用来解释,这一组影响因素与y的线性相关程度.8R值越接近1，说明因变量y与回归变量之间的函数关系越密切；

反之，则说明因变量y与回归变量之间的函数关系不密切或不存在线性函数关系.通常R大于0.8（或0.9）才认为相关关系成立.多重判定系数R2在多元线性回归分析是度量多元回归方程拟合程度的一个统计量，反映了在因变量y的变差中被估计的回归方程所解释的比例.9给定显著水平，

通常通过F的统计值对应的概率P

<来说明因变量y与全体回归变量之间的当

时，称回归方程高度显著；当

时，称回归方程显著；当

时，称回归方程不显著.（2）F检验法构造统计量F则,F

服从第一自由度为m，第二自由度为n-m-1的F分布，线性相关性:103.每个变量的显著性检验

即使回归方程的整体线性相关性显著，仍然要对关每个回归变量与因变量之间的线性相关的显著性进行检验.常用t-检验法.(1)可以根据软件计算出的t统计量的统计值对应的概率大小来定性的判断该回归变量的显著性.(2)也可以通过对回归系数的置信区间来定性的判断每个变量的影响显著.如果某个变量的置信区间包含0点,则说明该变量对因变量影响不显著.(3)若存在不显著的变量，剔除后再进行其余变量的回归，直至余下的变量全部显著为止.方法如下:114.残差分析与模型诊断观测值与回归值之差称为残差.在回归模型定义中，我们假设随机误差不服从正态分布，则说明建立回归模型不够好，自相关等其他因素所致，需要

对于通过检验的模型，残差图中置信区间不经过0直线的残差所对应的个别数据，可从原数据中删除后再重新进行回归，这一点有时候很重要，直接影响到模型的结构，如例5中建立的模型.

除了利用残差外，由于剩余标准差

,如果残差也可能是数据存在进一步改进模型.12显然，s越接近0，说明回归预测值与原始数据拟合的越好.表示观测值偏离回归直线的平均误差，利用它也可以判断回归方程拟合的效果.5.预测

当我们获得显著的回归方程，就可以运用该回归方程进行分析预测了.给出自变量的一组观测值，代入回归方程即可得到的回归预测值.给定置信度，我们还可以得y的的预测区间（置信区间）.13函数功能regress(y,x,alpha)计算回归系数及其区间估计，残差及其置信区间，并g给出检验回归模型的参数（决定系数，F统计量等），alpha缺省为0.05rcoplot(r,rint)画出残差及其置信区间fitlm(x,y,modelspec)以x为数据矩阵，以y为相应变量，用modelspec的方式建立一个线性回归模型.Modelspec方式见软件说明，可缺省.plotDiagnostics(mdl,plottype)以plottype选项的方式显示数据与回归模型的数据诊断图plotResiduals(mdl,plottype)以plottype指定选项的方式显示数据与回归模型的误差图predict(mdl,Xnew)返回(线性、非线性)模型mdl在Xnew的预测值和99%置信区间相关matlab的命令函数二、多元线性回归分析的matlab实现实验7.2多元线性回归分析14例3

某公司调查某种商品的两种广告费用1和广告费用2对该产品销售量的影响，得到如下数据，试建立线性回归模型并进行检验，诊断是否有异常点.表7-8数据表销量Y9690959295959494广告费1(x1)1.52.01.52.53.32.34.22.5广告费2(x2)5.02.04.02.53.03.52.53.015（1）输入数据，观察x与y是否线性关系X=[1.5 2.0 1.5 2.5 3.3 2.3 4.2 2.5；5.02.04.0 2.5 3.0 3.5 2.5 3.0];Y=[96 90 95 92 95 95 94 94]’;plot(X(:,1),Y,’r*’)↙plot(X(:,2),Y,’k*’)↙用fitlm函数执行实验实验方法如下：在命令窗口输入:(2)求回归方程首先观察数据的散点图，在命令行窗口输入：16

如上图所示,Y与x1,x2“大致”呈线性关系，我们可首先建立线性回归模型.17在命令行窗口输入：dlmf=fitlm(X,Y)↙dlmf=

线性回归模型:y~1+x1+x2估计系数:EstimateSEtStatpValue(Intercept)83.2121.713948.557.0048e-08x11.29850.349243.71790.013742x22.33720.331137.05820.00088245观测值数目:8,误差自由度:5均方根误差:0.7R-方:0.909,调整R-方0.872统计模型（常量模型）:24.9,p-值=0.0025118因此，回归方程为y=83.212+1.2985+2.3372(3)显著性检验方程显著性的整体检验：统计量R2=0.909的数值较大，说明方程拟合数据较好.F=24.94，对应于F的概率P=0.00251<0.01,总体上说明模型线性整体线性相关(4)诊断分析剩余标准差s=0.7相对因变量的值较小.性高度显著.所以回归方程与原数据拟合的较好.进一步观察是否有异常数据.19在命令窗口输入：plotDiagnostics(dlmf,'cookd')↙如图可见第一个数据的残差大于平均值.20[~,larg]=max(dlmf.Diagnostics.CooksDistance);↙dlmf2=LinearModel.fit(X,Y,'Exclude',larg)↙dlmf2=（略）剔除该异常数据，重新回归：

可见各项统计指标R2和F检验概率都显著增大,剩余标准差变得更小，模型拟合效果得到进一步改善.注：

如果使用Regress命令，还可以利用rcoplot(r,rint)得到残差及置信区间的图.21

两种广告费用之间会不会有影响呢?大家可以深入考虑模型应该如何改进.也可以利用

LinearModel.fit(X,Y,'interactions')试一试.这里‘interactions’选项表示上述模型的基础上引入x1*x2交叉项.22例4.

某产品的收率Y（%）与处理压强X1(1.0e+5Pa)及温度x2(摄氏度)有关，测的实验数据如下;

请检验产品收率Y与处理压强X1及温度X2之间是否存在显著的线性相关关系;如果存在，求Y关于X1及X2的线性回归方程.X1X2YX1X2Y6.8665409.1700657.2685499.3680587.6690559.5685598700639.7700678.26956510650568.46705710.3690728.66755810.5670688.86906223（1）输入数据观察散点图在命令窗口输入：X1=[6.87.27.68.08.28.48.68.89.19.39.59.710.010.310.5]’;X2=[665685690700695670675690700680685700650690670]’;Y=[404955636557586265585967567268]’;plot(X1,Y,'*')，plot(X2,Y,'*')↙24（2）计算回归系数在命令窗口输入:X=[ones(15,1)X1X2];[b,bint,r,rint,stats]=regress(Y,X)↙b=-200.45545.68340.3075bint=-290.7876-110.12313.99987.36700.17900.4360……stats=0.862137.50020.000010.3172此处省略了r，rint参数的显示。25（3）回归方程的显著性检验（4）残差分析在命令窗口输入:rcoplot(r,rint)↙检验拟合效果，以进一步优化和改进模型。必要时可以剔除异常数据。

回归系数的置信区间都不包含0，统计量R2=0.8621数值较大,F=37.5002，P<0.05.说明模型线性相关性显著.线性回归方程结果为+0.3075y=-200.4554+5.683426

从残差图7-12可以看出，残差分布在0直线附近，且残差的置信区间均包含零点，分布正常，但是有的数据的残差较大，而且误差方差的估计值s2=10.3172较大，残差分布呈现一定的趋势性，模型有待进一步改进.请大家思考.图7.1327

此外，一元多项式回归和多元多项式回归也是比较常用的线性回归模型，相关MATLAB命令有确定多项式系数的命令polyfit，一元多项式回归交互式命令polytool和求预测及预测误差估计的polyval,polyconf，以及曲线拟合工具箱cftool等.多元二项式回归命令有rstool，具体用法可参与其他资料.

多元回归分析建模是一个复杂的过程，多元回归分析有着丰富的数学理论，有兴趣的同学可以进一步深入学习相关知识.28三、逐步回归分析

为了得到“最优”的回归模型，我们要保留对因变量影响大的变量，剔除对因变量影响小的变量.这里最有效的方法是逐步回归法.(1)从一个自变量开始，根据对因变量y的影响程度，从大到小地依次逐个引入回归方程.但当引入的自变量由于后面的自变量的引入而变得不明显时,要将其除掉.(2)每引入或剔除一个自变量，都要对y进行一次检验，以确保每次引入新变量前回归方程中只包含对y作用显著的变量．(3)这个过程反复进行，直至没有显著影响变量引入，也没有不显著影响变量剔除为止．29

通过观察决定系数R2、F统计量和剩余标准差（RMSE）、回归系数的区间的变化来判断该判定变量对模型的影响的显著性.可用剩余标准差（RMSE）最小作为衡量变量选择的一个数量标准.引入或剔除变量的依据：命令函数用法介绍stepwise（x，y,inmodel,alpha）x是自变量数据矩阵，y是因变量数据矩阵；alpha是显著性水平（缺省时为0.05），inmode是自变量初始集合的指标.stepwiselm(x,y,modelspec)x和y的意义同stepwise，modelspec用来提供模型的类别，详见matlab2012的软件帮助.30例5

表中数据是某建筑公司去年20个地区是销售量（Y千元），推销开支、实际帐目数、同类商品竞争数和地区潜力分别是影响建筑材料销售量的因素，试分析哪些是主要的影响因素，并建立该因素的线性回归模型.地区I推销开支（x1）实际帐目数(x2)同类商品竞争数(x3)地区销售潜力(x4)销售量y15.53110879.322.55586200.138.067129163.2…..........203.038815146.031利用MATLAB中逐步回归命令解决问题：(1)输入数据Y=[79.3200.1163.2200.1146.0177.730.9291.9160.0339.4159.686.3237.5107.2155.0201.4100.2135.8223.3195.0]’;X=[5.531.010.08.02.555.08.06.0

…………..7.059.09.011.0]；此处省略了部分数据(2)逐步回归

为简单起见，先在模型中考虑全部变量，然后逐个加入对y影响不显著的变量.在命令窗口输入：stepwise(X,Y)↙32点击左侧的按钮nextstep,X3的置信区间线段变为蓝色，表示变量X3已经进入模型中.如下图33再点击左侧的按钮nextstep，X2的置信区间线段变为蓝色.34在以上过程中,决定系数R2增大，F的值增加，对应的概率P值明显减小，剩余标准差RMSE明显减小.35

也可以直接点击右侧的Allsteps按钮一次完成变量的自动选择.

如果再加入变量X1，X4时，决定系数几乎没变化，但F的值明显减小了，RMSE的值增大了.这两个变量不宜留在模型中，因此可以忽略X1和X4对y的影响！因此X2、X3是Y的主要影响因素.下面使用另一个命令函数再做这个实验.

如果使用stepwiselm(x,y)命令，并不出现交互的窗口，而是在命令窗口直接动态地给出结果.36在命令窗口输入：mdl=stepwiselm(X,Y)↙1.正在添加x3,FStat=29.0275,pValue=4.04289e-052.正在添加x2,FStat=49.6984,pValue=1.95221e-06mdl=线性回归模型:y~1+x2+x3估计系数:EstimateSEtStatpValue(Intercept)186.0535.8435.19067.3688e-05x23.09070.438417.04971.9522e-06x3-19.5142.3915-8.15962.7862e-07观测值数目:20,误差自由度:17均方根误差:24R-方:0.902,调整R-方0.891统计模型（常量模型）:78.6,p-值=2.56e-0937（3）变量Y和X2、X3的回归方程

回归模型的常数项Intercept为186.0484，X2的系数为3.0907，X3的系数为-19.5140，故模型为

以上计算的结果并没有把变量X1引入模型，和我们的直观感觉不一致。X1表示推销的开支，推销应该是有利于销售的.下面我们对模型做进一步的诊断.

用X2和X3建立的上述回归模型中，剩余标准差s=24相对y的值来说较大.下面分析残差：(4)模型的分析与诊断38在命令窗口输入：plotResiduals(mdl,'probability')↙可见有一个残差小于-80(图省略)且严重偏离拟合直线.找到对应的数据：在命令窗口输入：find(mdl.Residuals.Raw<-80)↙ans=16排除第16个数据重新拟合：

在命令窗口输入：39stepwiselm(X,Y,'Exclude',16)↙1.Addingx3,FStat=29.4601,pValue=4.5237e-052.Addingx2,FStat=587.0139,pValue=4.882293e-143.Addingx1,FStat=5.9635,pValue=0.0274684.Addingx1:x3,FStat=8.3808,pValue=0.011759mdl=线性回归模型:（模型为：y~1+x2+x1*x3估计系数:EstimateSEtStatpValue(Intercept)135.9317.0157.98891.3918e-06x19.56172.69763.54460.0032357x23.44060.1106531.0932.5454e-14x3-16.5521.8891-8.76184.6777e-07x1:x3-0.852220.29438-2.8950.011759观测值数目:19,误差自由度:14均方根误差:5.54R-方:0.996,调整R-方0.994统计模型（常量模型）:802,p-值=2.33e-1640所以回归方程为：结果分析：

变量x1可以进入模型（P=0.03208<0.05），这与我们的常识相符合：推销总是有利于销售的；X1*X3的交互项也进入模型，说明推销开支和同类商品的竞争的交互作用和销售量y有线性相关性，这也符合我们的常识——该公司的推销必然会影响同类商品的销售，最终会反映到该公司的销售量上；

销售潜力X4未进入了模型，说明潜力和销量没有显著地线性关系.进一步的分析，留给读者自己去完成.通过这个例子也再次说明了残差分析的重要性.41注意：

在（4）这一步如果用stepwise函数得不到X1与X3的交互项.如果想要得到含交叉项的模型，就需要重新设计数据矩阵X，在其中加入该交叉项的数据列后再做回归.我们将在最后一个例子来具体说明.42四、多元非线性回归分析多元非线性回归模型可以表示为是随机变量（是向量）.其中y是因变量（是向量），Ｘ是回归变量（是向量或矩阵），是待定系数（是向量），

通过一个例子主要来说明利用Matlab软件进行非线性回归分析的参数估计方法以及利用回归方程进行预测的方法.43非线性回归Matlab命令函数函数功能nlinfit(x,y,’model’beta0)计算非线性回归的系数，残差，估计预测误差的数据nlintool(x,y,’model’,beta0,alpha)产生拟合曲线和y的置信区间等信息的交互画面nlpredci(‘model’,x,beta,r,J)求回归函数在x处的预测值y及其置信区间nlparci(beta,r,J)计算回归系数的置信区间fitnlm(x,y,fun,beta0)与nlinfit函数采用相同算法的另一个非线性回归命令plotDiagnostics(mdl,plottype)以plottype选项的方式显示数据与回归模型的数据诊断图plotResiduals(mdl,plottype)以plottype指定选项的方式显示数据与回归模型的误差图predict(mdl,Xnew)返回(线性、非线性)模型mdl在Xnew的预测值和99%置信区间实验7.1一元回归分析44例6

在研究化学动力学反应过程中，建立了一个反应速度和反应物含量的数学模型.其中是未知的参数，是三种反应物（氢，n戊烷，异构戊烷）的含量，y是反应速度.今测得一组数据如下表，试由此确定参数区间.，并给出其置信参数的参考值为（0.1，0.05，0.02，1，2）.序号反应速度y氢x1n戊烷x2异构戊烷x318.554703001023.7928580104534.8247030012040.024708012052.754708010614.391001901072.54100806584.3547019065913.0010030054108.50100300120110.05100801201211.3228530010133.1328519012046下面用fitnlm来处理上述非线性回归问题：（1）首先以回归系数和自变量为输入变量，将要拟合的模型写成函数huaxue.m文件functionyhat=huaxue(beta,x);yhat=(beta(4)*x(:,2)-x(:,3)/beta(5))./(1+beta(1)*x(:,1)+...beta(2)*x(:,2)+beta(3)*x(:,3));（2）输入数据在命令窗口输入：x0=[18.5547030010;23.792858010;34.82470300120;40.0247080120;52.754708010;614.3910019010;72.541008065;84.3547019065;913.0010030054;108.50100300120;110.0510080120;1211.3228530010;133;47x=x0(:,3:5);y=x0(:,2);beta0=[0.1,0.05,0.02,1,2];%回归系数的初值（3）回归系数计算与分析nlfit=fitnlm(x,y,@huaxue,beta0)↙nlfit=Non线性回归模型:y~huaxue(b,X)估计系数:EstimateSEtStatpValueb10.0627760.043561.44120.18751b20.0400480.0308841.29670.23087b30.112420.0751551.49580.17308b41.25260.866991.44480.18652b51.19140.836711.42390.192348观测值数目:13,误差自由度:8均方根误差:0.193R-方:0.999,调整R-方0.998F-statisticvs.zeromodel:3.91e+03,p-值=2.54e-13

从以上结果可得非线性回归的系数，与剩余标准差RMSE=0.193.下面分析拟合的效果：

下面分析拟合的效果,plotResiduals(nlfit,'fitted')↙在命令窗口输入：49生右图.可见残差分布在0均值附近，且数值较小说明方程拟合效果较好.50（4）回归方程的预测：[yhat,yci]=predict(nlfit,x)；↙plot(x0(:,1),y,'r+',x0(:,1),yhat,'bo')图7-164可以看出，原始数据与回归预测值拟合的较好.（+~原始数据，o~拟合结果）51（5）交互方式工具图

用plotSlice(nlfit)函数得到一个交互式画面图调节变量的值可得到相应的结果，如y的预测值.左下方的Export可向