R语言实战笔记–第八章 OLS回归分析

R语言实战笔记–第八章 OLS回归分析

标签（空格分隔）： R语言 回归分析

myfit<-lm(I(Y^(a))~x+I(x^2)+I(log(x))+var+…+[-1],data=dataframe) 

Call: #拟合公式 lm(formula = Murder ~ Population + Illiteracy + Income + Frost, data = states) Residuals: #残差的图基五数 Min 1Q Median 3Q Max -4.7960 -1.6495 -0.0811 1.4815 7.6210 Coefficients: #回归系数及相关 Estimate Std. Error t value Pr(>|t|) #系数估计值 标准误 t统计量 P值 (Intercept) 1.235e+00 3.866e+00 0.319 0.7510 Population 2.237e-04 9.052e-05 2.471 0.0173 * Illiteracy 4.143e+00 8.744e-01 4.738 2.19e-05 * Income 6.442e-05 6.837e-04 0.094 0.9253 Frost 5.813e-04 1.005e-02 0.058 0.9541 --- Signif. codes: 0 ‘*’ 0.001 ‘’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 #上面一行标明在上面的回归系数后面跟着的*的含义，0-0.001之间为'*',0.001-0.001之间为''……，直观表示P值的范围 Residual standard error: 2.535 on 45 degrees of freedom #残差标准误 Multiple R-squared: 0.567, Adjusted R-squared: 0.5285 #决定系数R平方及修正R平方 F-statistic: 14.73 on 4 and 45 DF, p-value: 9.133e-08 #F统计量 

residplot <- function(fit, nbreaks=10) { #标准化残差z z <- rstudent(fit) #画直方图 hist(z, breaks=nbreaks, freq=FALSE,xlab="Studentized Residual",main="Distribution of Errors") #画轴须图 rug(jitter(z), col="brown") #画正态曲线图（以标准化残差的均值及标准差画正常的正态曲线） curve(dnorm(x, mean=mean(z), sd=sd(z)),add=TRUE, col="blue", lwd=2) #拟合核密度曲线 lines(density(z)$x, density(z)$y,col="red", lwd=2, lty=2) #右上角加入图例 legend("topright",legend = c( "Normal Curve", "Kernel Density Curve"),lty=1:2, col=c("blue","red"), cex=.7) } residplot(fit，nbreaks=n) #调用函数，fit为回归模型lm(),nbreaks为直方图分组数 

　　误差的独立性，car包主要提供的是Durbin-Watson检验（函数durbinWatsonTest()，一种检测样本是否存在自相关关系的方法，其原假设为独立，即不相关，当p
<显著水平时，显著拒绝原假设）。其它书上介绍了误差自相关时的解决方案：误差不独立，可以看作是自相关，也叫序列相关，一般来说，有三种原因导致其产生，1）模型的数学形式错误，可以通过利用残差对较高次幂的解释变量进行回归，然后做dw检验，如果序列相关性消失，则为数学形式错误，2）经济变量的惯性，比如时间序列，3）回归模型中略去了带有自相关的重要解释变量，可以通过利用残差对可能影响但未列入模型的解释变量进行回归，并做显著性检验，若显著，则则列入模型。若误差存在自相关，并且排除了上述因素1和因素3后，可以使用orcutt包的cochrane.orcutt(lmobject)函数来进行求解其回归系数。
　　

　　线性，成份残差图/偏残差图，使用crPlot来画出残差与各成份（即变量）的散点图，实线与虚线拟合度越高越好，即实线最好能与虚线重合。
　　
　　同方差性，ncvTest()函数将生成一个卡方检验量，零假设为方差不变，即P值小于显著水平时，显著拒绝方差不变的零假设。另外，car包还给出了一个Suggested power transformation（函数为spreadLevelPlot），建议幂次变换值 $\lambda$，使 $Y^{\lambda}$得到一个较为平稳的误差方差，可若建议为0，则使用对数变换，若为0.5，则使用 $\sqrt{y}$代替y。建议以0.5为一进位，spreadLevelPlot给出了一张图，若图形的线非水平，，则可以按照建议值进行y值的幂次变换。另外，异方差的克服也可以通过Glejser方法进行检验并得到残差绝对值与变量的回归函数关系式，然后对现有回归方程中所有变量均除这个关系式结果得到新解释变量观测值，利用新观测值进行回归拟合即可。
　　
　　综合验证，gvlma包给出了一个综合验证函数gvlma()上述假设的综合验证结果（acceptable/NOT satisfied）若出现不满意，则需要逐个排除。
　　
　　多重共线性，当自变量之间存在极高相关性时，将会出现整体的F检验非常显著，但各个回归系数却不显著的情况，要检测此值，可以使用方差膨胀因子VIF进行检测，一般来说 $\sqrt{vif}>2$就表明存在多重共线性问题。多重共线性的克服可以通过逐步回归方法（setp()函数）对变量的筛选来进行，但逐步回归方法会对变量进行删除，会有一定程度的失真，另一种方法叫 岭回归，它可以在不删除变量的基础上，选择适当的岭参数对回归系数进行修正，使岭估计相对稳定，函数为MASS包中的lm.ridge()函数，相对lm()函数多出了一个岭参数lambda的设置，岭迹法就是把lambda参数设置为序列，并画出来，选择其稳定点的具体参数进行代入计算，选择稳定点可以使用select()函数进行更精确的选择（有HKB、L-W和GCV三种方法，具体先哪个可以分别测试），方法如下：

library(MASS) lm.ridge(Y~X,data=data)#未设置岭参数时，岭回归与lm结果一样 plot(lm.ridge(Y~X,data=data,lambda=seq(0,1,0.0001)))#画岭迹图 select(lm.ridge(Y~X,data=data,lambda=seq(0,1,0.0001)))#选岭参数 lm.ridge(Y~X,data=data,lambda=0.0038)#根据选择的具体参数代入计算 

hat.plot <- function(fit){ p <- length(coefficients(fit)) n <- length(fitted(fit)) plot(hatvalues(fit), main = "Index Plot of Hat Values") abline(h = c(2, 3) * p/n, col = "red", lty = 2) identify(1:n, hatvalues(fit), names(hatvalues(fit))) } #调用 hat.plot(fit)