《多元统计分析及R语言》第4章多元相关与回归分析

问题：

如何确定身高和体重之间相关关系，如果存在关系性，如何建立数据模型对身高或体重进行预测？分析方法：相关分析和回归分析问题提出身高171175159155152158154164168166159164体重576441383544415157494746相关与回归相关关系平行关系相关分析依存关系回归分析函数关系有精确的数学表达式相关与回归一元相关分析多元相关分析复相关典型相关一元回归多元回归第4章多元相关与回归分析4.1变量间的关系分析4.2多元线性回归分析4.3多元线性相关分析4.4回归变量的选择方法4.1变量间的关系分析1.简单相关分析结论：相关系数|ρ|≦1,绝对值越接近1，相关性越密切，接近0则相关性弱。计算公式：Pearson相关系数：(1)两变量线性相关系数计算(2)相关系数检验结论：说明该组人群的身高与体重之间呈现正的线性相关关系。相关系数显著性需进行假设检验。身高171175159155152158154164168166159164体重576441383544415157494746lxy<-function(x,y){n=length(x);sum(x*y)-sum(x)*sum(y)/n}#建立离均差乘积和函数lxy(x1,x1)#x1的离均差平方和[1]556.9lxy(x2,x2)#x2的离均差平方和[1]813lxy(x1,x2)#x1与x2的离均差乘积和[1]645.5(r=lxy(x1,x2)/sqrt(lxy(x1,x1)*lxy(x2,x2)))#用离均差乘积和计算相关系数r=0.9593结论：说明该组人群的身高与体重之间呈现正的线性相关关系。相关系数显著性需进行假设检验。身高171175159155152158154164168166159164体重576441383544415157494746cor(x,y=NULL,method=c("pearson","kendall","spearman"))method计算方法，包括"pearson","kendall"或"spearman"三种，默认pearsonr=0.9593t检验统计量来自ρ=0的总体所有样本相关系数呈对称分布，故r的显著性可用t检验来进行。（2）相关系数的假设检验—t检验n=length(x1)#向量的长度tr=r/sqrt((1-r^2)/(n-2))#相关系数检验t统计量p=2*（1-pt(tr,df=10,0.05)）=8.21×10-7#双侧检验结论：p<0.05,拒绝原假设，接受备则假设，ρ不显著为0。身高和体重正相关。相关与回归相关关系平行关系相关分析依存关系回归分析函数关系有精确的数学表达式相关与回归一元相关分析多元相关分析复相关典型相关一元回归多元回归2.一元线性回归模型和参数估计实际值与估计值离差平方和：最小二乘法，对a，b求偏导：此处lxx表示的离差平方和，lxy表示与的离差积和。基本模型形式：参数估计（最小二乘法）(1)模型系数估计(2)模型和系数检验2.直线回归方程的建立步骤（2）计算斜率:b=lxy(x,y)/lxy(x,x)=-140.364（1）确定变量：x=x1#自变量,y=x2#因变量,（3）计算截距：a=mean(y)-b*mean(x)=1.159（4）建立回归方程（5）plot(x,y)；lines(x,a+bx)3.回归模型假设检验：方差分析离差平方和的分解xyy{}}残差平方和(SSE)扣除x对y线性影响后变异回归平方和(SSR)x对y的线性影响离差平方和(SST)与自由度结合：F统计量为其中结论：p<0.05,拒绝原假设，接受备则假设，模型参数不显著同时0。模型有意义。假设条件：原假设:b=0，模型无意义;

备择假设：b≠0,模型有意义p=1-pf(F,df1=1,df2=(n-2))=8.21×10-7自变量个数n-1-自变量个数结论：p<0.05,拒绝原假设，接受备则假设，模型参数不显著为0。模型有意义。t统计量为：其中Sb标准化为：3.回归系数假设检验：t检验假设条件：原假设:β=0;备择假设：β≠04实例分析

财政收入与税收有密切的依存关系。收集了我国1978年改革开放以来到2008年共31年的税收(x,百亿元)和财政收入(y，百亿元)数据。tyx197811.32625.1928197911.46385.3782198011.59935.7170198111.75796.2989198212.12337.0002………1.数据形式问题：分析税收与财政收入之间的依存关系。

