SAS讲义第三十二课多元线性回归分析报告

word第三十二课 多元线性回归分析一、 多元回归模型表示法通常，回归模型包括k个变量，即一个因变量和k个自变量包括常数项。由于具有N个方程来概括回归模型：(32.1)模型的相应矩阵方程表示为：Error! Reference source not found.(32.2)式中;(32.3)其中,Y为因变量观察的N列向量，X为自变量观察的N (k+1) 矩阵，为末知参数的(k+1) )列向量，e 为误差观察的N列向量。在矩阵X表达式中，每一个元素Xij都有两个下标，第一个下标表示相应的列变量，第二个下标表示相应的行观察。矩阵X的每一列表示相应的给定变量的N次观察的向量，与截矩有关的所有观察值都等于1。经典的线性回归模型的假设可以阐述如下：l 模型形式由(32.1)给定；l 矩阵X的元素都是确定的，X的秩为(k+1)，且k小于观察数N；l e 为正态分布，E(e )=0 和 ，式中I为NN单位矩阵。根据X的秩为(k+1) 的假定，可以保证不会出现共线性。如果出现完全共线性，矩阵X的一列将为其余列的线性组合，而X的秩将小于(k+1) )，关于误差的假设是最有用的假设，因为用它可以保证最小二乘法估计过程的统计性质。除了正态性外，我们还假定每一个误差项的平均值为0，方差为常数， 以与协方差为 0 。假假设我们按Y的分布来表示第三个假设，如此可写成下式：(32.4)二、 最小二乘法估计我们的目的是求出一个参数向量使得残差平方和最小，即：(32.5)式中：(32.6)(32.7)其中，表示回归残差的N列向量，而表示Y拟合值的N列向量，表示为估计参数的(k+1) 列向量，将式(32.6)和式(32.7)代入式(32.5)，如此得：(32.8)为了确定最小二乘法估计量，我们求ESS对进展微分，并使之等于0，即：(32.9)所以：(32.10)被称为“交叉乘积矩阵，即Error! Reference source not found.矩阵能够保证逆变换，这是因为我们假设X的秩为(k+1),该假设直接导致了的非奇异性。最小化的二阶条件是， 是一个正定矩阵。最小二乘法残差有一个有益的特性，即：(32.11)这个结果说明自变量和残差的交叉乘积的总和为O，这个公式在一些推导中是非常有用的。 现在可以考虑最小二乘估计量的性质。首先可以证明它们是无偏估计量。因为：(32.12)设式中，且是常数，这样：(32.13)根据式(32.13) ，可以看到，只要遗漏变量都是随机分布的，与X无关，并且具有 0 均值，如此最小二乘法估计量将是无偏的。(32.14)我们看到，最小二乘法估计量为线性和无偏估计量。事实上，为的最优线性无偏估计量，也就是说，它在全部无偏估计量中方差最小，这就是著名的高斯马尔可夫定理。为了证明高斯马尔可夫定理，我们需要证明，任何其他线性估计量b的方差比的方差大。请注意=AY。为了不失去一般性，我们可写成：(32.15)假设b是无偏的，如此：(32.16)式(32.16)成立的一个必要和充分的条件是，这样就可以研究矩阵。由于，所以有：(32.17)由于：因为，所以，即：(32.18)我们可以看出，为一半正定矩阵。该矩阵的二次型为0，只有当0所有元素为0时才出现。当时，另外的估计量就是普通最小二乘法估计量，这样，我们的定理就得到证明。三、 的估计和t检验为了计算估计参数的方差协方差矩阵，我们需要给出的估计量，该估计量自然选为：(32.19)证明为 的一个无偏估计量，虽很单调冗长，但不困难。因此，是Var()的估计。当为时，可用正态分布假设检验。当用近似时，我们不得不用t假设检验。为此，我们利用以下的统计结果：l 假设 ，如此服从分布，具有Nk1个自由度；l Error! Reference source not found.服从分布，具有Nk1个自由度；l Error! Reference source not found.，当i=0，1，2,，k时，服从正态分布，平均值为0，方差为，其中vi为的第i个对角线元素；l Error! Reference source not found.和 相互独立。由此得出：(32.20)该式为t分布，具有(N-k-1)个自由度。这就使我们能按照与前面所述一样的方式确定各个回归参数的置信区间。假设t值的绝对值相当大，就可以在适当选定的置信水平上否认原假设，参数的置信区间可由下式得出：(32.21)其中，为与显著水平有关的t分布临界值。四、 R2和F检验我们可将Y的总变差分成两局部，一局部代表已说明变差，另一局部代表未说明变差。为了简化公式推导过程，首先我们假定Y变量具有0平均值，即 =0，如此有：(32.22)由于和，因此：(32.23)式中为总平方和，为回归已说明平方和，为残差未说明平方和，归纳成回归方差分析表，见表32.1。表32.1 回归方差分析表变异来源source离差平方和SS自由度df均方MSF统计量FP概率值P回归RP误差E总变异T从而(32.24)假设因变量不具有0平均值，我们必须改良一下的定义。这样，由此可以得出：(32.25)和(32.26)注意到一个数学上的事实：随着模型中增添新的变量，必定会增加，从而只要给模型增添越来越多的新因素，就可能使得人为地增大。在一元回归时已经指出较大常指模型与数据拟合得较好，在多元回归时很容易错误地去寻找一个极大化的回归模型。我们应该知道一个好的多元回归模型，应是具有合理个数的有意义自变量的简单模型。为了解决这个问题，提出了修正，使得只有当新增变量确实对因变量有所作用时修正才会增加。我们定义为修正的，它校正拟合优度对自由度的依赖关系，如下式所示：(32.27)现在就可以考虑对回归系数集的统计检验。最通常利用的检验是, 这个联合假设的检验。适宜的F统计量为：(32.28)具有k和Nk1自由度。较大的值，可使我们否认原假设。五、 reg回归过程 在SAS/STAT中有多个进展回归的过程，如reg、glm等，常用于进展一般线性回归模型分析的为reg过程。1. proc reg过程Reg过程一般由如下语句控制：proc reg data=数据集集名 ;model 因变量=自变量名列 ；var 变量列表；output out=数据集名 ;plot 绘图表达式 ；print 关键字列；weight 变量；freq 变量；by 变量；restrict 方程1，方程2， ；test 方程1，方程2，；run ;其中，model语句是必需要有的，其他语句都是可选的。2. proc reg 语句中的l outest=SAS数据集将有关模型的参数估计和选择的统计量输出到指定的SAS数据集中。l outsscp=SAS数据集要求把平方和与叉积矩阵输出到type=sscp的数据集中。l all屏幕输出所有内容。l usscp对用在该过程中的所有变量输出平方和与叉积矩阵。l noprint不在屏幕输出任何内容。3. model语句中的（1） 确定变量筛选方法的选择项l selection =none | forward | backward | stepwise | maxr | minr | rsquare | cp | adjrsq依次表示全部变量进入法none、前进法forward、后退法backward、逐步筛选法stepwise前进法与后退法的结合、最大R2增量法maxr、最小R2增量法minr、R2选择法rsquare、Mallows Cp 选择法cp、修正R2选择法adjrsq。（2） 其他选择项表3.2所示的是可在model语句中选用的其他选项。表32.2 model语句中的其他选项acovxpxspecpcorr1slentrydetailsaiccovbistbpcorr2slstaylackfitsbccorrbpcliscorr1startcollinss1mserclmscorr2bestcollinointss2ssebjpadjrsqincludeinfluencevifseqbdwrmsegmsepstoppartialtolallpcspnointsigmanoprintbic其中一些选择项的意义如下：l acov存在异方差时，输出参数估计量的渐近协方差阵的估计。l spec进展关于方差异性的检验。l slentry | sle =显著性水平规定入选变量进人方程的显著性水平。l slstay | sls=剔除水平规定从方程中剔除变量的显著性水平。l includen强迫前n个自变量进入模型。l starts以含有model语句中前3个自变量的模型开始，进展比拟、选择过程仅用于maxr或minr方法。l stops当找到最优的s个变量模型之后，逐步回归便停止(仅用于maxr或minr方法)。l p要求计算各观测点上因变量的预测值。l r作残差分析，同时给出因变量的预测值。l cli给出各自变量x0所对应的因变量y0的95置信上、下限。l clm给出各自变量所对应的因变量预测值均数Eyi i的95置信上、下限。l noint指明回归方程不带截距项常数项。l stb要求输出标准回归系数。l covb要求输出回归系数估计的协方差阵估计。l corrb要求输出回归系数估计的相关矩阵估计。l mse要求输出随机扰动项方差的估计。l rmse要求输出。l collin在对截距未进展校正的情形下，诊断多重共线性，条件数越大越可能存在共线性。l collinoint在对截距进展校正的情形下，诊断多重共线性。l tol表示共线性水平的容许值。对于某个变量容许值定义为1，其中是由这个变量和模型中所有其他回归变量建立的回归模型所得到的。tol越小说明其可用别的自变量解释的局部越多，自然就越可能与别的自变量存在共线性关系，tol与vif互为倒数。l vif输出变量间相关性的方差膨胀系数，vif越大，说明由于共线性的存在，使方差变大。l influence要求对异常点进展诊断。对每一观测点，输出如表32.3所示的统计量：表32.3 诊断异常点的统计量名称统计量含义“异常的判别准如此Leverage(hi)杠杆率hi，第i次观测自变量的取值在模型中作用的量度 0hi1hi越大，如此第i 次观测在模型中的作用就越大Cooks DCOOKD统计量，对某一观测点引起回归影响大小的度量。用于诊断异常点。假设D 50，如此可认为该观测点对模型的拟合有强的影响covratio协方差矩阵的行列式之比去掉某一观测点后、前比照假设| covratio | 3(自变量个数+i),如此第i个观测点值得引起注意defits此值大于2，明确该点影响较大debetas此值大于2，明确该点影响较大l i要求打印其中X 为设计矩阵。l xpx输出模型的叉积矩阵。l ss1要求打印第一类的模型参数估计的顺序平方和。l ss2要求打印第二类的模型参数估计的偏平方和。l all要求输出SAS所分析的以下选择项的特性：xpx，ss1，ss2，stb，covb，corrb，seqb，p，r，cli，clm，spec，acov，tol, pcorr1，pcor，r2，scorr1，scorr2。l partial给出每一回归变量的偏回归残差图。l dw一阶自相关检验的Durbin-Watson统计量。4. 其他选择语句l output语句用于把一些计算结果输出到指定的数据集中。有关的关键字与其意义如表32.4所示。表32.4 reg过程的output语句中的关键字关键字意义关键字意义关键字意义predicted预测值l95m95clm下限stdpclm的标准差residual残差u95m95clm上限stdr残差的标准差press残差/(1hi)l9595cli下限stdicli的标准差rstudent刀切残差u9595cli上限cookedCook D统计量student学生氏残差h杠杆点统计量hil var语句列出叉积矩阵中的变量，仅当具有outsscp=sasdataset这个选择时才使用。l plot语句绘制两变量的散点图。语句格式为：plot x*y / 选项。其中x和y变量，可以是原始数据集中的变量，也可以是统计量关键字。假设变量是统计量关键字时，需要在其后加上一个小圆点“。l restrict语句要求计算线性等式约束的最小二乘估计，其中的方程就是关于回归系数用自变量表示的等式，方程与方程间用逗号分隔。例如，对于模型model y=a1 a2 b1 b2，可以用restrict a1+a2=1语句，表示参数估计是在a1+a2=1的条件下，求最小二乘估计。l test语句要求进展线性等式约束的显著性检验，即Tintner检验，其中的方程就是关于回归系数用自变量表示的等式，方程与方程间用逗号分隔；test语句一般不与restrict语句同用。例如，对于模型model y=a1 a2 b1 b2，可以用test a1+a2=1语句，表示在a1+a2=1原假设条件下作F检验。5. 交互式语句下面的这局部语句可以用在proc reg过程中，但常用在reg过程激活后，以交互方式运行。l add 变量名列表向模型中增加变量。l delete 变量名列表删除原拟合模型中的有关变量。l refit重新拟合模型。l print输出有关模型的相关信息。reg过程其详细用法可参阅SAS/STAT的用户手册。六、 实例分析 表32.5列举了一个班级的学生情况的调查数据，试分析身高对体重的影响。表32.5 bclass记录数据name某某age年龄sex性别height身高(厘米)weight体重(公斤)name某某age年龄Sex性别height身高(厘米)weight体重(公斤)KATE12女145FREDRICK14男154LOUISE12女149ALFRED14男157JANE12女135HENRY14男159JACLYN12女162LEWIS14男157LILLIE12女127EDWARD14男167TIM12男147CHRIS14男157JAMES12男149JEFFERY14男169ROBERT12男125MARY15女152BARBARA13女147AMY15女157ALICE13女149ROBERT15男164SUSAN13女137WILLIAM15男159JOHN13男159CLAY15男162JOE13男154MARK15男152MICHAEL13男142DANNY15男162DAVID13男145MARTHA16女159JUDY14女149MARIAN16女147ELIZABET14女152PHILLIP16男167LESLIE14女159LINDA17女152CAROL14女154KIRK17男167PATTY14女152LAWRENCE17男172分析和操作步骤过程如下。1. 建立数据文件首先要将表32.5中的数据输入到SAS数据集中，可调用SAS的数据步data过程，建立我们所需的bclass数据集。程序如下：data study.bclass ;input name $ age sex $ height weight ;cards ;KATE 12 F 145 43.1 LOUISE 12 F 149 55.8 ;run ;2. 制作变量的散点图建立完SAS数据集bclass后，一般需要对数据集中要分析的变量weight与height制作散点图，以便能从图示中反映学生的身高与体重的关系。一般的处理操作有菜单操作方法和编程方法2种。如果用菜单操作方法，在SAS/Assist环境中，从Primary Menu主菜单中选择Graphics/High resolution/ Plots/Simple x*y plot菜单命令，再选择Active data set为study.bclass，Vertical axis为weight，Horizontal axis为height，可以在additional options选项菜单中通过Line and Symbol子选项选定所需要的连线类型和点的符号等，最后选择Locals/Run菜单命令，提交运行即可显示图形。如果用编程方法，程序如下：goptions reset=global gunit=pct cback=white border htitle=6 htext=3 ftext=swissb colors=(back);proc gplot data=study.bclass ;plot weight * height ;run ;运行后，在Graph窗口得到如图32.1所示的结果。 体重与身高(weight与height)的散点图3. 相关系数计算如果用菜单操作方法，可选择Globals/SAS/Assist/Data Analysis/Elementary/Correlation命令，再选择Active data set为study.bclass，Columns to be correlated为weight和height，然后提交运行。直接编写调用相关系数计算的程序为：proc corr data= study.bclass ;var weight height ;run ;运行后，在Output窗口得到如表32.6所示的结果。表32.6 身高与体重(weight与height)的相关系数Correlation Analysis2 VAR Variables: WEIGHT HEIGHTSimple StatisticsVariable N Mean Std Dev Sum Minimum MaximumPearson Correlation Coefficients / Prob |R| under Ho: Rho=0 / N = 40 WEIGHT HEIGHT 从输出表32.6可以看出，身高与体重之间的相关系数为0.70844。4. 回归分析如果用菜单操作方法，可选择Globals/SAS/Assist/Data Analysis/Regression/Linear regression命令，再选择Active data set为study.bclass，Dependent为weight，Independent为height，然后提交运行。编程实现回归方法为：proc reg data= study.bclass ;model weight = height /r clm cli dw;run ;其中，模型参数r表示要输出残差分析，包括因变量的观察值、由输入数据和估计模型来计算的预测值、残差值、标准误差、学生化残差、COOKD统计量。模型参数clm表示对每个观察输出因变量期望值的95%置信上界和下界，仅考虑到参数估计的偏差，没有考虑误差项的偏差。模型参数cli表示对因变量的各个预测值输出95%置信上界和下界，这个置信界反映了误差的偏差以与参数估计的偏差。模型参数dw表示要进展误差项的独立性检验，计算Durbin-Watson统计量。运行后，在Output窗口得到如表32.7所示的结果。表32.7 回归分析结果 Model: MODEL1 Dependent Variable: WEIGHT Analysis of Variance(方差分析) Sum of Mean Source DF Squares Square F Value ProbF Parameter Estimates(参数估计) Parameter Standard T for H0: Variable DF Estimate Error Parameter=0 Prob |T| INTERCEP 1 HEIGHT 1 误差项的独立性检验(For Number of Obs.) 40置信区间 Dep Var Predict Std Err Lower95% Upper95% Lower95% Upper95% Std Err Obs WEIGHT Value Predict Mean Mean Predict Predict Residual Residual残差分析 Student Cooks Obs Residual -2-1-0 1 2 DSum of Residuals 0回归分析根据所选择的模型参数的输出，分为假设干段，下面逐个段地给以说明：方差分析表提供关于拟合模型的一般信息。总观察数N=40，自变量个数k=1，回归模型带有截距i=1。回归模型的离差平方和RSS=1986.48457，自变量的个数k=1，所以自由度df=k=1，计算公式见式(31.29)。因变量的样本离差平方和TSS=3958.05375，自由度为df=N1=401=39，计算公式见式(31.34)。误差项的样本离差平方和ESS=1971.56918，自由度df=Nk1=40-1-1=38，计算公式见式(31.32)。注意TSS=RSS+ESS，即3958.05375=1986.48457+1971.56918。回归模型的离差平方和平均值MSR=RSS/df=1986.48457/1=1986.48457，误差项的离差平方和平均值MSE=ESS/df=1971.56918/38=51.88340。在原假设所有自变量的回归系数都为0的情况下，本例只有一个自变量height，即H0：，F(1,38)=MSR/MSE=1986.48457/51.88340=38.287，查F分布表，p值为0.0001小于显著水平0.05，明确可拒绝原假设，并有足够的证据断定回归线的斜率不为零。所以，这一模型拟合数据比基线模型好。无偏的误差估计标准值Root MSE=7.20301，因变量weight平均值Dep Mean=47.66250，变异系数或称方差系数CV=Root MSE/Dep Mean100=15.11254，它表示与单位无关的方差。RSquare是01之间的值，它表示贡献给模型而不是贡献给拟合残差的总方差的那局部，它也称为决定系数或拟合优度，用于判断回归模型拟合好坏。R2=1ESS/TSS=RSS/TSS=1986.48457/3958.05375=0.5019，调整 R2=1ESS/TSS(N-i)/N-k-i39/38=0.4888，R2越是接近1说明模型拟合得越好，等于1如此说明完全拟合，没有任何信息丢失，本例的R2值明确有一半信息丢失没有被回归模型表示出来，通常R2应该超过0.7以上才比拟好。参数估计表给出截距和斜率的估计值，方程明确截距的估计值为56.748575，斜率的估计值为0.681312，计算公式见式(31.17)和式(31.19)。估计截距的标准误差计算公式见式(31.37)，其中，自变量height的平均值=153.25，自变量height的离差平方和=4279.5，估计误差51.88340，所以估计截距的标准误差=16.912396。在截距等于零的原假设下，计算出的t(38)=56.748575/16.912396=3.355，大于此临界点绝对值出现的概率为0.0018，远远地小于5%，有充足的理由否决截距为零的原假设。估计斜率的标准误差计算公式见式(31.38)，估计斜率的标准误差=0.1101077，在斜率等于零的原假设下，计算出的t0.1101077=0.9042135。置信区间分析，输出了weight因变量Dep Var的40条原始观察值和回归模型的预测均值Predict Value，与预测均值的标准差Std Err Predict、预测均值的置信区间下界Lower95% Mean和上界Upper95% Mean、预测值的置信区间下界Lower95% Predict和上界Upper95% Predict、残差Residual、残差的标准差Std Err Residual。我们以第一条观察Obs=1为例来说明计算过程，第一条的观察=43.1，=145，根据回归模型最小二乘法计算出的估计参数，可以得到预测均值为145=42.0417。第一条观察的杠杆率计算公式见式(31.42)，=0.040904311，所以预测均值的标准差=1.457=44.9908。预测值的方差除了要考虑参数估计的偏差，还要考虑误差项的偏差，所以要在预测均值的偏差上加上一个误差项的偏差，计算公式见式(31.44)，预测值的标准差=7.34885394=56.9187，我们从上面的置信区间计算中可以发现两个知识点，第一个知识个点，预测值的置信区间要大于预测均值的置信区间，第二个知识点，越是接近自变量heightheight观察值，它的因变量weightheight观察值，它的因变量weight预测均值和预测值的置信区间越是宽，从图形上直观地看置信区间为中间窄，两头形成喇叭口。残差分析，我们仍然以第一条观察为例来说明计算过程。残差=43.100042.0417= 1.0583。标准残差的计算公式见(31.46)式，标准残差=7.054，学生化残差Student Residual=残差/标准残差=1.0583/7.054=0.150。由于学生化残差服从标准正态分布，将学生化残差画在残差图上，我们可以清楚地看到大约68%的学生化残差值落在一个标准差1到1之间，而大约95%学生化残差值落在两个标准差2到2之间。根本上认为模型的误差项服从正态分布与满足同方差假设，在诊断上没有太大问题。残差之和=0，残差的平方和=1971.5692。COOKD统计量用于预测每个观察点是否为强影响点或称异常点，它是通过删除这个观察点后重新用最小二乘估计求解参数值，来分析这个观察点。观察点的COOKD统计量小于50%，我们认为不存在异常情况。PRESS统计量是预测残差的平方和，第i个观察的残差定义为，其中，为删除第i个观察后从余下的组数据中重新用最小二乘法求出的参数估计而计算出的第i个观察的预测值。第i个观察的预测残差为。误差的独立性检验，它是回归模型的三大假设之一。我们采用针对残差一阶自相关性进展计算的Durbin-Watson统计量来检验，计算公式见式(31.48)，相邻残差之差的平方和=2899.603，DW=2899.603/1971.56918=1.471，DW值靠近2说明误差根本上是独立的，小于2说明是正相关。残差一阶自相关系数=0.185，接近0也说明了误差根本上是独立的。残差一阶自相关系数的计算方法与一般的相关系数计算公式类似，残差值的第一个序列数据为第1个残差到第39个残差，第二个序列数据为第2个残差到第40个残差，第一、二个序列残差数据的平均值为0，标准化时公式的分母取1到40个残差值，即。5. 输出带有回归线的散点图如果我们需要输出带有回归线的散点图，菜单操作方法是通过在additional options选项菜单中选择Regression Plots/Plots of dependent by independent columns命令，重新再提交一次。注意，此时还可以同时选择输出残差图。程序的方法是在proc reg过程里增加plot语句，要注意SAS的关键字使用在plot语句中时要加小圆点，这里是预测值p关键字，增加的plot语句如下：plot weight * height=+ p.* height=*/ overlay ;如果我们需要输出高分辨率的回归线图形，可以先在reg过程中将拟合的预测值p输出到一个SAS数据集如bclassg中，再调用gplot过程绘制图形。增加的output语句如下：outputout=study.bclassg p=predict l95=clil95 u95=cliu95;绘制高分辨率的带有回归线的散点图程序如下：goptions reset=global gunit=pct cback=white border htitle=6 htext=3 ftext=swissb colors=(back);proc gplot data=bclassg ;plot weight*height predict*height clil95*height cliu95*height/overlay;symbol1 v=plus c=red i=none h=2.5;symbol2 i=spline v=none c=blue;symbol3 i=spline v=none c=red l=3;symbol4 i=spline v=none c=black l=3;run ;注意，我们也可以用图形自带i=rlcli95选项，直接绘制预测值的置信区间上下界。运行后，在Graph窗口得到如图32.2所示的结果。图32.2 带有回归线、95%置信线的体重与身高(weight与height)散点图从绘制出的带有回归线的图形可形象地看出模型拟合数据比只用均值预测的基线模型好。仔细观察图形，两条95%的上下预测值置信带呈现两头喇叭口。此外，还可用性别来分组，分别对男生和女生进展回归分析，分别建立男生和女生的回归模型。 研究耗氧量模型。这是有关身体适应性测试的例子，肺活量与一些简单的锻炼测试数据的拟合，目的是为了在锻炼测试的根底上而不是在昂贵笨重的氧气消耗测试的根底上得到方程来预测适应性。由于回归是相关的，因此，理论上还应该请求共线性诊断。该数据名为fitness，这是一个对31位成年人心肺功能的调查结果，它包含的变量见表32.8，测试的各项数据见表32.9。表32.8 fitness数据集的变量名变量名含 义age年龄weight体重oxygen耗氧量runtime跑15英哩的时间分rstpulse休息时每分钟心跳次数runpulse跑步时每分钟心跳次数maxpulse每分钟心跳次数最大值表32.9 fitness数据集中的测试数据ageweight oxygen runtime rstpulse runpulse maxpulse 4462178182406218518544451561684240166172385517818047581761764070176180436416217044631741763848170186444516816845561861924551176176474716216454501661704944180185515716817251481621684848162164497616816857581741765462156165524816416650481461555148172172544416817251591861885749148155495618618848521701765253170172在这个锻炼测试数据里，我们感兴趣的是耗氧量是如何依赖于其他变量的。1. 建立数据文件程序如下：data fitness ;input age weight oxygen runtime rstpulse runpulse maxpulse;cards ; 44 89.47 44.609 11.37 62 178 18240 75.07 45.313 10.07 62 185 18552 82.78 47.467 10.50 53 170 172;run ;2. 制作变量的散点图fitness数据集中的变量较多，我们需要制作每两个不同变量oxygen、age 、weight、runtime、rstpulse、runpulse和maxpulse之间的所有散点图，即散点图矩阵，共有76=42个散点图。我们可以通过在SAS/Insight软件中绘制散点图矩阵，操作步骤为：在SAS命令框中键入insight后按Enter，在SAS/Insight：Open对话单中，选择work.fitness数据集后单击Open按钮，将在屏幕的窗口中显示当前打开的数据集work.fitness内容，再选择菜单上的Analyze/Scatter PlotY X命令，在出现的Scatter PlotY X对话单中，把fitness数据集中的7个变量依上面的次序全部参加Y轴和X轴的列表框中，最后单击OK。如图32.3所示。图32.3 fitness变量间的散点图矩阵散点图矩阵图中第一行的六个散点图分别表示oxygen变量作为y轴,其他六个变量作为x轴的散点图,第一列的六个散点图分别表示oxygen变量作为x轴,其他六个变量作为y轴的散点图。对角线上第一个oxygen小方格里的左下角数字37.388和右上角数字60.055,表示变量oxygen的数值变化X围。其他行和列中的散点图所代表的两个变量类同。绘制散点图矩阵图是为了研究变量间的相关性，从图32.3中，我们发现变量runpulse与maxpulse之间存在有较强的共线性，如果在回归模型中增加方差膨胀系数vif，共线性水平的容许值tol，条件数collin选项对回归进展诊断，也会得到一样的结论。另外，我们从图中还发现耗氧量oxygen与变量runtime有较强的负相关。从下面的相关系数计算中也能得到这些一样认识。3. 相关系数计算编写的相关系数计算程序如下：proc corrdata= fitness;var oxygen age weight runtime rstpulse runpulse maxpulse;label oxygen = Oxygen consumptionage = Age in years weight = weight in kg runtime = Min. to run 1.5 miles rstpulse= Heart rate while resting runpulse= Heart rate while running maxpulse= Maximum heart rate;run ;运行后，在Output窗口得到表32.10所示的结果。表32.10 fitness变量间的相关系数Correlation Analysis 7 VAR Variables: OXYGEN AGE WEIGHT RUNTIME RSTPULSE RUNPULSE MAXPULSESimple StatisticsVariable N Mean Std Dev Sum Minimum Maximum LabelOXYGEN 31 47.3758 5.3272 1468.7 37.3880 60.0550 Oxygen consumption AGE 31 47.6774 5.2114 1478.0 38.0000 57.0000 Age in years WEIGHT 31 77.4445 8.3286 2400.8 59.0800 91.6300 Weight in kg RUNTIME 31 10.5861 1.3874 328.2 8.1700 14.0300 Min. to run 1.5 miles RSTPULSE 31 53.7419 8.2944 1666.0 40.0000 76.0000 Heart rate while resting RUNPULSE 31 169.6452 10.2520 5259.0 146.0000 186.0000 Heart rate while runningMAXPULSE 31 173.7742 5387.0 155.0000 192.0000 Maximum heart ratePearson Correlation Coefficients / Prob |R| under Ho: Rho=0 / N = 31 OXYGEN AGE WEIGHT RUNTIME RSTPULSE RUNPULSE MAXPULSE OXYGEN 1.00000 -0.30459 -0.16275 -0.86219 -0.34641 -0.39797 Oxygen consumption 0.0 0.0957 0.3817 0.0001 0.0563 0.0266 AGE Age in years 0.0957 0.0 WEIGHT Weight in kg RUNTIME Min. to run 1.5 miles 0.0001 0.3092 0.4412 0.0 RSTPULSE Heart rate while resting RUNPULSE Heart rate while running 0.0266 0.0630 0.3284 0.0858 0.0813 0.0 MAXPULSE Maximum heart rate 0.1997 0.0150 0.1761 0.2213 0.1620 0.0001 0.0 4. 回归分析编写的回归分析程序如下：proc regdata= fitness;model oxygen=age maxpulse rstpulse runpulse runtime weight/ss1 ss2 ;run ;Model语句中增加ss1和ss2选项，回归模型将计算每个自变量两种不同类型的平方和，其中，ss1是按model语句中自变量的排列顺序依次计算每个自变量的平方和，也称为第一类平方和或称顺序平方和，ss2是把model语句中每个自变量排到变量列表的最后，所计算的第一类平方和，也称为第二类平方和或称局部平方和。通过分析每个自变量的这两类平方和，能知道回归模型总的平方和的构成和各个自变量所贡献的平方和，进而能知道哪些自变量是最重要的回归变量，哪些回归变量可能是无关紧要的，配合参数估计的t检验，最终为缩减回归变量提供依据，达到简化模型的目的。运行后，在Output窗口得到表32.11所示的结果。表32.11 fitness数据的回归分析结果 Model: MODEL1Dependent Variable: OXYGEN Analysis of Variance Sum of Mean Source DF Squares Square F Value ProbF Parameter Estimates Parameter Standard T for H0:Variable DF Estimate Error Parameter=0 Prob |T| Type I SS Type II SSINTERCEP 1 102.238339 12.45304719 8.210 0.0001 69578 AGE 1 -0.219916 0.09959154 -2.208 0.0370 MAXPULSE 1 0.304735 0.13722472 2.221 0.0361 RSTPULSE 1 -0.000844 0.05863130 -0.014 0.9886 RUNPU