R语言方法总结要点.pdf

计算描述性统计量: 1、summary(): 例：summary(mtcarsvars) summary()函数提供了最小值、最大值、四分位数和数值型变量的均值，以及因子向量和逻 辑型向量的频数统计。 2、apply()函数或 sapply()函数 计算所选择的任意描述性统计量。mean、 sd、 var、 min、 max、 median、 length、 range 和 quantile。函数 fivenum() 可返回图基五数总括（Tukeys five -number summary，即最小值、 下四分位数、中位数、上四分位数和最大值）。 sapply() 例：mystats 2 就表明存在多重共线性问题 例： vif(fit) sqrt(vif(fit) 2 异常观测值 1、outlierTest() ：car 包 ： 求得最大标准化残差绝对值Bonferroni 调整后的p 值 例： library(car) outlierTest(fit) 高杠杆值点 1、hat.plot() : 观测点的帽子值大于帽子均值的2 或 3 倍，即可以认定为高杠杆值点 例： 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) 强影响点：Cook s D值大于 4/(n-k -1) ，则表明它是强影响点，其中 n 为样本量大小，k 是 预测变量数目。 例： cutoff - 4/(nrow(states) - length(fit$coefficients) - 2) plot(fit, which = 4, cook.levels = cutoff) abline(h = cutoff, lty = 2, col = “red“) 1、influencePlot()：car 包： 离群点、杠杆值和强影响点的信息整合到一幅图形中 例： influencePlot(fit, id.method = “identify“, main = “Influence Plot“, sub = “Circle size is proportial to Cook's Distance“) 纵坐标超过 +2 或小于2 的州可被认为是离群点， 水平轴超过0.2 或 0.3 的州有高杠杆值（通常为预测值的组合）。 圆圈大小与影响成比例，圆圈很大的点可能是对模型参数的估计造成的不成比例影响的强影 响点 变量变换 1、powerTransform():car 包 :函数通过 的最大似然估计来正态化变量x。 例： library(car) summary(powerTransform(states$Murder) 2、boxTidwell():car 包： 通过获得预测变量幂数的最大似然估计来改善线性关系 例： library(car) boxTidwell(Murder Population + Illiteracy, data = states) 模型比较 1、anova()：基础包： 比较两个嵌套模型的拟合优度 例： fit1 - lm(Murder Population + Illiteracy + Income + Frost, data = states) fit2 - lm(Murder Population + Illiteracy, data = states) anova(fit2, fit1) 2、AIC() ：AIC 值越小的模型 (可以不嵌套 )要优先选择，它说明模型用较少的参数获得了足 够的拟合度。 例： fit1 - lm(Murder Population + Illiteracy + Income + Frost, data = states) fit2 - lm(Murder Population + Illiteracy, data = states) AIC(fit1, fit2) 变量选择 1、stepAIC() ：MASS 包： 逐步回归模型 例： library(MASS) fit1 - lm(Murder Population + Illiteracy + Income + Frost, data = states) stepAIC(fit, direction = “backward“) 2、regsubsets() ：leaps包：全子集回归 例： library(leaps) leaps - regsubsets(Murder Population + Illiteracy + Income + Frost, data = states, nbest = 4) plot(leaps, scale = “adjr2“) 交叉验证 1、crossval() 函 数： bootstrap 包 ：实 现 k 重 交 叉 验 证 例： shrinkage - function(fit, k = 10) require(bootstrap) # define functions theta.fit - function(x, y) lsfit(x, y) theta.predict - function(fit, x) cbind(1, x) %*% fit$coef # matrix of predictors x - fit$model, 2:ncol(fit$model) # vector of predicted values y - fit$model, 1 results - crossval(x, y, theta.fit, theta.predict, ngroup = k) r2 - cor(y, fit$fitted.values)2 r2cv - cor(y, results$cv.fit)2 cat(“Original R-square =“, r2, “n“) cat(k, “Fold Cross-Validated R-square =“, r2cv, “n“) cat(“Change =“, r2 - r2cv, “n“) 2、shrinkage()：交叉验证;R 平方减少得越少，预测则越精确。 例： fit - lm(Murder Population + Income + Illiteracy + Frost, data = states) shrinkage(fit) 相对重要性 1、scale()：将数据标准化为均值为0、标准差为1 的数据集，这样用R 回归即可获得标准 化的回归系数。注意，scale()函数返回的是一个矩阵，而lm() 函数要求一个数据框 例： zstates - as.data.frame(scale(states) zfit - lm(Murder Population + Income + Illiteracy + Frost, data = zstates) coef(zfit) 2、relweights() ： 相对权重 例： relweights - function(fit, .) R - cor(fit$model) nvar - ncol(R) rxx - R2:nvar, 2:nvar rxy - R2:nvar, 1 svd - eigen(rxx) evec - svd$vectors ev - svd$values delta - diag(sqrt(ev) # correlations between original predictors and new orthogonal variables lambda - evec %*% delta %*% t(evec) lambdasq - lambda2 # regression coefficients of Y on orthogonal variables beta - solve(lambda) %*% rxy rsquare - colSums(beta2) rawwgt - lambdasq %*% beta2 import - (rawwgt/rsquare) * 100 lbls - names(fit$model2:nvar) rownames(import) - lbls colnames(import) - “Weights“ # plot results barplot(t(import), names.arg = lbls, ylab = “% of R-Square“, xlab = “Predictor Variables“, main = “Relative Importance of Predictor Variables“, sub = paste(“R-Square = “, round(rsquare, digits = 3), .) return(import) # using relweights() fit - lm(Murder Population + Illiteracy + Income + Frost, data = states) relweights(fit, col = “lightgrey“) 方差分析 1、aov() =lm() 单因素方差分析 2、plotmeans():绘制带置信区间的图形 例： library(multcomp) attach(cholesterol) table(trt) aggregate(response, by = list(trt), FUN = mean) aggregate(response, by = list(trt), FUN = sd) fit - aov(response trt) summary(fit) library(gplots) plotmeans(response trt, xlab = “Treatment“, ylab = “Response“, main = “Mean Plotnwith 95% CI“) detach(cholesterol) 多重比较 1、TukeyHSD() ：对各组均值差异的成对检验 例： TukeyHSD(fit) par(las = 2) par(mar = c(5, 8, 4, 2) plot(TukeyHSD(fit) par(opar) 2、glht()：multcomp 包： 多重均值比较 例： library(multcomp) par(mar = c(5, 4, 6, 2) tuk - glht(fit, linfct = mcp(trt = “Tukey“) plot(cld(tuk, level = 0.05), col = “lightgrey“) par(opar) 评估检验的假设条件 1、正态检验：library(car) qqPlot(lm(response trt, data = cholesterol), simulate = TRUE, main = “QQ Plot“, labels = FALSE) 2、方差齐性检验：bartlett.test(response trt, data = cholesterol) 3、检测离群点：outlierTest() car 包 library(car) outlierTest(fit) 单因素协方差分析 例： data(litter, package = “multcomp“) attach(litter) table(dose) aggregate(weight, by = list(dose), FUN = mean) fit - aov(weight gesttime + dose) summary(fit) 1、effects() ： effects包 ：计算调整的均值 例：library(effects) effect(“dose“, fit) 2、ancova() ：HH 包 ：绘制因变量、协变量和因子之间的关系图 例： library(HH) ancova(weight gesttime + dose, data = litter) 3、interaction.plot()：函数来展示双因素方差分析的交互效应 例： interaction.plot(dose, supp, len, type = “b“, col = c(“red“, “blue“), pch = c(16, 18), main = “Interaction between Dose and Supplement Type“) 4、plotmeans()：gplots包 ：展示交互效应 例： library(gplots) plotmeans(len interaction(supp, dose, sep = “ “), connect = list(c(1, 3, 5), c(2, 4, 6), col = c(“red“, “darkgreen“), main = “Interaction Plot with 95% CIs“, xlab = “Treatment and Dose Combination“) 5、interaction2wt() ：HH 包 ：可视化结果 例:library(HH) interaction2wt(len supp * dose) 6、colMeans():计算每列的平均值 7、nrow()/ncol : 计算数组额行数和列数 8、mahalanobis()：用协方差来计算两点之间距离的方法 稳健多元方差分析 Wilks.test() : 稳 健 单 因 素 MANOV A