实验8_多元线性回归分析与非线性回归分析
SAS 软件与统计应用实验徐雅静、汪远征编著实验 8多元线性回归分析与非线性回归分析多元线性回归分析研究多个变量的数量伴随关系,内容主要包括模型的假定与检验、参 数的估计与检验、回归诊断与预测。 很多非线性回归问题都可以转化为线性回归问题处理,如多项式回归、指数回归、对数 回归、幂函数回归等。8.1实验目的掌握使用 SAS 多元线性回归分析与非线性回归分析的方法。8.2实验内容一、用“分析家”作多元线性回归分析 二、用 INSIGHT 模块作多项式回归 三、使用 REG 过程作回归分析 四、一元非线性回归分析8.3实验指导一、用“分析家”作多元线性回归分析【实验 8-1】某研究人员需要分析我国固定资产投资状况的影响因素,选取 5 个可能的影 响因素:国内生产总值、商品房屋销售额、财政支出、社会消费品零售总额、进出口总额,统 计 1987~2001 共 15 年的各项指标如表 8-1 所示(sy8_1.xls)所示。试在 0.05 的显著性水平下 进行多元回归分析,判断哪些因素对固定资产投资有着显著影响,给出回归方程。表 8-1年度 1987 1988 1989 1990 1991 1992 1993 固定投资总额 3791.7 4753.8 4410.4 4517 5594.5 8080.1 13072.3 国内生产总值 11962.5 14928.3 16909.2 18547.9 21617.8 26638.1 34634.415 年的统计数据财政支出 2262.18 2491.21 2823.78 3083.59 3386.62 3742.2 4642.3 社会消费品零售总额 5820 7440 8101.4 8300.1 9415.6 10993.7 12462.1 进出口总额 3084.2 3821.8 4155.9 5560.1 7225.8 9119.6 11271商品房屋销售额 1100967 1472164 1637542 2018263 2378597 4265938 86371411SAS 软件与统计应用实验1994 1995 1996 1997 1998 1999 2000 2001 17042.1 20019.26 22913.55 24941.11 28406.17 29854.71 32917.73 37213.49 46759.4 58478.1 67884.6 74462.6 78345.2 82067.46 89442.2 95933.3 10184950 12577269 14271292 17994763 25133027 29878734 39354423 48627517 5792.62 6823.72 7937.55 9233.56 10798.18 13187.67 15886.5 18902.58 16264.7 20620 24774.1 27298.9 29152.5 31134.7 34152.6 37595.2徐雅静、汪远征编著20381.9 23499.9 24133.8 26967.2 26857.7 29896.3 39274.2 42193.31. 生成数据集 在“分析家”中直接打开上面的 Excel 数据表(sy8_1.xls),选择编辑状态,修改每个变量 的属性,将变量名分别改为:年度:n、固定投资总额:y、国内生产总值:x1、商品房屋销 售额:x2、财政支出:x3、社会消费品零售总额:x4、进出口总额:x5。 以数据集 Mylib.sy8_1 存盘。 2. 全回归分析 1) 选 择 主 菜 单 “ Statistics ( 统 计 ) ”→ “Regression(回归) ”→“Linear(线性) ” ,打开 “Linear Regression(线性回归) ”对话框。 2) 选择变量列表中的变量 y, 单击 “Dependent” 按钮,选定响应变量,选择变量列表中的变量 x1、 x2、x3、x4、x5,单击“Explanatory”按钮,选定 图8-1 Linear Regression对话框 解释变量,如图 8-1 所示。 3) 单击“OK”按钮,得到分析结果如图 8-2 所示。图8-2多元回归分析结果分析结果包括方差分析表、拟合的汇总信息以及回归系数估计值与显著性检验。方差分 析表中显示模型的作用是显著的(F 统计量的值为 1567.35,p 值SAS 软件与统计应用实验徐雅静、汪远征编著不都是显著的。例如 x3、x4、x5 的作用就不显著。因此有必要适当选择变量建立一个“最优” 的回归方程。 3. 逐步回归分析 1) 重复上面 2 中 1),在“Linear Regression(线 性回归) ”对话框(图 8-1)中,单击“Model”按钮, 打开 “Linear Regression: Model” 对话框。 在 “Method” 选项卡中选择“Stepwise selection(逐步选择法) ” ,如 图 8-3 所示。 两次单击“OK”按钮,得到分析结果。 2) 在显示结果中, 第 1 步记录了只有x1 进入回归 图8-3 选择逐步回归法 方程的回归分析结果,其中回归方程和系数的检验均 2 为显著,此时R =0.9911,C(p)=58.5161;接着第 2 步 是自变量x1 和x2 进入回归方程后的回归分析结果, 回归方程及x1 和x2 的系数检验均为显著, 但常数项检验不显著。接着第 3 步是自变量x1、x2 和x3 进入回归方程后的回归分析结果。其 C(p)=5.5226 中回归方程及所有系数检验均为显著, 常数项检验也显著。 且R2=0.9984 提高了, 减少了。图8-4逐步回归第1、2步、3步及最后结果在图 8-4 右下中指出在 0.05 的检验水平下,不能再有其它变量进入模型。比较R2和C(p) 的值(图 8-4 右) ,应取包含变量x1、x2和x3的第三个模型作为较优的模型,对应的回归方程 是:y = 3023.27814 + 0.36911x1 + 0.00078157 x 2 − 2.09048 x33SAS 软件与统计应用实验徐雅静、汪远征编著4. 回归诊断 进行回归诊断的步骤如下: 1) 重复上面 2 中 1),在打开的“Linear Regression (线性回归) ” 对话框中, 单击 “Plots” 按钮。在打开的“Linear Regression:Plots”对 话框中,选择“ Residual ”选项卡,按图 8-5 所示选择有关复选框。 2) 两次单击“OK”按钮,得到回归诊断 结果,在“分析家”窗口的项目管理器中依次 图8-5 Linear Regression:Plots对话框 双击 “Residual Plots” 下的 “Plot of STUDET vs PRED”和“Plot of RESIDUAL vs NQQ”得到 标准化后的残差图(图 8-6 左)和残差的 QQ 图(图 8-6 右) 。图8-6残差图和残差的QQ图从标准化后的残差图(图 8-6 左)看出,数据点随机地散布在零线附近,表明模型中误 差等方差、独立性的假设没有问题。残差的 QQ 图(图 8-6 右)近似一条直线,可以初步判 定残差来自正态分布总体,所建回归模型是有效的。 3) 对残差作进一步检验: 在上述操作打开的 “Linear Regression (线 性回归) ”对话框中,单击“Save Data”按钮。 在打开的“Linear Regression:Save Data”对 话框中, 选中 “Crate and save diagnostics data” 复选框,并将列表中的第二项“ RESIDUAL Residuals” 添加到左边方框内, 如图 8-7 所示。 两次单击“OK”后得到分析结果。 4) 在“分析家”窗口的项目管理器中双 击“Diagnostics”下的“Diagnostics Table” 图8-7 Linear Regression:Save Data对话框 可以看到在数据集中生成了残差数据,如图 8-8 所示。4SAS 软件与统计应用实验徐雅静、汪远征编著图8-8生成残差数据将“Diagnostics Table”存盘(sy8_1_r)后在“分析家”中打开。 5) 选择主菜单“Statistics(统计) ”→“Descriptive(描述性统计) ”→“Distributions… (分布) ” ,打开“Distributions”对话框,选择变量列表中的_RESID,单击“Analysis”按钮, 选定分析变量,如图 8-9 左所示。图8-9设置选项6) 单击“Fit(拟合) ”按钮,在打开的对话框中选择拟合的分布类型:Normal,使用样 本估计量(Sample estimates) ,如图 8-9 右所示。 7) 两次单击“OK”按钮,并在分析家窗口的项目管理器中双击“Fitted Distributions of sy8_1_r”项,得到对残差_RESID 的正态分布检验结果,如图 8-10 所示。图8-10残差分布检验结果5SAS 软件与统计应用实验徐雅静、汪远征编著三种检验均有p值>0.05,因此不能拒绝残差来自正态总体的假定。 5. 预测 通过回归诊断得知模型:y = 3023.27814 + 0.36911x1 + 0.00078157 x 2 − 2.09048 x3是合适的,可以用于预测。 1) 假定 02,03 年国内生产总值(x1) 、商品房屋销售额(x2) 、财政支出(x3)的数据 已存入数据集 Mylib.sy8_1_new 中,如图所示。图8-11数据集Mylib.sy8_1_new2) 重复上面逐步回归步骤,并在图 8-1 所示的“Linear Regression(线性回归) ”对话框 中,单击“predictions”按钮,打开“Linear Regression:predictions”对话框。按图 8-12 所 示进行预测的 Input(输入) 、Output(输出)设置。图8-12“Linear Regression:predictions”对话框3) 两次单击“OK” ,得到结果。在分析家的项目管理器中点击“predictions”可以看到 预测结果,如图 8-13 所示。图8-13预测结果二、用 INSIGHT 模块作多项式回归【实验 8-2】某研究人员统计了房地产行业 2003 年主营业务收入和净利润的关系,从中 6SAS 软件与统计应用实验徐雅静、汪远征编著随即抽取 20 家上市公司,统计后的主营业务收入和净利润如表 8-2(sy8_2.xls)所示。试采用 多项式回归方程给出净利润y与主营业务收入x的关系。表 8-2股票代码 0000029 0000031 0000046 0000511 0000558 0000608 0600064 0600322 0600376 主营业务收入 96.24 33.09 57.60 36.77 20.10 45.53 110.07 72.05 39.8220 家公司的净利润和主营业务收入(百万元)净利润 11.35 85.11 36.61 27.26 6.72 36.37 80.13 20.88 37.15 股票代码 600533 600603 600638 600639 600641 600648 600663 600675 600736 主营业务收入 48.05 0.85 51.66 61.07 64.54 35.50 226.55 188.60 81.32 净利润 70.50 2.57 63.70 127.13 33.59 25.01 664.02 175.31 108.681. 生成数据集 将表 8-2 在 Excel 中修改后导入成 SAS 数据集 Mylib.sy8_2,在 INSIGHT 中打开后如图 8-14 所示。其中变量 n、x、y 分别表示股票代码、主营业收入和净利润。 2. 回归分析 为了大致地了解y与x的关系,首先利用INSIGHT作y与x的带有曲线拟合的散点图,具体 方法是: 1) 选择菜单“Analyze(分析)”→“Fit(Y X) (拟合)”打开“Fit(Y X)”对话框,选择变量 如图 8-15 所示。图8-14数据集Mylib.sy8_2图8-15选择变量单击“OK” ,得到拟合线性模型:y = β 0 + β1 x + ε主要部分如图 8-16 所示。 7SAS 软件与统计应用实验徐雅静、汪远征编著图8-16一元线性模型的拟合结果从图 8-16 可以发现,虽然模型检验显著,但R2只有 0.6494,且常数项未通过检验。为改 进模型,移动图 8-16 右上第一张表中“Degree(Polynomial)”栏下的滚动条,可以做高次多 项式拟合试验,其中R-Square、F Stat可以说明拟合的效果。 作二次和三次多项式拟合试验的结果如图 8-17 所示。图8-172调整拟合模型拟合三次多项式的R 已接近 94%,不需要再高的阶次了。 重新回到“Fit(Y X)”对话框,按图 8-18 左设置分析变量,可以拟合下面的三次多项式 模型:y = β 0 + β1 x + β 2 x 2 + β 3 x 3 + ε8SAS 软件与统计应用实验徐雅静、汪远征编著图8-18设置拟合变量单击“Output”按钮,按图 8-18 右选择输出结果,两次单击“OK”后,得到拟合结果, 如图 8-19 所示。图8-19三次多项式拟合结果虽然模型检验通过,但常数项未通过检验,进一步改进模型,再次回到“Fit(Y X)”对话 框,按图 8-18 进行有关设置后,取消常数项复选框,如图 8-20 所示。图8-20取消常数项复选框单击“OK”按钮,得到拟合模型 y = 所示。β1 x + β 2 x 2 + β 3 x 3 + ε ,主要部分如图 8-21、8-229SAS 软件与统计应用实验徐雅静、汪远征编著图8-21不含常数项的三次多项式拟合结果其中模型和系数检验均通过,且R2高达 84.98%,所以该模型为比较理想的模型,最后的 回归方程具体形式为:y = 3.3145 x − 0.0527 x 2 + 0.0002 x 33. 回归诊断 从图 8-22 可以看出残差数据随机独立并大致服从均值为零的正态分布。图8-22残差图和残差QQ图三、使用 REG 过程作回归分析【实验 8-3】某种水泥在凝固时放出的热量y(cal/g)与水泥中四种化学成分x1,x2,x3,x4 有关,现测得 13 组数据,如表 8-3(sy8_3.xls)所示。试从中选出主要的变量,建立y关于它 们的线性回归方程。表 8-3x1 7 1 11 11 7 11 x2 26 29 56 31 52 55热量 y 与四种化学成分的实测数据x3 6 15 8 8 6 9 x4 60 52 20 47 33 22 y 78.5 74.3 104.3 87.6 95.9 109.210SAS 软件与统计应用实验3 1 2 21 1 11 10 71 31 54 47 40 66 68 17 22 18 4 23 9 8 6 44 22 26 34 12 12徐雅静、汪远征编著102.7 72.5 93.1 115.9 83.8 113.3 109.41. 建立数据集 输入以下代码建立数据集 sy8_3 并显示:data mylib.sy8_3; input x1 x2 cards; 7 1 11 11 7 11 3 1 2 21 1 11 10 ; Title '数据集 sy8_3'; Proc print ; run; 26 29 56 31 52 55 71 31 54 47 40 66 68 6 15 8 8 6 9 17 22 18 4 23 9 8 60 52 20 47 33 22 6 44 22 26 34 12 12 78.5 74.3 104.3 87.6 95.9 109.2 102.7 72.5 93.1 115.9 83.8 113.3 109.4 x3 x4 y;运行结果如图所示。 2. 向后逐步剔出法进行回归 执行以下代码:proc var reg y data = Mylib.sy8_3; x1 - x4;model y = x1 - x4/selection=backward; plot residual. * predicted.; run;输出结果如下:图8-23 数据集sy8_311SAS 软件与统计应用实验徐雅静、汪远征编著图8-24向后逐步剔除的第0步(全回归)图8-25向后逐步剔除的第1步图8-26向后逐步剔除法第2步12SAS 软件与统计应用实验徐雅静、汪远征编著图8-27向后逐步剔除法结果汇总向后逐步剔除法的分析结果给出回归模型: Y = 52.57735 + 1.46831x1 + 0.66225x2 残差对预测值的散点图显示如下:图8-28残差散点图3. 结果分析 采用向后逐步剔除法回归的第 0 步是做全回归,结果如图 8-24 所示,所有系数均未通过 检验(P值均大于 0.05) ,向后逐步剔除法第 1 步将变量x3 剔除,结果如图 8-25 所示,其中 x2 和x4 的系数仍不能通过检验,接下来第 2 步将变量x4 剔除,结果如图 8-26 所示,此时的 回归方程及x1 和x2 的系数均能通过检验,残差对预测值的散点图(图 8-28)基本正常符合 模型假定,所以方程Y = 52.57735 + 1.46831x1 + 0.66225x2为有效回归方程。四、一元非线性回归分析【实验 8-4】已知数据如表 8-4(sy8_4.xls)所示。试分别采用指数回归、对数回归、幂函 数回归和倒幂函数回归 4 种非线性回归方法进行回归分析,并选择一个较好的回归方程。 13SAS 软件与统计应用实验 表 8-4X Y 1.1 109.95 1.2 40.45 1.3 20.09 1.4 24.53 1.5 11.02 1.6 7.39徐雅静、汪远征编著实验数据1.7 4.95 1.8 2.72 1.9 1.82 2 1.49 2.1 0.82 2.2 0.3 2.3 0.2 2.4 0.221. 生成数据集 运行下面程序生成并显示数据集 sy8_4,如图 8-29 所示。data sy8_4; input x y; cards; 1.1 1.2 1.3 1.4 1.5 1.6 1.7 1.8 1.9 2 2.1 2.2 2.3 2.4 ; run; title '数据集 sy8_4'; proc print; run; 109.95 40.45 20.09 24.53 11.02 7.39 4.95 2.72 1.82 1.49 0.82 0.3 0.2 0.22 图8-29 数据集sy8_42. 对 x 和 y 作相关分析 执行如下代码:/*画 x 和 y 的散点图*/ goptions ftext='宋体'; proc gplot data = sy8_4; plot y*x; title 'x 和 y 的散点图'; symbol v=dot i=none cv=orange ; run; /*求 x 和 y 的相关系数*/ proc var run; x corr y; data = sy8_4;14SAS 软件与统计应用实验徐雅静、汪远征编著运行上面程序,得到散点图(图 8-30 左)以及 x 与 y 的相关系数(图 8-30 右) :图8-30x与y的散点图与相关系数由图可见x和y有一定的非线性关系,根据散点分布的形状考虑用下面几种非线性回归方 法建立非线性回归方程,并从中选出较为合适的回归方程。 3. 倒幂函数 y = a + b1 回归 x首先考虑倒幂函数拟合,执行如下代码:goptions ftext='宋体'; data new1; set sy8_4; u = 1/x; v = y; run; /*画 u 和 v 的散点图*/ title 'u 和 v 的散点图'; proc gplot data = new1; 图8-31 u和v的散点图 plot v*u; symbol v=dot i=none cv=red ; run;运行结果得到散点图 8-31,由图可见,u 和 v 有着较弱的线性关系。做线性回归:proc var reg v data = new1; u;model v = u; print cli; title '残差图'; plot residual. * predicted.; run;运行结果如图 8-32 和图 8-33 所示。15SAS 软件与统计应用实验徐雅静、汪远征编著图8-32倒幂函数回归结果倒幂函数回归结果(图 8-32) :方差分析表中显示模型的作用是显著的(F 统计量的值为 24.95,p 值y = −78.56560 + 156.538871 x残差对预测值的散点图(图 8-33)表明,残差有一定趋势,不符合模型的假定,以上回归 方程无效。图8-33残差对预测值的散点图4. 幂函数 y = ax b 回归 考虑幂函数拟合,执行如下代码:data new2; set sy8_4; u = log(x);16SAS 软件与统计应用实验 v = log(y); run; /*画 u 和 v 的散点图*/ title 'u 和 v 的散点图'; proc gplot data = new2; plot v*u; symbol v=dot i=none cv=red ; run; title '残差图'; proc var reg v data = new2; u; 图8-34徐雅静、汪远征编著model v = u; print cli; plot residual. * predicted.; run; u与v的散点图得到散点图如图 8-34 所示: 幂函数回归的结果见图 8-35 左:图8-35幂函数回归结果与残差对预测值的散点图得回归方程:v = 5.51053 – 7.93588u 即: y = 247.2822 x −7.93588 残差对预测值的散点图(如图 8-35 右)表明,残差有微弱趋势,不符合模型的假定,上面 回归方程不佳。 5. 指数函数 y = ae bx 回归 考虑指数函数拟合,执行如下代码:data new3; set sy8_4; u = x;17SAS 软件与统计应用实验 v=log(y); run; /*画 u 和 v 的散点图*/ title 'u 和 v 的散点图'; proc gplot data = new3; plot v*u; symbol v=dot i=none cv=red ; run; title '残差图'; proc var reg v data = new3; u; 图8-36徐雅静、汪远征编著model v = u; print cli; plot residual. * predicted.; run; u与v的散点图得到散点图如图 8-36 所示: 指数函数回归结果见图 8-37 左:图8-37指数函数回归结果与残差对预测值的散点图得回归方程:v = 9.58399 – 4.73895u 即: y = 14530.28e −4.73895 x 从残差对预测值的散点图(如图 8-37 右)可以看出,残差基本符合模型的假定,上面回归 方程有效。 6. 对数 y = a + b ln x 回归 考虑对数函数拟合,执行如下代码:data new4; set sy8_4; u = log(x);18SAS 软件与统计应用实验 v = y; run; /*画 u 和 v 的散点图*/ title 'u 和 v 的散点图'; proc gplot data = new4; plot v*u; symbol v=dot i=none cv=red ; run; title '残差图'; proc var reg v data = new4; u; 图8-38徐雅静、汪远征编著model v = u; print cli; plot residual. * predicted.; run; u与v的散点图得到散点图如图 8-38 所示。 对数函数回归结果见图 8-39 左。图8-39对数函数回归结果得回归方程:v = 64.58847 – 91.11730u 即:y = 64.58847 – 91.11730lnx 从残差对预测值的散点图(如图 8-39 右)可以看出, 残差有二次趋势, 不符合模型的假定, 上面回归方程无效。 7. 结论 比较上述 4 个回归方程,第三种指数函数回归的 Root MSE (均方残差平方根)最小 ,效果最好。 (0.25991) 、R-Square(判定系数R2)最大(0.9844) 执行下述代码,得到模型 y = 14530.28e −4.73895 x 的拟合图形如图 8-40 所示。data new5; set new1;19SAS 软件与统计应用实验 y1 = 14530.28*exp(-4.73895*x); run; title '回归图'; proc gplot data = new5; plot y*x=1 y1*x=2/overlay ; symbol v=dot i=none cv=red ; symbol2 i=sm color=blue; run;徐雅静、汪远征编著图8-40指数函数拟合图形8.4上机演练【练习 8-1】某学校 20 名一年级女大学生体重(公斤) 、胸围(厘米) 、肩宽(厘米)及 肺活量(升)实测值如表 8-5(lx8_1.xls)所示。试对影响女大学生肺活量的有关因素作多元 回归分析。表 8-5编号 1 2 3 4 5 6 7 8 9 10 11 1220 名一年级女大学生肺活量及有关变量测量结果胸围 X2/厘米 73.6 83.9 78.3 77.1 81.7 74.8 73.7 79.4 72.6 79.5 83.8 88.4 肩宽 X3/厘米 36.4 34.0 31.0 31.0 30.0 32.0 36.5 37.0 30.1 37.1 33.9 38.0 肺活量 Y/升 2.9 3.11 1.91 2.63 2.86 1.91 2.98 3.28 2.52 3.27 3.10 3.28体重 X1/公斤 51.3 48.9 42.8 55.0 45.3 45.3 51.4 53.8 49.0 53.9 48.8 52.620SAS 软件与统计应用实验13 14 15 16 17 18 19 20 42.7 52.5 55.1 45.2 51.4 48.7 51.3 45.2 78.2 88.3 77.2 81.6 78.3 72.5 78.2 74.7 30.9 38.1 31.1 30.2 36.5 30.0 36.4 32.1徐雅静、汪远征编著1.92 3.27 2.64 2.85 3.16 2.51 3.15 1.92【练习 8-2】在光刻工艺过程中,要求找出国产光致抗蚀剂显影的腐蚀速率与显影时间的 关系,实验中观测的数据经整理如表 8-6(lx8_2.xls)所示。表 8-6x y 5 94.40 10 54.90 15 30.80 20 26.10腐蚀速率与显影时间的关系25 20.28 30 19.96 35 14.99 40 14.30 45 14.70 50 13.90试分别采用指数回归、对数回归、幂函数回归和倒幂函数回归四种非线性回归方法分别 给出回归方程,并选择一个较为合适的回归方程。8.5实验报告请按练习内容写出包括如下内容的实验报告: 一、实验目的; 二、实验内容及程序与结果分析; 三、实验体会(问题、评价、感想与建议等)。21