yx=read.table(“clipboard”,header=T)#加载数据(fm=lm(y~x,data=yx))#显示一元线性回归模型lm(formula=y~x)Coefficients:(Intercept)x-1.1971.1162.模型估计线性回归方程为：plot(y~x,data=yx)#做散点图abline(fm)#添加回归线3.模型检验结论：F统计量对应的p<0.05,则拒绝原假设，模型系数不显著为0，模型有统计学意义假设条件：原假设:β=0;备择假设：β≠0

anova(fm)#模型方差分析 DfSumSqMeanSqFvaluePr(>F)x1712077712077

=27428<2.2e-16***Residuals29753264.模型系数检验结论：t统计量对应的p<0.05,则拒绝原假设，模型系数不显著为0，模型有统计学意义注意：t2=F，当自变量只有1个时，方差分析与检验的结果是等价的。但在多元分析中，方差分析与检验的结果并不等价假设条件：原假设:β=0;备择假设:β≠0summary(fm)#回归系数t检验变量回归系数标准误t值P值截距-1.19656/1.1612

=-1.0304

0.3113x1.1162/0.0067=165.6144<2e-16***F-statistic:27428.133on1and29DF,p-value:<2.22045e-164.2多元线性回归分析1.多元线性回归模型建立模型基本公式：n组观测值回归模型：矩阵形式：2.基本假设4、n>p.即要求样本容量个数多于解释变量的个数。1、解释变量（自变量）一般说来是非随机变量2、误差等方差及不相关假定（G-M条件）：3、误差正态分布的假定条件为：3.多元线性回归最小二乘估计模型基本公式：对参数β求偏导，使其实际值与估计值离差平方和最小：β估计值为：结论：在正态性假定条件下，β的最小二乘估计与极大似然估计值完全相同。4.模型检验结论：F统计量对应的p<0.05,则拒绝原假设，模型有意义，否则，接受原假设，模型无意义。其中：SSR为回归平方和，SSE为残差平方和。假设条件：原假设:β1=β2=…=βn=0;备择假设：β1,β2,…,βn

不全为0离差平方和分解：检验F统计量：5.模型估计参数检验结论：t统计量对应的p<0.05,则拒绝原假设，系数不显著为0，有意义，否则，接受原假设，系数无意义。假设条件：原假设:βj=0;备择假设:βj≠0检验t统计量：其中：为估计值为第j个偏回归系数的标准误差4实例分析财政收入:反映一国经济实力的重要标志.分析财政收入(y百亿元)和国内生产总值(x1百亿元)、税收(x2百亿元)、进出口贸易总额(x3百亿元)、经济活动人口(x4百万人)之间的关系。tyx1x2x3x4197811.326236.2415.19283.55406.82197911.463840.3825.37824.12415.92198011.599345.1785.71705.70429.03198111.757948.6036.29898.904441.65198212.123353.0187.000212.801456.74198318.669559.5747.555915.903467.071.数据形式结论：为了消除单位数量级影响，需要对数据进行标准化或归一化处理标准化方法：

yX=read.table("clipboard",header=T)#加载例4.4数据

YX=data.frame(scale(yX))#标准化(fm=lm(y~x1+x2+x3+x4,data=yX))#显示多元线性回归模型lm(formula=y~x1+x2+x3+x4)Coefficients:(Intercept)x1x2x3x423.532109-0.0033871.1641150.000292-0.0437422.模型估计多元线性回归方程为：结论：由于前期未对数据进行标准化处理，此时需要对系数进行后期标准化处理。3.模型估计系数标准化结论：税收对财政收入线性影响最大标准化公式为：其中：si和sy分别是各自变量和因变量的标准差。library(mvstats)coef.sd(fm)#标准化偏回归系数结果$coef.sdx1x2x3x4-0.01741.04230.0009-0.03714.模型检验结论：F统计量对应的p<0.05,则拒绝原假设，模型系数不同时为0，模型有统计学意义假设条件：原假设:β1=β2=…=βn=0;备择假设：β1,β2,…,βn

不全为0来源离均差平方和自由度均方F值P值回归712627.36/4=178156.84=22890.80

