R语言常用统计方法实现.ppt

常用统计方法用R实现 描述性统计 位置的度量: 均值、顺序统计量、中位数、百分位数。 均值计算： 若x是向量、矩阵，则mean(x)返回其全部元素均值。 若要返回数组某一维的均值：apply(x,dim,mean); dim=1计算行 均值，dim=2计算列均值。 若x是数据框，则mean(x)返回各列的均值 Mean的一般用法： mean(x,trim=0,na.rm=FALSE) trim指定去掉x两端数的比例；na.rm=TRUE允许有缺失值。 类似有sum(x)函数可求x的和。 顺序统计量 将n个数据(观测值)按从小到大的顺序排列 后,称其为顺序统计量. 函数sort(x)给出了样本x的顺序统计量 order ( )给出排序后的下标 rank( )给出了样本x的秩次统计量 x|t|) (Intercept) -62.96336 16.99976 -3.704 0.004083 * X1 2.13656 0.17534 12.185 2.53e-07 * X2 0.40022 0.08321 4.810 0.000713 * - Signif. codes: 0 * 0.001 * 0.01 * 0.05 . 0.1 1 Residual standard error: 2.854 on 10 degrees of freedom Multiple R-squared: 0.9461, Adjusted R-squared: 0.9354 F-statistic: 87.84 on 2 and 10 DF, p-value: 4.531e-07 预测结果如下： fit lwr upr 1 123.9699 117.2889 130.6509 回归诊断 par(mfrow=c(2,2) #设置画图为2x2的格式 plot(lm.sol,which=c(1:4) #模型检验4张图，包括残差图、QQ图和Cook距离图 数据太少，上面诊断结果并不理想。 library(car) #载入程序包Car,vif()函数在其内 round(vif(lm.sol),2) #计算模型的方差膨胀因子，用2位小数点的格式展示 各变量的方差膨胀因子情况如下： X1 X2 1.96 1.96 可以看到所有参数估计的VIFj=1/(1-Rj2)值都远远小于10， 并且接近1。因此这里我们不用担心多重共线性的问题。 二项选择模型 当我们考虑多个连续解释变量对某个取0-1值的响应变量的影响时， R中常用probit或logit回归来分析。 probit: -1(PY=1)=0+X logit: logit(PY=1)=log(PY=1/(1-PY=1)=0+X 对二项选择的probit/logit回归,R软件可用glm()处理. fm-glm(formula,family=binomial(link=probit),data=data.frame) fm-glm(formula,family=binomial(link=logit),data=data.frame) 在用glm()函数处理二项选择模型时，公式中响应变量y的输入形式 有两种：1、y中第一列为对应自变量的响应次数，第2列是不响应的 次数；2、y是只由0、1构成的向量，分别表示对应自变量取值是不 响应还是相应。 二项选择案例1 研究小电流对农场动物的影响.选择了7头牛,6种电 击强度0,1,3,4,5毫安.给出每种电击强度70次试验 中牛发生响应的总次数.试分析电击对牛的影响。 案例1的程序 norell-data.frame(x=0:5,n=rep(70,6),success=c(0,9,21,47,60,63) norell$Ymat-cbind(norell$success,norell$n-norell$success) glm.sol-glm(Ymatx,family=binomial,data=norell) # logit回归 #glm.sol-glm(Ymatx,family=binomial(link=probit),data=norell) summary(glm.sol) 预测: pre-predict(glm.sol,data.frame(x=3.5) p-exp(pre)/(1+exp(pre);p d-seq(0,5,len=100) pre-predict(glm.sol,data.frame(x=d) p-exp(pre)/(1+exp(pre) norell$y-norell$success/norell$n plot(norell$x,norell$y);lines(d,p) 二项选择案例2 50位急性林巴细胞性白血病病人，在入院治疗时取得了外 辕血中的细胞数X1、林巴结浸润等级X2(分为0,1,2,3级);出 院后有无巩固治疗X3(“1”表示有巩固治疗，“0”表示无巩固 治疗).并取得病人的生存时间,Y=0表示生存时间在1年以内 ，Y=1表示生存时间在1年或1年以上.试分析病人生存时间 长短的概率与X1,X2,X3的关系. 案例2的程序 life-data.frame( X1=c(2.5,173,119,10,502,4,14.4,2,40,6.6, 21.4,2.8,2.5,6,3.5,62.2,10.8,21.6,2,3.4, 5.1,2.4,1.7,1.1,12.8,1.2,3.5,39.7,62.4,2.4, 34.7,28.4,0.9,30.6,5.8,6.1,2.7,4.7,128,35, 2,8.5,2,2,4.3,244.8,4,5.1,32,1.4), X2=rep(c(0,2,0,2,0,2,0,2,0,2,0,2,0,2,0,2,0,2,0,2,0,2,0), c(1,4,2,2,1,1,8,1,5,1,5,1,1,1,2,1,1,1,3,1,2,1,4), X3=rep(c(0,1,0,1,0,1,0,1,0,1,0,1,0,1,0,1,0,1), c(6,1,3,1,3,1,1,5,1,3,7,1,1,3,1,1,2,9), Y=rep(c(0,1,0,1),c(15,10,15,10) glm.sol-glm(YX1+X2+X3,family=binomial,data=life) summary(glm.sol) 定序回归模型 当我们考察多个连续解释变量对某个取有序离散值的响应 变量的影响时，适用定序回归，R中常用MASS包中polr函 数来实现。 Probit定序回归: -1(PYk)=k-kX Logit定序回归: logit(PY=1)= k-kX polr(formula, data