第四课 回归分析.docx

上传人:b****6 文档编号:6461785 上传时间:2023-01-06 格式:DOCX 页数:23 大小:230.23KB
下载 相关 举报
第四课 回归分析.docx_第1页
第1页 / 共23页
第四课 回归分析.docx_第2页
第2页 / 共23页
第四课 回归分析.docx_第3页
第3页 / 共23页
第四课 回归分析.docx_第4页
第4页 / 共23页
第四课 回归分析.docx_第5页
第5页 / 共23页
点击查看更多>>
下载资源
资源描述

第四课 回归分析.docx

《第四课 回归分析.docx》由会员分享,可在线阅读,更多相关《第四课 回归分析.docx(23页珍藏版)》请在冰豆网上搜索。

第四课 回归分析.docx

第四课回归分析

回归分析

一、回归分析:

是建立因变量Y与自变量X之间关系的模型,一元线性回归使用一个自变量X,多元线性回归使用超过一个自变量,如:

.

1线性回归模型

进行n次独立观测:

,i=1,2,…,n

称为残差向量

2回归系数显著性检验:

判断某一自变量

的系数

是否为0,如果统计量(

表示)的

(通常取0.05),认为

3回归方程显著性检验:

检查是否可用线性方程来处理数据。

即检验方程系数

是否全为0,如果统计量(

表示)的

,认为可以用线性方程来处理问题。

4相关性检验:

检验因变量

与自变量

相关性程度。

用相关系数的平方

来衡量,

的值接近1表示相关,接近0表示不相关。

二、线性回归模型的计算

1、lm()函数:

完成多元线性回归系数的估计、回归系数的检验、回归方程的检验等工作,返回值称为拟合结果的对象,本质上是一个具有类属性值lm的列表.

格式:

lm(formula,data,subset,weights,na.action,method=”qr”,model=TRUE,x=FALSE,y=FALSE,qr=TRUE,singular.ok=TRUE,contrasts=NULL,offset,…)

formula为模型公式,形如y~1+x1+x2的形式,表示常数项、

的系数和

的系数,去掉

公式中的1,其意义不变。

拟合成齐次线性模型:

y~0+x1+x2或y~-1+x1+x2的形式;

data为数据框,由样本数据构成;

subset为可选项,表示所使用的样本子集;

weights为可选向量,表示对应样本的权重;

na.action为函数,表示当数据中出现缺失数据(NA)的处理方法;

method为估计回归系数的计算方法,默认值为“qr”;

model、x、y、qr为逻辑变量,如果取TRUE,函数的返回值将给出模型的框架、模型矩阵、

响应变量,及QR分解;

singular.ok为逻辑变量,取FALSE表示奇异值拟合是错误的;

contrasts为可选列表;

offset为NULL,或者是数值向量。

…为附加参数。

2、summary()函数

lm()函数的返回值称为拟合结果的对象,是一个有类属性值lm的列表,有model,coefficients,residuals(残差)等成员,为了获得更多信息,summary()常与lm()函数一起使用,格式:

summary(object,correlation=FALSE,symbolic.cor=FALSE,…)

参数意义:

object为lm()函数生成的对象,correlation为逻辑变量,取TRUE表示给出估计参数的相关矩

阵,symbolic.cor为逻辑变量,取TRUE表示用符号形式给出估计参数的相关矩阵,此参数只有当correlation=TRUE时才有效。

例1年龄相等情况下,建立血压的收缩压

与体重

(kg),年龄

(岁数)有关,建立

的线性回归方程。