<.0001误差202.3541/26=7.7829总计712829.718630残差正态性检验：shapiro.test(fm$residuals)#原假设正态性4.模型系数检验结论：税收x2和经济活动人口x4的p<0.05，拒绝原假设，说明系数不显著为0，这两个因素对财政收入有显著影响，国内生产总值x1和进出口贸易x3未通过显著性检验，未对财政收入有显著影响，可以去掉。另：结果与经济现实不符，需要消除共线性或去除部分变量假设条件：原假设:βj=0;备择假设:βj≠0变量回归系数标准误t值P值标准回归系数x023.53214.5995.122.5e-05……x1-0.0033870.0081-0.420.68-0.01745x21.164110.040528.75<2-e-161.04235x30.0002920.00850.030.950.00096x4-0.043740.0092-4.727.0e-05-0.037104.3多元线性相关分析1.矩阵相关分析n个变量：相关系数：相关系数公式：财政收入:反映一国经济实力的重要标志.分析财政收入(y百亿元)和国内生产总值(x1百亿元)、税收(x2百亿元)、进出口贸易总额(x3百亿元)、经济活动人口(x4百万人)之间的关系。2.实例分析问题：分析财政收入与各变量间的相关性tyx1x2x3x4197811.326236.2415.19283.55406.82197911.463840.3825.37824.12415.92198011.599345.1785.7175.70429.03198111.757948.6036.29898.904441.65198212.123353.0187.000212.801456.74198318.669559.5747.555915.903467.07结论：财政收入与和国民生产总值及税收、进出口贸易总额关系显著，与经济活动人口关系不够显著。yx1x2x3x4y10.98710.99950.99120.6957x10.987110.99070.98680.7818x20.99950.990710.99170.7154x30.99120.98680.991710.7074x40.69570.78180.71540.70741cor(yX)#多元数据相关系数矩阵3.复相关分析定义：研究多个变量同时和某个变量间的相关关系，采用复相关系数对相关性进行度量回归模型决定系数：研究因变量与一组自变量间的相关程度。计算公式：结论：R2称为决定系数，体现模型拟合程度好坏，一般大于0.8就表示模型拟合程度好。计算公式：分析财政收入(y百亿元)和国内生产总值(x1百亿元)、税收(x2百亿元)、进出口贸易总额(x3百亿元)、经济活动人口(x4百万人)之间的关系。2.实例分析问题：分析财政收入与各变量间的回归模型拟合程度。tyx1x2x3x4197811.326236.2415.19283.55406.82197911.463840.3825.37824.12415.92198011.599345.1785.7175.70429.03(R2=summary(fm)$r.sq)#多元线性回归模型决定系数[1]0.9997(R=sqrt(R2))#多元数据复相关系数[1]0.9999结论：模型可决系数大于0.9，拟合效果好，但随着变量增加，决定系数会增加，故不能仅仅根据其判断模型好坏4.4回归变量选择方法一.变量选择标准选择标准常用准则平均残差平方和最小准则误差平方根MSE最小准则校正复相关系数平方准则Cp准则AIC准则和BIC准则全局择优法残差平方和(RSS)最小复相关系数平方(R2)和最大准测1.全局择优法定义:从所有可能回归模型中，按一定准则选择最优模型

