
《第三讲 回归分析》由会员分享,可在线阅读,更多相关《第三讲 回归分析(90页珍藏版)》请在文档大全上搜索。
1、第第3章章 回归分析回归分析 MATLAB数据分析方法数据分析方法 第第3章章 回归分析回归分析 回归分析是最常用的数据分析方法之一。它回归分析是最常用的数据分析方法之一。它是根据已得的试验结果以及以往的经验来建立统是根据已得的试验结果以及以往的经验来建立统计模型,并研究变量间的相关关系,建立起变量计模型,并研究变量间的相关关系,建立起变量之间关系的近似表达式即经验公式,并由此对相之间关系的近似表达式即经验公式,并由此对相应的变量进行预测和控制等应的变量进行预测和控制等.3.1一元回归模型一元回归模型 3.1.1一元线性回归模型一元线性回归模型1.一元线性回归的基本概念一元线性回归的基本概念第
2、第3章章 回归分析回归分析 通常,我们对总体通常,我们对总体(x,Y)进行进行n次的独立观测,获得次的独立观测,获得n组数据(称为样本观测值)组数据(称为样本观测值) (x1,y1),(x2,y2),(xn,yn)利用最小二乘法可以得到回归模型参数利用最小二乘法可以得到回归模型参数 0, 1的最的最小二乘估计小二乘估计 设设Y是一个可观测的随机变量,它受到一个非随机变是一个可观测的随机变量,它受到一个非随机变量因素量因素x和随机误差和随机误差 的影响。若的影响。若Y与与x有如下线性关系:有如下线性关系:,10 xY(3.1.1)且且E =0,D = 2,则称则称(3.1.1)为一元线性回归模型
3、为一元线性回归模型.其中其中 0, 1为回归系数,为回归系数,x为自变量,为自变量,Y为因变量为因变量.10,第第3章章 回归分析回归分析 (3.1.2) 其中其中 于是建立经验公式模型:于是建立经验公式模型:(3.1.3)一元线性回归分析的主要任务:一是利用样本观测值一元线性回归分析的主要任务:一是利用样本观测值对回归系数对回归系数 0, 1和和 作点估计;二是对方程的线性关作点估计;二是对方程的线性关系即系即 1作显著性检验;三是在作显著性检验;三是在x=x0处对处对Y作预测等作预测等.以下举例说明建立经验公式(以下举例说明建立经验公式(3.1.3)的方法。)的方法。.,110 xxxyL
4、Lxy,niixnx11,niiyny11,niixxxxL12)()( )(1yyxxLiniixy01 yx第第3章章 回归分析回归分析 例例3.1.1 近近10年来,某市社会商品零售总额与职工工年来,某市社会商品零售总额与职工工资总额(单位:亿元)数据如下表资总额(单位:亿元)数据如下表3.1。表表3.1 商品零售总额与职工工资表商品零售总额与职工工资表 (单位:亿元)(单位:亿元)建立社会商品零售总额与职工工资总额数据的回归模型建立社会商品零售总额与职工工资总额数据的回归模型工资总额工资总额23.827.631.632.433.734.943.252.863.873.4零售总额零售总额
5、41.451.861.767.968.777.595.9137.4155.0175.0解:解:% 首先输入数据首先输入数据x=23.80,27.60,31.60,32.40,33.70,34.90,43.20,52.80,63.80,73.40;y=41.4,51.8,61.70,67.90,68.70,77.50,95.90,137.40,155.0,175.0;第第3章章 回归分析回归分析 % 然后作散点图然后作散点图plot(x,y,*) %作散点图作散点图xlabel(x(职工工资总额职工工资总额) %横坐标名横坐标名ylabel(y(商品零售总额商品零售总额) %纵坐标名纵坐标名20
6、304050607080406080100120140160180 x(职 工 工 资 总 额 )y(商品零售总额)图图3.1商品零售总额与职工工资总额数据散点图商品零售总额与职工工资总额数据散点图第第3章章 回归分析回归分析 % 计算最佳参数计算最佳参数Lxx=sum(x-mean(x).2);Lxy=sum(x-mean(x).* (y-mean(y);b1=Lxy/Lxx;b0=mean(y)-b1*mean(x);运行后得到:运行后得到:b1 = 2.7991,b0 = -23.5493所以,回归模型为所以,回归模型为 5493.237991. 2xy问题问题1:当:当x=0,得到,得
7、到y=-23.5493亿元如何理解?亿元如何理解?问题问题2:如何检验:如何检验E =0? D = 2?第第3章章 回归分析回归分析 2. 一元多项式回归模型一元多项式回归模型1110.nnnnya xaxa xa在一元回归模型中,如果变量在一元回归模型中,如果变量y与与x的关系是的关系是n次多次多项式,即项式,即其中其中 是随机误差,服从正态分布是随机误差,服从正态分布N(0, 2)a0,a1,an为回归系数为回归系数,则称则称(3.1.4)为多项式回归模型为多项式回归模型. (3.1.4)(1)多项式曲线拟合多项式曲线拟合在在MATLAB7的统计工具箱中,有多项式曲线拟合的的统计工具箱中,
8、有多项式曲线拟合的命令命令polyfit,其调用格式有以下三种:,其调用格式有以下三种:第第3章章 回归分析回归分析 p=polyfit(x,y,n) p,S=polyfit(x,y,n) p,S,mu=polyfit(x,y,n)其中,输入其中,输入x,y分别为自变量与因变量的样本观测数据分别为自变量与因变量的样本观测数据向量;向量;n是多项式的阶数,对于一元线性回归则取是多项式的阶数,对于一元线性回归则取n=1;输出输出p是按照降幂排列的多项式的系数向量,是按照降幂排列的多项式的系数向量,S是一个是一个矩阵,用于估计预测误差或供矩阵,用于估计预测误差或供MATLAB的其它函数的的其它函数的
9、调用调用 。例例3.1.2 某种合金中的主要成分为某种合金中的主要成分为A,B两种金属,经过两种金属,经过试验发现:这两种金属成分之和试验发现:这两种金属成分之和x与合金的膨胀系数与合金的膨胀系数y有如下关系,建立描述这种关系的数学表达式有如下关系,建立描述这种关系的数学表达式.第第3章章 回归分析回归分析 表表3.2 合金的膨胀系数表合金的膨胀系数表解:解:%首先输入数据首先输入数据x=37:0.5:43; y=3.4,3,3,2.27,2.1,1.83,1.53,1.7,1.8,1.9,2.35,2.54,2.9; %其次做散点图其次做散点图plot(x,y,*)xlabel(x(两种合金
10、之和两种合金之和) %横坐标名横坐标名ylabel(y(合金膨胀系数合金膨胀系数) %纵坐标名纵坐标名 %然后根据散点图猜测曲线类别然后根据散点图猜测曲线类别(2.1.7) x3737.5 3838.5 3939.5 4040.5 4141.5 4242.5 43y3.4332.27 2.11.83 1.53 1.71.81.92.35 2.54 2.9第第3章章 回归分析回归分析 由于散点图呈抛物线,故选择二次函数曲线进行拟合由于散点图呈抛物线,故选择二次函数曲线进行拟合.p = polyfit(x,y,2) %注意取注意取n=2运行得到回归系数:运行得到回归系数:p=0.1660 -13.
11、3866 271.6231 即二次回归模型为:即二次回归模型为:20.16613.3866271.6231yxx第第3章章 回归分析回归分析 多项式曲线拟合预测的命令多项式曲线拟合预测的命令polyval,其调用格,其调用格式有以下两种:式有以下两种:Y=polyval(p,x0)Y,Delta=polyconf(p,x0,S,alpha)其中,输入其中,输入p,S是由多项式拟合命是由多项式拟合命p,S=polyfit(x,y,n)的输出,的输出,x0是要预测的自变量的值是要预测的自变量的值.输出输出Y是是polyfit所所得的回归多项式在得的回归多项式在x处的预测值。处的预测值。 (2) 多
12、项式回归的预测与置信区间多项式回归的预测与置信区间如果输入数据的误差相互独立,且方差为常数,则如果输入数据的误差相互独立,且方差为常数,则YDelta至少包含至少包含95%的预测值;的预测值;alpha缺省时为缺省时为0.05。(Y-Delta, Y+Delta)即即95%的置信区间的置信区间第第3章章 回归分析回归分析 (3) 多项式回归的多项式回归的GUI界面命令界面命令多项式回归的多项式回归的GUI界面命令界面命令polytool,其典型调用格式,其典型调用格式 polytool(x,y,n,alpha)其中,输入其中,输入x,y分别为自变量与因变量的样本观测数据分别为自变量与因变量的样
13、本观测数据向量;向量;n是多项式的阶数;置信度为是多项式的阶数;置信度为(1-alpha)%,alpha缺省时为缺省时为0.05。 该命令可以绘出总体拟合图形以及该命令可以绘出总体拟合图形以及(1-alpha)上、下置信区间的直线(屏幕上显示为红色)上、下置信区间的直线(屏幕上显示为红色).此此外,用鼠标拖动图中纵向虚线,就可以显示出对于外,用鼠标拖动图中纵向虚线,就可以显示出对于不同的自变量数值所对应的预测状况,与此同时图不同的自变量数值所对应的预测状况,与此同时图形左端数值框中会随着自变量的变化而得到的预报形左端数值框中会随着自变量的变化而得到的预报数值以及数值以及(1-alpha) 置信
14、区间长度一半的数值。置信区间长度一半的数值。第第3章章 回归分析回归分析 第第3章章 回归分析回归分析 例例3.1.3为了分析为了分析X射线的杀菌作用,用射线的杀菌作用,用200千伏的千伏的X射线来照射细菌,每次照射射线来照射细菌,每次照射6分钟用平板计数法估分钟用平板计数法估计尚存活的细菌数,照射次数记为计尚存活的细菌数,照射次数记为t,照射后的细,照射后的细菌数菌数y如表如表3.3所示。所示。t123456789101112131415y3522111971601421061046056383632211915表表3.3 X射线照射次数与残留细菌数射线照射次数与残留细菌数试求:试求: 给出
15、给出y与与t的二次函数回归模型;的二次函数回归模型; 在同一坐标系内做出原始数据与拟合结果的散点图在同一坐标系内做出原始数据与拟合结果的散点图 预测预测t=16时残留的细菌数;时残留的细菌数; 根据问题实际意义选择多项式函数是否合适?根据问题实际意义选择多项式函数是否合适?数据来源:数据来源:http/www.ilr.cornell.edu/hadi/RABE第第3章章 回归分析回归分析 解:解:% 输入原始数据输入原始数据t=1:15;y=352,211,197,160,142,106,104,60,56,38,36,32,21,19,15;p=polyfit(t,y,2); % 作二次多项
16、式回归作二次多项式回归y1= polyval(p,t); % 模型估计与作图模型估计与作图plot(t,y,-*,t,y1,-o); legend(原始数据原始数据,二次函数二次函数) xlabel(t(照射次数照射次数) ylabel(y(残留细菌数残留细菌数)t0=16;yc1= polyconf(p,t0) % 预测预测t0=16时残留的细菌数时残留的细菌数运行结果为运行结果为p =1.9897 -51.1394 347.8967,yc1 =39.0396即二次回归模型为即二次回归模型为21 1.9897t -51.1394t+347.8967y 第第3章章 回归分析回归分析 yc1 =
17、39.0396,表明照射表明照射16次后,用二次函数计算出次后,用二次函数计算出细菌残留数为细菌残留数为39.0396,显然与实际不相符合。,显然与实际不相符合。调用多项式回归的调用多项式回归的GUI界面命令界面命令polytool,如图如图3.4原始数据与拟合结果的散点图如图原始数据与拟合结果的散点图如图3.3所示,从所示,从图形可知拟合效果较好图形可知拟合效果较好.图图 3.3 原始数据与拟合结果的散点图原始数据与拟合结果的散点图第第3章章 回归分析回归分析 根据实际问题的意义可知:尽管二次多项式拟合效根据实际问题的意义可知:尽管二次多项式拟合效果较好,但是用于预测并不理想。因此如何根据原
18、果较好,但是用于预测并不理想。因此如何根据原始数据散点图的规律,选择适当的回归曲线是非常始数据散点图的规律,选择适当的回归曲线是非常重要的,因此有必要研究非线性回归分析重要的,因此有必要研究非线性回归分析. 图图 3.4 二次函数预测交互图二次函数预测交互图第第3章章 回归分析回归分析 3.1.2一元非线性回归模型一元非线性回归模型 为了便于正确地选择合适的函数进行回归分析为了便于正确地选择合适的函数进行回归分析建模,我们给出通常选择的六类曲线如下所示:建模,我们给出通常选择的六类曲线如下所示:1. 非线性曲线选择非线性曲线选择(1)双曲线)双曲线1/y=a+b/x(见图见图3.5)。图图3.
19、5双曲线双曲线图图3.5双曲线双曲线第第3章章 回归分析回归分析 (2) 幂函数曲线幂函数曲线y=axb, 其中其中x0,a0(图图3.6)。图图3.6 幂函数曲线幂函数曲线(3)指数曲线)指数曲线y=aebx,其中参数,其中参数a0(见图见图3.7)。图图3.7 指数曲线指数曲线第第3章章 回归分析回归分析 xbaey/(4)倒指数曲线)倒指数曲线 ,其中a0(图3.8)。图图3.8 倒指数曲线倒指数曲线(5)y=a+blnx (见图见图3.9)。图图3.9 对数曲线对数曲线第第3章章 回归分析回归分析 xbeay1(6)S型曲线型曲线 (见图3.10)。图图3.10 S型曲线型曲线 对于非
20、线性回归建模通常有两种方法:一是通过适对于非线性回归建模通常有两种方法:一是通过适当的变换转化为线性回归模型,例如双曲线模型当的变换转化为线性回归模型,例如双曲线模型(图图3.5)。如果无法实现线性化,可以利用最小二乘法直接建立。如果无法实现线性化,可以利用最小二乘法直接建立非线性回归模型,求解最佳参数。非线性回归模型,求解最佳参数。第第3章章 回归分析回归分析 2.非线性回归的非线性回归的MATLAB命令命令MATLAB统计工具箱中实现非线性回归的命令有统计工具箱中实现非线性回归的命令有nlinfit、nlparci、lpredci和和nlintool。下面逐一介绍。下面逐一介绍调用格式。调
21、用格式。 非线性拟合命令非线性拟合命令nlinfit,调用格式:,调用格式:beta,r,J = nlinfit(x,y,model,beta0)其中,输人数据其中,输人数据x,y分别为分别为nm矩阵和矩阵和n维列向量,维列向量,对一元非线性回归,对一元非线性回归,x为为n维列向量,维列向量,model是事先用是事先用M文件定义的非线性函数,文件定义的非线性函数,beta0是回归系数的初值是回归系数的初值(需要通过解方程组得到需要通过解方程组得到),beta是估计出的最佳回归系是估计出的最佳回归系数,数,r是残差,是残差,J是是Jacobian矩阵,它们是估计预测误矩阵,它们是估计预测误差需要
22、的数据。差需要的数据。第第3章章 回归分析回归分析 非线性回归预测命令非线性回归预测命令nlpredci,调用格式:,调用格式: ypred = nlpredci(FUN,inputs,beta,r,J)其中,输入参数其中,输入参数beta,r,J是非线性回归命令是非线性回归命令nlinfit的输的输出结果出结果, FUN 是拟合函数,是拟合函数,inputs是需要预测的自变是需要预测的自变量;输出量量;输出量ypred是是inputs的预测值。的预测值。非线性回归置信区间命令非线性回归置信区间命令nlparci,调用格式:,调用格式: ci = nlparci(beta,r,J,alpha)
23、输入参数输入参数beta,r,J就是非线性回归命令就是非线性回归命令nlinfit输出的输出的结果,输出结果,输出ci是一个矩阵,每一行分别为每个参数的是一个矩阵,每一行分别为每个参数的(1-alpha)% 的置信区间,的置信区间,alpha缺省时默认为缺省时默认为0.05.第第3章章 回归分析回归分析 非线性回归的非线性回归的GUI界面命令界面命令nlintool,典型调用格式,典型调用格式 nlintool(x,y,fun,beta0)其中参数其中参数x,y,fun,beta0与命令与命令nlinfit中的参数含义相同中的参数含义相同. 例例3.1.4. 在在M文件中建立函数文件中建立函数
24、y=a(1-be-cx),其中,其中a,b,c为待定的参数。为待定的参数。解:解:fun=inline(b(1)*(1-b(2)*exp(-b(3)*x),b,x);此处,将此处,将b看成参变量,看成参变量,b(1),b(2),b(3)为其分量为其分量.例例3.1.5 炼钢厂出钢时所用盛钢水的钢包,由于钢水对炼钢厂出钢时所用盛钢水的钢包,由于钢水对耐火材料的侵蚀,容积不断增大,我们希望找出使用耐火材料的侵蚀,容积不断增大,我们希望找出使用次数与增大容积之间的函数关系次数与增大容积之间的函数关系.实验数据如表实验数据如表3.4。 第第3章章 回归分析回归分析 使用次数使用次数(x)2345678
25、9增大容积增大容积(y)6.428.29.589.59.7109.939.99使用次数使用次数(x)10111213141516增大容积增大容积(y)10.4910.5910.610.810.610.910.76表表3.4 钢包使用次数与增大容积钢包使用次数与增大容积(1)建立非线性回归模型)建立非线性回归模型1/y=a+b/x;(2)预测钢包使用)预测钢包使用x0=17次后增大的容积次后增大的容积y0;(3)计算回归模型参数的)计算回归模型参数的95%的置信区间。的置信区间。MATLAB脚本程序如下脚本程序如下:x=2:16;y=6.42,8.2,9.58,9.5,9.7,10,9.93,9
26、.99,10.49,10.59,10.6,10.8,10.6,10.9,10.76;第第3章章 回归分析回归分析 %建立非线性双曲线回归模型建立非线性双曲线回归模型b0=0.084,0.1436; % 初始参数值初始参数值fun=inline(x./(b(1)*x+b(2),b,x); beta,r,J=nlinfit(x,y,fun,b0);beta % 输出最佳参数输出最佳参数y1=x./(0.0845*x+0.1152); % 拟合曲线拟合曲线plot(x,y,*,x,y1,-or)legend(原始数据原始数据,拟合曲线拟合曲线)注意:初始值要先计算后,才能得到上面程序中的注意:初始值
27、要先计算后,才能得到上面程序中的b0,由于确定两个参数值,因此我们选择已知数据,由于确定两个参数值,因此我们选择已知数据中的两点(中的两点(2,6.42)和()和(16,10.76)代入设定方程,)代入设定方程,得到方程组得到方程组第第3章章 回归分析回归分析 上述方程组有两种解法:手工方法与上述方程组有两种解法:手工方法与Matlab方法。方法。下面用下面用Matlab方法解方程组:方法解方程组:a,b=solve(6.42*(2*a+b)=2,10.76*(16*a+b)=16)输出为输出为a =.83961597702347450462657355615004e-1b =.1436032
28、843460839152740622358104926.426.42(2)221610.76(16)1610.7616abababab第第3章章 回归分析回归分析 图图3.11钢包使用次数与增大容积的非线性拟合图钢包使用次数与增大容积的非线性拟合图在例在例3.1.5中,预测钢包使用中,预测钢包使用17次后增大的容积,可次后增大的容积,可在执行上面的程序中,继续输入命令在执行上面的程序中,继续输入命令ypred=nlpredci(fun,17,beta,r,J)得到:得到:ypred =10.9599即钢包使用即钢包使用17次后增大的容积次后增大的容积10.9599。第第3章章 回归分析回归分析
29、 求回归模型参数的求回归模型参数的95%的置信区间,只要继续添加程序的置信区间,只要继续添加程序ci = nlparci(beta,r,J)运行后得到运行后得到ci =0.0814 0.0876 0.0934 0.1370即回归模型即回归模型 中参数中参数a,b的的95%的置信区间分别为的置信区间分别为(0.0814 ,0.0876) 与(与(0.0934,0.1370).我们求出的最佳参数分别为我们求出的最佳参数分别为 a=0.0845,b=0.1152均属均属于上述置信区间。于上述置信区间。第第3章章 回归分析回归分析 图图3.12给出钢包使用次数与增大容积的非线性拟合的给出钢包使用次数与
30、增大容积的非线性拟合的交互图形,图中的的圆圈是实验的原始数据点,两条交互图形,图中的的圆圈是实验的原始数据点,两条虚线为虚线为95%上、下置信区间的曲线(屏幕上显示为红上、下置信区间的曲线(屏幕上显示为红色),中间的实线(屏幕上显示为绿色)是回归模型色),中间的实线(屏幕上显示为绿色)是回归模型曲线,纵向的蓝色虚线显示了在自变量为曲线,纵向的蓝色虚线显示了在自变量为8.9502,横,横向的虚线给出了对应的预测值为向的虚线给出了对应的预测值为10.2734.图图3.12 钢包使用次数与增大容积的非线性拟合交互图钢包使用次数与增大容积的非线性拟合交互图第第3章章 回归分析回归分析 例例3.1.6
31、对例题对例题3.1.3进行非线性回归,并预测照射进行非线性回归,并预测照射16次后细菌残留数目,给出模型参数的次后细菌残留数目,给出模型参数的95%的置信的置信区间,绘出模型交互图形区间,绘出模型交互图形.解:我们选取函数解:我们选取函数y=aebt进行非线性回归,该方程的两进行非线性回归,该方程的两个参数具有简单的物理解释,个参数具有简单的物理解释,a表示实验开始时的细菌表示实验开始时的细菌数目,数目,b表示细菌死亡(或衰变)的速率。表示细菌死亡(或衰变)的速率。MATLAB脚本程序如下:脚本程序如下:t=1:15;y=352 211 197 160 142 106 104 60 56 38
32、 36 32 21 19 15;fun=inline(b(1)*exp(b(2)*t),b,t) % 非线性函数非线性函数beta0=148,-0.2; % 参数初始值参数初始值beta,r,J=nlinfit(t,y,fun,beta0); % 非线性拟合非线性拟合beta % 输出最佳参数输出最佳参数y1=nlpredci(fun,t,beta,r,J); % 模型数值计算模型数值计算第第3章章 回归分析回归分析 plot(t,y,*,t,y1,-or),legend(原始数据原始数据,非线性回归非线性回归)xlabel(t(照射次数照射次数) ylabel(y(残留细菌数残留细菌数)yp
33、red = nlpredci(fun,16,beta,r,J) % 预测残留细菌数预测残留细菌数ci = nlparci(beta,r,J) % 参数参数95%区间估计区间估计nlintool(t,y,fun,beta0) % 作出交互图形作出交互图形运行后结果如下运行后结果如下:beta = 400.0904 -0.2240即即,最佳参数为最佳参数为:a=400.0904,b=-0.2240故非线性回归模型为故非线性回归模型为-0.224400.0904tye第第3章章 回归分析回归分析 预测为:预测为:ypred =11.1014即,照射即,照射16次后细菌残留数目为次后细菌残留数目为11
34、.1014,该预测符,该预测符合实际,显然比例合实际,显然比例3.1.3中多项式回归的结果合理。中多项式回归的结果合理。ci =355.2481 444.9326 -0.2561 -0.1919即参数即参数a置信度为置信度为95%的置信区间的置信区间 (ci的第一行的第一行)为:为: 355.2481 , 444.9326参数参数b的置信度为的置信度为95%的置信区间的置信区间 (ci的第二行的第二行)为为 -0.2561 -0.1919显然,最佳参数显然,最佳参数a=400.0904,b=-0.2240,均属于各,均属于各自自置信度为置信度为95%的置信区间。的置信区间。第第3章章 回归分析
35、回归分析 图图3.13原始数据与非线性回归图形原始数据与非线性回归图形图图3.14 原始数据与非线性回归原始数据与非线性回归GUI图形图形第第3章章 回归分析回归分析 从交互图形从交互图形3.14可以看出:圆圈为原始数据,两条虚线可以看出:圆圈为原始数据,两条虚线(屏幕上显示红色)是置信区间曲线;两条虚线内的实(屏幕上显示红色)是置信区间曲线;两条虚线内的实线(屏幕上显示绿色)是回归模型曲线;纵向虚线指示线(屏幕上显示绿色)是回归模型曲线;纵向虚线指示照射照射8次,此时对应的水平虚线表示模型得到的残留细次,此时对应的水平虚线表示模型得到的残留细菌数为:菌数为:66.6451。图图3.14 原始
36、数据与非线性回归原始数据与非线性回归GUI图形图形第第3章章 回归分析回归分析 3.1.3一元回归建模实例一元回归建模实例例例3.1.7在四川白鹅的生产性能研究中,得到如下一组在四川白鹅的生产性能研究中,得到如下一组关于雏鹅重(关于雏鹅重(g g)与)与7070日龄重日龄重(g)(g)的数据,试建立的数据,试建立7070日日龄重龄重( (y) )与雏鹅重与雏鹅重( (x) )的直线回归方程,计算模型误差的直线回归方程,计算模型误差平方和以及可决系数,当雏鹅重分别为:平方和以及可决系数,当雏鹅重分别为:85,9585,95 ,115,115时预测其时预测其7070日龄重,以及置信区间。日龄重,以
37、及置信区间。 表表3.5 5 四川白鹅重与四川白鹅重与7070日龄重测定结果日龄重测定结果 (单位:(单位:g g)编号编号123456789101112雏鹅重雏鹅重(x)80869890120102958311310511010070日龄重日龄重(Y)235024002720250031502680263024003080292029602860第第3章章 回归分析回归分析 解:(解:(1)作散点图。以雏鹅重()作散点图。以雏鹅重(x)为横坐标,)为横坐标,70日龄重(日龄重(y)为)为纵坐标作散点图,如图纵坐标作散点图,如图2-14。在在MATLAB命令窗口中输入命令窗口中输入:x=80
38、86 98 90 120 102 95 83 113 105 110 100; % 雏鹅重雏鹅重y=2350 2400 2720 2500 3150 2680 2630 2400 3080 2920 2960 2860; %70日龄重日龄重plot(x,y,*) %作散点图作散点图xlabel(x(雏鹅重雏鹅重) %横坐标名横坐标名ylabel(y(70日龄重日龄重) %纵坐标名纵坐标名图图3.15 四川白鹅的雏鹅重与四川白鹅的雏鹅重与70日龄重散点图和回归直线图日龄重散点图和回归直线图第第3章章 回归分析回归分析 由图形由图形3.15可见白鹅的可见白鹅的70日龄重与雏鹅重间存在直线日龄重与雏
39、鹅重间存在直线关系,且关系,且70日龄重随雏鹅重的增大而增大。因此,日龄重随雏鹅重的增大而增大。因此,可认为可认为y与与x符合一元线性回归模型。符合一元线性回归模型。(2)建立直线回归方程。在建立直线回归方程。在MATLAB中调用命令中调用命令polyfit,从而求出参数,从而求出参数 0, 1的最小二乘估计的最小二乘估计. 在在MATLAB命令窗口中继续输入命令窗口中继续输入:n= size(x,1) % 计算样本容量计算样本容量p,s=polyfit(x,y,1); % 调用命令调用命令polyfit计算回归参数计算回归参数y1=polyval(p,x); % 计算回归模型的函数值计算回归
40、模型的函数值hold onplot(x,y1) % 作回归方程的图形,结果如图作回归方程的图形,结果如图3.15p % 显示参数的最小二乘估计结果显示参数的最小二乘估计结果p=582.1850 21.7122第第3章章 回归分析回归分析 即参数即参数 的最小二乘估计为的最小二乘估计为),(107122.21,1850.58210所以所以70日龄重日龄重(y)与雏鹅重与雏鹅重(x)的直线回归经验方程为的直线回归经验方程为xy7122.211850.582(3)误差估计与决定系数。在)误差估计与决定系数。在MATLAB命令窗口中继续输入命令窗口中继续输入:TSS=sum(y-mean(y).2)
41、%计算总离差平方和计算总离差平方和RSS=sum(y1-mean(y).2) %计算回归平方和计算回归平方和ESS=sum(y-y1).2) %计算残差平方和计算残差平方和R2=RSS/TSS; %计算样本决定系数计算样本决定系数R2.第第3章章 回归分析回归分析 输出输出:TSS =8.314917e+005RSS =7.943396e+005ESS =3.715217e+004R2= 0.9553TSS =8.314917e+005RSS =7.943396e+005ESS =3.715217e+004R2=0.9553由于样本决定系数由于样本决定系数R2=0.9553接近于接近于1,因此
42、模型的,因此模型的拟合的效果较好。拟合的效果较好。第第3章章 回归分析回归分析 (4)回归方程关系显著性的)回归方程关系显著性的F检验。在检验。在MATLAB命令窗口中继续输入命令窗口中继续输入:F=(n-2)*RSS/ESS %计算的计算的F统计量统计量F1=finv(0.95,1,n-2) %查查F统计量统计量0.05的分位数的分位数F2=finv(0.99,1,n-2) %查查F统计量统计量0.01的分位数的分位数输出结果输出结果:F=2.138e+002,F1 =4.9646,F2 =10.0442为了方便,将以上的计算结果列成表为了方便,将以上的计算结果列成表3.6。表表3.6 四川
43、白鹅四川白鹅70日龄重与雏鹅重回归关系方差分析表日龄重与雏鹅重回归关系方差分析表自由度(自由度(df) 平方和平方和(SS) 均方和(均方和(MS) F值值F0.05F0.01回归回归1794339.60794339.60213.81*4.9610.04残差残差1037152.073715.21总离差总离差 11831491.67第第3章章 回归分析回归分析 因为因为 表明四川白鹅表明四川白鹅70日龄日龄重与雏鹅重间存在显著的线性关系。重与雏鹅重间存在显著的线性关系。(5)回归关系显著性的)回归关系显著性的t检验。在检验。在MATLAB命令窗口中命令窗口中继续输入继续输入:T=p(2)/sqr
44、t(ESS/(n-2)*sqrt(sum(x-mean(x).2)%计算计算T统计量统计量T1=tinv(0.975,n-2) %t统计量统计量0.05的分位数的分位数T2=tinv(0.995,n-2) %t统计量统计量0.01的分位数的分位数,04.1010, 181.21301. 02)(FFF输出输出:T =14.622, T1 =2.228, T2 =3.169因为因为T=14.62t0.01(10),否定,否定H0,接受,接受H1即四川白鹅即四川白鹅70日龄重日龄重(y)与雏鹅重与雏鹅重(x)的线性回归系数是显著的的线性回归系数是显著的,可用所可用所建立的回归方程进行预测和控制。建
45、立的回归方程进行预测和控制。第第3章章 回归分析回归分析 (6)预测预测x1=85,95,115; % 输入自变量输入自变量yc=polyval(p,x1) % 计算预测值计算预测值Y,Delta=polyconf(p,x1,s); I1=Y-Delta,Y+Delta % 置信区间置信区间输出输出:yc = 2427.72 2644.84 3079.08I1 = 2279.47 2575.96 2503.01 2786.67 2927.55 3230.62所以当雏鹅重分别为所以当雏鹅重分别为85,95,115时,时,白鹅白鹅70日龄重分别为日龄重分别为2427.72, 2644.84, 30
46、79.08;且且95%的的置信区间分别为:置信区间分别为:2279.47 ,2575.96,2503.01,2786.67, 2927.55,3230.62.第第3章章 回归分析回归分析 在程序中加入:在程序中加入:polytool(x,y) % 交互功能交互功能bar(x,y-y1), % 残差图残差图legend(残差残差)h=lillietest(y-y1) % 残差正态性检验残差正态性检验输出输出h = 0得到交互图形如图得到交互图形如图3.16所示,可以看出当雏鹅重为所示,可以看出当雏鹅重为100时,模时,模型给出型给出70日龄鹅重为日龄鹅重为2753.4016.图图3.16 四川白
47、鹅四川白鹅70日龄重与雏鹅重线性模型交互图日龄重与雏鹅重线性模型交互图第第3章章 回归分析回归分析 3.2多元线性回归模型多元线性回归模型3.2.1多元线性回归模型及其表示多元线性回归模型及其表示第第3章章 回归分析回归分析 );,(21YXXXp);, 2 , 1)(;,(21pnniyxxxiipii对于总体对于总体的的n组观测值组观测值它应满足式(它应满足式(3.2.1),即),即nnppnnnppppxxxyxxxyxxxy2211022222211021112211101其中其中 i (i=1,2,n)相互独立,且设相互独立,且设 记记), 2 , 1)(, 0(2niNinyyyY
48、21npnnppxxxxxxxxxX212222111211111p10n21, , , 第第3章章 回归分析回归分析 则模型则模型(3.2.2)可用矩阵形式表示为可用矩阵形式表示为 Y=X + (3.2.3)其中其中Y称为观测向量,称为观测向量,X称为设计矩阵,称为设计矩阵, 称为待称为待估计向量,估计向量, 是不可观测的是不可观测的n维随机向量,它的分维随机向量,它的分量相互独立,假定量相互独立,假定 .), 0(2nIN2. 多元线性回归建模的基本步骤多元线性回归建模的基本步骤(1)对问题进行直观分析,选择因变量与解释变量,作对问题进行直观分析,选择因变量与解释变量,作出与因变量与各解释
49、变量的散点图,初步设定多元线出与因变量与各解释变量的散点图,初步设定多元线性回归模型的参数个数;性回归模型的参数个数;(2)输入因变量与自变量的观测数据输入因变量与自变量的观测数据(y,X)调用命令调用命令 b,bint,r,rint,s=regress(y,X,alpha),计算参数的估计。计算参数的估计。(3)调用命令调用命令rcoplot(r,rint),分析数据的异常点情况。,分析数据的异常点情况。第第3章章 回归分析回归分析 (4)作显著性检验,若检验通过,则用模型作预测。作显著性检验,若检验通过,则用模型作预测。(5)对模型进一步研究:如残差的正态性检验,残差的对模型进一步研究:如
50、残差的正态性检验,残差的异方差检验,残差进行自相关性的检验等。异方差检验,残差进行自相关性的检验等。3.2.2 MATLAB的回归分析命令的回归分析命令在在MATLAB7.0的统计工具箱中,与多元回归模型的统计工具箱中,与多元回归模型有关的命令有多个,下面逐一介绍。有关的命令有多个,下面逐一介绍。1.多元回归建模命令多元回归建模命令regeress,其调用格式有以下三,其调用格式有以下三种:种:(1)b = regress(y,X)(2)b,bint,r,rint,stats = regress(Y,X)(3)b,bint,r,rint,stats = regress(Y,X,alpha) 第
51、第3章章 回归分析回归分析 三种方式的主要区别在输出项参数多少上,第三种方式的主要区别在输出项参数多少上,第3种方种方式可称为全参数方式。以第式可称为全参数方式。以第3种为例来说明种为例来说明regeress命令的输入与输出参数的含义。命令的输入与输出参数的含义。输入参数:输入量输入参数:输入量Y表示模型(表示模型(3.1.1)中因变量的)中因变量的观测向量;观测向量;X是一个的矩阵,其中第一列元全部是是一个的矩阵,其中第一列元全部是数数“1”,第,第j列是自变量列是自变量Xj的观测向量,即的观测向量,即,111,21222211121121npnnppnxxxxxxxxxXyyyY对一元线性
52、回归,取对一元线性回归,取p=1即可;即可;alpha为显著性水平为显著性水平输出参数:输出向量输出参数:输出向量b为回归系数估计值,为回归系数估计值,bint为回归为回归系数的系数的(1-alpha)置信区间;输出向量置信区间;输出向量r表示残差列向量表示残差列向量第第3章章 回归分析回归分析 输出输出rint为模型的残差的为模型的残差的 (1- )的置信区间;输出的置信区间;输出stats是用于检验回归模型的统计量,有是用于检验回归模型的统计量,有4个分量值:第一个个分量值:第一个是是R2,其中,其中R是相关系数,第二个是是相关系数,第二个是F统计量值,第三统计量值,第三个是与统计量个是与
53、统计量F对应的概率对应的概率P,当,当P 时拒绝时拒绝H0,即认,即认为线性回归模型有意义,第四个是方差为线性回归模型有意义,第四个是方差 2的无偏估计的无偏估计.例例3.2.1某销售公司将库存占用资金情况、广告投入的费某销售公司将库存占用资金情况、广告投入的费用、员工薪酬以及销售额等方面的数据作了汇总用、员工薪酬以及销售额等方面的数据作了汇总,该公司该公司试图根据这些数据找到销售额与其他变量之间的关系,试图根据这些数据找到销售额与其他变量之间的关系,以便进行销售额预测并为工作决策提供参考依据。以便进行销售额预测并为工作决策提供参考依据。(1)建建立销售额的回归模型;立销售额的回归模型;(2)
54、如果未来某月库存资金额为如果未来某月库存资金额为150万元,广告投入预算为万元,广告投入预算为45万元,员工薪酬总额为万元,员工薪酬总额为27万万元,试根据建立的回归模型预测该月的销售额。元,试根据建立的回归模型预测该月的销售额。 第第3章章 回归分析回归分析 表表3.7 占用资金、广告投入、员工薪酬、销售额(单位:万元)占用资金、广告投入、员工薪酬、销售额(单位:万元)月份月份库存资金额库存资金额(x1)广告投入广告投入(x2)员工薪酬总额员工薪酬总额(x3)销售额销售额(y)175.230.621.11090.4277.631.321.41133380.733.922.91242.1476
55、29.621.41003.2579.532.521.51283.2681.827.921.71012.2798.324.821.51098.8867.723.621826.397433.922.41003.31015127.724.71554.61190.845.523.2119912102.342.624.31483.113115.64023.11407.11412545.829.11551.315137.851.724.61601.216175.667.227.52311.717155.26526.52126.718174.365.426.82256.5第第3章章 回归分析回归分析 解:为
56、了确定销售额与库存占用资金、广告投入、员解:为了确定销售额与库存占用资金、广告投入、员工薪酬之间的关系,分别作出工薪酬之间的关系,分别作出y与与x1,x2,x3的散点图,的散点图,若散点图显示它们之间近似线性关系,则可设定若散点图显示它们之间近似线性关系,则可设定y与与x1,x2,x3的关系为三元线性回归模型的关系为三元线性回归模型01 12233yxxx%输入数据并作散点图(图输入数据并作散点图(图3.18)A=75.2 30.6 21.1 1090.4;77.6 31.3 21.4 113380.7 33.9 22.9 1242.1;76 29.6 21.4 1003.279.5 32.5
57、 21.5 1283.2;81.8 27.9 21.7 1012.298.3 24.8 21.5 1098.8;67.7 23.6 21 826.374 33.9 22.4 1003.3;151 27.7 24.7 1554.690.8 45.5 23.2 1199;102.3 42.6 24.3 1483.1115.6 40 23.1 1407.1;125 45.8 29.1 1551.3137.8 51.7 24.6 1601.2;175.6 67.2 27.5 2311.7155.2 65 26.5 2126.7;174.3 65.4 26.8 2256.5;第第3章章 回归分析回归分析
58、 m,n=size(A);subplot(3,1,1),plot(A(:,1),A(:,4),+),xlabel(x1(库存资金额库存资金额) ylabel(y(销售额销售额)subplot(3,1,2),plot(A(:,2),A(:,4),*),xlabel(x2(广告投入广告投入) ylabel(y(销售额销售额)subplot(3,1,3),plot(A(:,3),A(:,4),x),xlabel(x3(员工薪酬员工薪酬) ylabel(y(销售额销售额)所得图形如图所得图形如图3.18所示,可见销售额所示,可见销售额y与库存资与库存资金、广告投入、员工薪酬具有线性关系,因此可金、广告
59、投入、员工薪酬具有线性关系,因此可以建立三元线性回归模型以建立三元线性回归模型.第第3章章 回归分析回归分析 图图3.18销售额与库存、广告、薪酬散点图销售额与库存、广告、薪酬散点图% 调用命令调用命令regress建立三元线性回归模型建立三元线性回归模型x=ones(m,1), A(:,1), A(:,2), A(:,3);y=A(:,4)b,bint,r,rint,stats=regress(y,x);b,bint,stats, % 输出结果输出结果 第第3章章 回归分析回归分析 程序运行结果程序运行结果b =162.0632 7.2739 13.9575 -4.3996bint =-58
60、0.3603 904.4867 4.3734 10.1743 7.1649 20.7501 -46.7796 37.9805 stats =0.9574804050 105.0866520891 0.0000000008 10077.9867891125输出结果说明:输出结果说明:b就是模型中的参数就是模型中的参数 0 , 1 , 2 ,因此回归模型为,因此回归模型为123162.06327.273913.95794.3996yxxxb就是模型中的参数就是模型中的参数 0 , 1 , 2 ,因此回归模型为,因此回归模型为123162.06327.273913.95794.3996yxxxbin
61、t的各行分别为参数的各行分别为参数 0 , 1 , 2的的95%的置信区间。的置信区间。stats的第一列为模型可决系数,第二列为的第一列为模型可决系数,第二列为F统计量的观统计量的观测值,第三列得到概率测值,第三列得到概率p,最后一列为模型的残差平方和,最后一列为模型的残差平方和第第3章章 回归分析回归分析 2.多元回归辅助图形命令多元回归辅助图形命令(1)残差图命令rcoplot,其调用格式rcoplot(r,rint)其中,输入参数r,rint 是多元回归建模命令regress 输出的结果,运行该命令后展示了残差与置信区间的图形。该命令有助于对建立的模型进行分析,如果图形中出现红色的点,
62、则可以认作异常点,此时可删除异常点,重新建模,最终得到改进的回归模型。在上面的程序中加入在上面的程序中加入 rcoplot(r,rint)得到如下图形得到如下图形第第3章章 回归分析回归分析 图图3.19 残差与置信区间图残差与置信区间图 从图形中可以看到第五个点为异常点,实际上从从图形中可以看到第五个点为异常点,实际上从表表3.7可以发现第可以发现第5个月库存占用资金、广告投入、员个月库存占用资金、广告投入、员工薪酬均比工薪酬均比3月份少,为何销售额反而增加?这就可月份少,为何销售额反而增加?这就可以促使该公司的经理找出原因,寻找对策。下面的例以促使该公司的经理找出原因,寻找对策。下面的例题
63、介绍如何删除异常点,对模型进行改进的方法。题介绍如何删除异常点,对模型进行改进的方法。第第3章章 回归分析回归分析 例例3.2.2 葛洲坝机组发电耗水率的主要影响因素为库葛洲坝机组发电耗水率的主要影响因素为库水位水位,出库流量。数据如表出库流量。数据如表3.8所示,利用多元线性回所示,利用多元线性回归分析方法建立耗水率与出库流量、库水位的模型。归分析方法建立耗水率与出库流量、库水位的模型。表表3.8 某天耗水率与出库流量、库水位的数据某天耗水率与出库流量、库水位的数据 时间年时间年-月月-天天-时时 库水位(米库水位(米)出库流量(立方米)出库流量(立方米) 机组发电耗水率机组发电耗水率 (
64、(立方米立方米/ /万千瓦万千瓦) ) 2005-10-15:0065.081560760.462005-10-15:0265.101556560.282005-10-15:0465.121554060.102005-10-15:0665.171550759.782005-10-15:0865.211543259.442005-10-15:1065.371561959.252005-10-15:1265.381553658.912005-10-15:1465.391551458.762005-10-15:1665.401551958.732005-10-15:1865.431551058.63
65、2005-10-15:2065.471548958.482005-10-15:2265.531543758.312005-10-16:0065.621635557.962005-10-16:0265.581470857.062005-10-16:0465.701439356.432005-10-16:0665.841429655.83第第3章章 回归分析回归分析 解:解:% 输入原始数据输入原始数据A=65.08 15607 60.4665.10 15565 60.2865.12 15540 60.1065.17 15507 59.7865.21 15432 59.4465.37 15619
66、59.2565.38 15536 58.9165.39 15514 58.7665.40 15519 58.7365.43 15510 58.6365.47 15489 58.4865.53 15437 58.3165.62 16355 57.9665.58 14708 57.0665.70 14393 56.4365.84 14296 55.83;% 做散点图做散点图subplot(1,2,1),plot(A(:,1),A(:,3),+)xlabel(x1(库水位库水位) ylabel(y(耗水率耗水率)subplot(1,2,2),plot(A(:,2),A(:,3),o)xlabel(x2(出库流量出库流量) ylabel(y(耗水率耗水率)运行后得到的图形如图运行后得到的图形如图3.20所示,所示,从图中可以看到无论是库水位还从图中可以看到无论是库水位还是出库流量都与机组发电耗水率是出库流量都与机组发电耗水率具有线性关系,因此,可以建立具有线性关系,因此,可以建立机组发电耗水率与库水位和出库机组发电耗水率与库水位和出库流量的二元线性回归模型。流量的二元线性回归模型。第第3章章