blood<-data.frame(

X1=c(76.0,91.5,85.5,82.5,79.0,80.5,74.5,79.0,85.0,76.5,82.0,95.0,92.5),

X2=c(50,20,20,30,30,50,60,50,40,55,40,40,20),

Y=c(120,141,124,126,117,125,123,125,132,123,132,155,147)

lm.sol<-lm(Y~1+X1+X2,data=blood)

summary(lm.sol)

运行结果:

Call:

lm(formula=Y~1+X1+X2,data=blood)###调出函数使用的模型

Residuals:

(残差)###列出残差

Min1QMedian3QMax

-4.0404-1.01830.46400.69084.3274

Coefficients:

##列出估计值、估计值的标准差、t统计量、对应t统计量的P值

EstimateStd.ErrortvaluePr(>|t|)

(Intercept)-62.9633616.99976-3.7040.004083**

X12.136560.1753412.1852.53e-07***

X20.400220.083214.8100.000713***

---

极为显著***,高度显著**,显著*

Signif.codes:

0‘***’0.001‘**’0.01‘*’0.05‘.’0.1‘’1

Residualstandarderror:

2.854on10degreesoffreedom###残差的标准差

MultipleR-squared:

0.9461,AdjustedR-squared:

0.9354##相关系数的平方,修正相关系数的平方

F-statistic:

87.84on2and10DF,p-value:

4.531e-07###F统计量,F统计量对应的P值

回归方程与回归系数的检验都是显著的,回归方程为:

3预测区间与置信区间

如果回归效果显著,就可以利用回归方程进行预测,即对给定的回归自变量值,预测因变量可能的取值范围,是一个区间估计问题。

假设

=

回归方程的真实值(未知):

相应的估计值:

predict()函数计算

的估计值、预测区间和置信区间。

格式:

predict(object,newdata,se.fit=FALSE,scale=NULL,df=Inf,

interval=c(“none”,”confidence”,”prediction”),level=0.95,type=c(“response”,”terms”),terms=NULL,

na.action=na.pass,pred.var=res.var/weights,weights=1,…)

参数意义:

object为lm()函数得到的对象;

newdata为数据框,由预测点构成,如取默认值,将计算已知数据的回归值;

se.fit为逻辑变量,取TRUE表示输出预测值的标准差、自由度和残差尺度信息;

scale为计算标准差的尺度参数;

df为尺度参数的自由度;

interval为计算的区间类型,取“none”(默认值)表示不计算,取”confidence”表示计算置信区

间,取”prediction”表示计算预测区间;

level为置信水平,默认值为0.95,在计算置信区间和预测区间时用到;

type为预测类型,或者取“response”(默认值),或者取”terms”,当type=”terms”时,为全部

选择项;

terms为选择项,只能取NULL,1,2等;

na.action为函数,表示处理newdata中有缺失数据(NA)时的处理方法,默认预测为NA.pred.var

预测区间的方差;

weights为数值向量或单侧模型公式,用于预测时方差的权;

…附加参数。

例2继续上面的例1,设

,求

的估计值、预测区间和置信区间(置信水平为0.95)。

blood<-data.frame(

X1=c(76.0,91.5,85.5,82.5,79.0,80.5,74.5,79.0,85.0,76.5,82.0,95.0,92.5),

X2=c(50,20,20,30,30,50,60,50,40,55,40,40,20),

Y=c(120,141,124,126,117,125,123,125,132,123,132,155,147)

lm.sol<-lm(Y~1+X1+X2,data=blood)

newdata<-data.frame(X1=80,X2=40)

predict(lm.sol,newdata,interval="prediction")

predict(lm.sol,newdata,interval="confidence")

例分析合金强度

与合金中含碳量

之间的关系:

(它们之间的一组数据如下)

(1)完成一元线性回归的计算;

(2)计算自变量

在区间[0.10,0.23]内的回归方程的预测估计值、预测区间和置信区间(取

),并将数据点、预测估计曲线、预测区间曲线和置信区间曲线画在同一张图上。

##回归估计

x<-c(0.10,0.11,0.12,0.13,0.14,0.15,0.16,0.17,0.18,0.20,0.21,0.23)

y<-c(42.0,43.5,45.0,45.5,45.0,47.5,49.0,53.0,50.0,55.0,55.0,60.0)

lm.sol<-lm(y~1+x)

summary(lm.sol)

##计算预测值,并绘图

new<-data.frame(x=seq(0.10,0.24,by=0.01))##取x的范围值

pp<-predict(lm.sol,new,interval="prediction")

pc<-predict(lm.sol,new,interval="confidence")

par(mai=c(0.8,0.8,0.2,0.2))

matplot(new$x,cbind(pp,pc[,-1]),type="l",######绘出矩形各列pp和pc[,-1]共5列

xlab="X",ylab="Y",lty=c(1,5,5,2,2),

col=c("blue","red","red","brown","brown"),lwd=2)

points(x,y,cex=1.4,pch=21,col="red",bg="orange")###真实值x,y

legend(0.1,63,

c("Points","Fitted","Prediction","Confidence"),

###如果是字符串,用c(“”),而不用expr1数学表达式

pch=c(19,NA,NA,NA),lty=c(NA,1,5,2),

col=c("orange","blue","red","brown"))

savePlot("predict",type="eps")

三、回归诊断

为什么要做回归诊断?

因为有些数值例子得到的回归方程尽管通过了t检验和F检验,但做线性回归方程还是有问题的。

如Anscombe在1973年构造的数据集给出4组数据拟合简单的线性模型:

,得到的结果:

[[1]]

EstimateStd.ErrortvaluePr(>|t|)

(Intercept)3.00009091.12474682.6673480.025734051

x10.50009090.11790554.2414550.002169629

[[2]]

EstimateStd.ErrortvaluePr(>|t|)

(Intercept)3.0009091.12530242.6667580.025758941

x20.5000000.11796374.2385900.002178816

[[3]]

EstimateStd.ErrortvaluePr(>|t|)

(Intercept)3.00245451.12448122.6700800.025619109

x30.49972730.11787774.2393720.002176305

[[4]]

EstimateStd.ErrortvaluePr(>|t|)

(Intercept)3.00172731.12392112.6707630.025590425

x40.49990910.11781894.2430280.002164602

并且进一步检验4组数据的

、F值检验全合格,但绘出图分析4组数据做线性回归模型不是全合格:

只有第一组数据适合做线性回归模型。

1残差检验

检验模型误差是否满足正态性和方差齐性。

1)用残差图检验,残差为纵坐标,拟合值或数据观测序号或数据观测时间做横坐标,做出的散点图称为残差图。

(a)(b)(c)

正常情况异方差情况非线性情况

有异常值点的情况

检查1)有无有异常值点;2)画出残差的Q-Q散点图,若这些点位于一条直线上,则残差服从正态分布,否则不然。