残差平方和(RSS)最小复相关系数平方(R2)和最大准测：模型总个数：library(leaps)#加载leaps包varsel=regsubsets(y~x1+x2+x3+x4,data=yX)#线性回归变量选择模型result=summary(varsel)#变量选择方法结果data.frame(result$outmat,RSS=result$rss,R2=result$rsq)#结果展示x1x2x3x4RSSR21(1) *752.880.998942(1)**203.880.999713(1)***202.350.999724(1)****202.340.99972结论:自变量子集越多,RSS越小,R2越大,总是会选择自变量最多,模型,且当变量较多时,产生模型个数会很多(2p),计算量大.2.变量选择常用标准结论：随着模型变量个数增加，采用(n-p)-1方式对其进行惩罚，故RMSp越小越好。（1）平均残差平方和最小准则：结论：MSE即为剩余标准差sy.x，故MSEp越小越好。（2）误差均方根最小准则：结论：等价误差平方根最小，故adjR2越大越好。（3）校正复相关系数：结论：按AIC和BIC越小越好准则进行模型选择。结论：选择对应点（p,Cp）最接近第一象限角平分线，且Cp最小（4）Cp准则：（5）AIC和BIC准则：实例data.frame(result$outmat,adjR2=result$adjr2,Cp=result$cp,BIC=result$bic)#调整决定系数,Cp和BIC准则结果展示x1x2x3x4adjR2CpBIC1(1)*0.998969.745-205.62(1)**0.99971.199-242.63(1)***0.99973.001-239.44(1)****0.99975.000-236.0结论：按最小标准，选择x2和x4两个变量的模型较好。4.4.2.逐步回归分析思想：当自变量之间相关性很强时，会使模型方程无解或病态，需要将这部分变量删除，采用逐步筛选变量，寻找较优子空间的变量选择方法可以实现这一过程。三种方法向前引入法（forward）向后剔除法（backward）逐步筛选法（stepwise）1.向前引入法：只引入不剔除向前引入法x1,yx2,yxp,y分别建立p个一元回归模型，设置引入门槛(F进=3.87)分别对p个模型回归系数进行F检验选择最大F统计量且大于F进值引入模型继续引入操作，直到F统计量<F进为止…缺点：没有考虑引入变量引入方程后变为不显著的情况2.向后剔除法：只剔除不引入向后剔除法分别对模型p个回归系数进行F检验选择最小F统计量且小于F出值剔出模型继续剔除操作，直到F统计量>F出为止对所有自变量和因变量建立一个多元回归模型，设置剔除门槛(F出)缺点：计算量大3.逐步筛选法逐步筛选法计算方程与其他变量偏回归平方和，选择最大F进值引入方程检验引入方程模型，选择变量系数最小F统计量且小于F出值剔出模型选择相关性最强的1个变量建立一元回归方程，设置引入（F进）和剔除门槛(F出)说明:由于变量间存在相关性，单一引入或剔除没有考虑这种相关性，影响模型变量质量。重复操作，直到既不能引入也不能剔除为止4.实例分析fm=lm(y~x1+x2+x3+x4)#多元数据线性回归模型fm.step=step(fm,direction="forward")Start:AIC=68.15y~x1+x2+x3+x4（1）向前引入法变量选择（2）向后剔除法变量选择fm.step=step(fm,direction=“backward”,test=“F”)Start:AIC=68.15y~x1+x2+x3+x4DfSum

ofSqRSS

AICFPr(>F)-x310.009

202

660.00120.9730-x111

204

660.17590.6784<none>

202

68-x41174

376

8522.29547.005e-05***-x216433

6635174826.6460<2.2e-16***（2）向后剔除法变量选择fm.step=step(fm,direction=“backward”,test=“F”)Step:AIC=66.16y~x1+x2+x4

DfSumofSqRSSAICFvaluePr(>F)-x11

2

204640.20530.6541<none>

20266-x41

197

4008526.32502.143e-05***-x21

7382

7585176985.0464<2.2e-16***Step:AIC=64.39y~x2+x4

Df

SumofSq

RSSAICFvaluePr(>F)<none>

20464-x41

549

75310375.3961.972e-09***-x21

36765536785929550491.307<2.2e-16***（3）逐步筛选法变量选择fm.step=step(fm,direction=“both”,test=“F”)Start:AIC=68.15y~x1+x2+x3+x4DfSumofSqRSSAICFvaluePr(>F)-x310.009202660.00120.9730-x111204660.17590.6784<none>20268x411743768522.29547.005e-05***-x2164336635174826.6460<2.2e-16***（3）逐步筛选法变量选择fm.step=step(fm,direction=“both”,test=“F”)Step:AIC=66.16y~x1+x2+x4DfSumofSqRSSAICFvaluePr(>F)-x112204640.20530.6541<none>20266+x310.009202680.00120.9730-x411974008526.32502.143e-05***-x213827585176985.0464<2.2e-16***Step:AIC=64.39y~x2+x4DfSumofSqRSSAICFvaluePr(>F)<none>20464+x112202660.20530.6541+x310.18204660.02370.8789-x4154975310375.39621.972e-09***-x2136765536785929550491.3069<2.2e-16***4.5非线性回归模型可直线化曲线方程多项式曲线:y=a+bx+cx2+…+pxp对数曲线:y=a+blogx指数函数:y=aebx(a>0)幂函数:y=axb(a>0)双曲线函数:y=a+b/x曲线回归:通过两个相关变量x与y的实际观测值数据建立曲线回归方程，揭示x与y曲线关系类型，此处指可直线化非线性回归。经济增长模型:Y=A0emtLaKbCobb-Douglas生产函数:Q=aLb1Kb2+e消费函数:C=b1+b2Yb3+e不能精确直线化经典经济非线性模型非线性回归方程1、多项式曲线直线化多项式方程：直线化后方程：方法：令x1=x;x2=x2,…，xp=xpx1=x；x2=x^2lm.2=lm(y~x1+x2)#二次函数回归模型summary(lm.2)$coef#二次函数回归模型系数EstimateStd.ErrortvaluePr(>|t|)(Intercept)6.914690.33198720.8286.346e-09x1-0.465630.056969-8.1731.864e-05x20.010760.0020095.3534.604e-04决定系数；summary(lm.2)$r.sq=0.9513实例：分析百货商店销售额x与流通费率y之间的关系。结论：系数和模型都通过显著性检验，模型拟合效果好2、对数曲线直线化对数曲线方程：y=a+blogx直线化后方程：y=a+bx’方法：令x’=logx;lm.log=lm(y~log(x))#对数函数回归模型summary(lm.log)$coef#对数函数回归模型系数EstimateStd.ErrortvaluePr(>|t|)(Intercept)7.3640.168843.649.596e-13log(x)-1.7570.0677-25.951.660e-10summary(lm.log)$r.sq=0.9854实例：分析百货商店销售额x与流通费率y之间的关系。结论：系数和模型都通过显著性检验，模型拟合效果好.