注意:

上面图中只有(a)图是正常的,(b)图是喇叭口,回归值的大小与残差的波动有关,方差的假设有问题,(c)图表明线性模型不合适。

residuals()(或resid()函数)计算回归模型的残差,rstandard()计算回归模型的标准化残差,rstudent()计算回归模型的化学生化残差(删除第i个样本数据后得到的标准化残差)

格式:

residuals(object,...)resid(object,...)

rstandard(model,infl=lm.influence(model,do.coef=FALSE),

sd=sqrt(deviance(model)/df.residual(model)),...)

rstudent(model,infl=lm.influence(model,do.coef=FALSE),res=infl$wt.res,...)

参数object或model为lm生成的对象;infl为lm.influence返回值得到的影响结构;sd为模型的标准差;res为模型残差。

plot()函数绘出各种残差的散点图,格式:

plot(x,which=c(1:

3,5),...)),参数x由lm生成的对象,which是开关变量,其值1~6,默认子集(1,2,3,5),即绘出1、2、3和5号散点图,该函数可以绘出6张诊断图,第一张残差与预测值散点图,第二张残差的正态Q-Q图,第三张标准差的平方根与预测值的散点图等。

例20~60成年女性的血压与年龄之间的回归关系,并作残差分析

##

(1)回归

rt<-read.table("blood.dat",header=TRUE)

lm.sol<-lm(Y~X,data=rt);lm.sol

summary(lm.sol)

##

(2)残差图

pre<-fitted.values(lm.sol)###拟合值

res<-residuals(lm.sol)###残差

rst<-rstandard(lm.sol)####标准差

par(mai=c(0.9,0.9,0.2,0.1))

plot(pre,res,xlab="FittedValues",ylab="Residuals")

savePlot("resid-1",type="eps")

plot(pre,rst,xlab="FittedValues",ylab="StandardizedResiduals")

savePlot("resid-2",type="eps")

##(3)对残差作回归

rt<-read.table("blood.dat",header=TRUE)

rt$res<-res###########给rt数据框增加一列res

lm.res<-lm(abs(res)~X,data=rt);################利用残差的绝对值与自变量做回归

lm.res

summary(lm.res)

##(4)修正回归方程

s<-lm.res$coefficients[1]+lm.res$coefficients[2]*rt$X#####计算残差的标准差

lm.weg<-lm(Y~X,data=rt,weights=1/s^2);#####用标准差平方的倒数作权重

################权重的选取可减少非齐性方差带来的影响

lm.weg

summary(lm.weg)

2影响分析

影响分析就是检查是否有对估计有异常大数据(高杠杆点)情况。

影响分析的方法有:

DFFITS准则、Cook距离、COVRATIO准则、帽子值和帽子矩阵。

相关的R函数有:

influence.measures(),dffits(),dfbeta(),dfbetas(),cooks.distance(),covratio(),hatvalues(),hat();

格式:

influence.measures(model)

dffits(model,infl=,res=)

dfbeta(model,infl=lm.influence(model,do.coef=TRUE),...)

dfbetas(model,infl=lm.influence(model,do.coef=TRUE),...)

covratio(model,infl=lm.influence(model,do.coef=FALSE),res=weighted.residuals(model))

cooks.distance(model,infl=lm.influence(model,do.coef=FALSE),res=weighted.residuals(model),sd=sqrt(deviance(model)/df.residual(model)),hat=infl$hat,...)

hatvalues(model,infl=lm.influence(model,do.coef=FALSE),...)

hat(x,intercept=TRUE)

参数model为lm生成的对象,infl为lm.influence返回值得到的影响结构,res为模型残差,x为设计矩阵。

例通过智力测试数据建立智力随年龄变化的关系(在exam0606.R文件中)。

第一步:

运行,计算回归系数,并做回归系数与回归方程检验

##1.回归分析

intellect<-data.frame(

x=c(15,26,10,9,15,20,18,11,8,20,7,9,10,11,11,10,12,42,17,11,10),######年龄(月)

y=c(95,71,83,91,102,87,93,100,104,94,113,96,83,84,102,100,105,57,121,86,100)

##智力

lm.sol<-lm(y~1+x,data=intellect)

summary(lm.sol)

第二步做回归诊断

##2.回归诊断

op<-par(mfrow=c(2,2),mar=0.4+c(4,4,1,1),oma=c(0,0,2,0))

plot(lm.sol,1:

4)#####通过绘图分析有无异常值点

par(op)

savePlot("diagnoses1",type="eps")

influence.measures(lm.sol)#####通过函数计算分析有无异常值点,下面是运行结果:

dfb.1_dfb.xdffitcov.rcook.dhatinf

10.016640.003280.041271.1668.97e-040.0479

20.18862-0.33480-0.402521.1978.15e-020.1545

.....................

180.83112-1.11275-1.155782.9596.78e-010.6516*

190.143480.273170.853740.3962.23e-010.0531*

20-0.207610.10544-0.263851.0433.45e-020.0567

210.02791-0.016220.032981.1875.74e-040.0628

influence.measures()函数得到回归诊断结果共7列,1,2列(dfbetas)对应常数和变量x值,3列是DFFITS准则值,4列是COVRATIO准则值,5列是Cook距离,6列高杠杆值,7列是影响点记号。

(18,19号样本是强影响点记号)。

第三步处理强影响点

##3.回归修正

n<-length(intellect$x)

weights<-rep(1,n);weights[18]<-0.5

lm.correct<-lm(y~1+x,data=intellect,subset=-19,weights=weights)

summary(lm.correct)

第四步验证

##4.验证####在脚本里做

intellect<-data.frame(

x=c(15,26,10,9,15,20,18,11,8,20,7,9,10,11,11,10,12,42,17,11,10),######年龄(月)

y=c(95,71,83,91,102,87,93,100,104,94,113,96,83,84,102,100,105,57,121,86,100)

##智力

lm.sol<-lm(y~1+x,data=intellect)

attach(intellect)######把intellect调入内存,可操作它的元素

par(mai=c(0.8,0.8,0.2,0.2))

plot(x,y,cex=1.2,pch=21,col="red",bg="orange")########真值

abline(lm.sol,col="blue",lwd=2)###########修正前,lm.sol是拟合值

text(x[c(19,18)],y[c(19,18)],label=c("19","18"),adj=c(1.5,0.3))

###adj()调整符号坐标位置

detach()

abline(lm.correct,col="red",lwd=2,lty=5)##

展开阅读全文
相关资源
猜你喜欢
相关搜索

当前位置:首页 > 高中教育 > 数学

copyright@ 2008-2022 冰豆网网站版权所有

经营许可证编号:鄂ICP备2022015515号-1