模型为:y=7.364-1.757logx3、指数曲线直线化指数曲线方程：y=aebx(a>0)直线化后方程：y’=a+bx方法：两边取对数logy=loga+bxlm.exp=lm(log(y)~x)#指数函数回归模型summary(lm.exp)$coef#指数函数回归模型系数EstimateStd.ErrortvaluePr(>|t|)(Intercept)1.759660.07510123.434.543e-10x-0.048810.004697-10.391.116e-06summary(lm.exp)$r.sq=0.9153实例：分析百货商店销售额x与流通费率y之间的关系。结论：系数和模型都通过显著性检验，模型拟合效果好.

模型为:y=5.81e-0.049x4、幂函数曲线直线化幂函数曲线方程：y=axb(a>0)直线化后方程：y’=a’+bx’方法：两边取对数logy=loga+blogxlm.pow=lm(log(y)~log(x))#幂函数回归模型summary(lm.pow)$coef#幂函数回归模型系数

EstimateStd.ErrortvaluePr(>|t|)(Intercept)2.19070.0295174.234.806e-15log(x)-0.47240.01184-39.902.337e-12summary(lm.pow)$r.sq=0.9938实例：分析百货商店销售额x与流通费率y之间的关系。结论：系数和模型都通过显著性检验，模型拟合效果好.

模型为:y=8.942x-0.47245、不能精确直线化非线性回归模型可直线化曲线方程多项式曲线:y=a+bx+cx2+…+pxp对数曲线:y=a+blogx指数函数:y=aebx(a>0)幂函数:y=axb(a>0)双曲线函数:y=a+b/x经济增长模型:Y=A0emtLaKbCobb-Douglas生产函数:Q=aLb1Kb2+e消费函数:C=b1+b2Yb3+e不能精确直线化经典经济非线性模型非线性回归方程基本形式：（1）非线性回归模型参数计算Gauss-Newton迭代法求解:基本思想是使用泰勒级数展开式去近似地代替非线性回归模型，然后通过多次迭代，多次修正回归系数，使回归系数不断逼近非线性回归模型的最佳回归系数，最后使原模型的残差平方和达到最小。残差平方和：最小二乘法：在β=β1处Taylor展示：结论：对于给定初始值β=β1

，最小二乘原理，β2可求出Gauss-Newton迭代法求解对于给定初始值β=β1

，β2可求出结论：在给定初始值条件下，一直迭代,βn+1≈βn,即收敛，βn为最优解，残差平方和达到最小。该方法需要确定迭代初始值，通常取其线性模型结果作为其初始值β2估计值：估计拟线性模型：β3估计值：βn+1估计值：(S4=nls(y~a*(x^b),start=list(a=5,b=-0.1)))#幂函数回归Nonlinearregressionmodelmodel:y~a*(x^b)data:parent.frame()ab8.609-0.452residualsum-of-squares:0.164Numberofiterationstoconvergence:5Achievedconvergencetolerance:6.07e-07使用非