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 c 75 64 47 4 66 9 62 2 62 2 58 7 63 5 sort x order x 中位数 中位数描述数据中心位置的数字特征 大体上比中位数大或小的数据个数为整个数据的一半 对于对称分布的数据 均值与中位数比较接近 对于偏态分布的数据 均值与中位数不同 中位数的又一显著特点是不受异常值的影响 具有稳健性 因此它是数据分析中相当重要的统计量 在R软件中 函数median 给观测量的中位数 如x c 75 64 47 4 66 9 62 2 62 2 58 7 63 5 median x median x na rm TRUE 若数据中有缺失值 百分位数 百分位数 percentile 是中位数的推广 将数据按从小到大的排列后 0 p 1 它的p分位点定义为 在R软件中 quantile 函数计算观测量的百分位数 如w c 75 0 64 0 47 4 66 9 62 2 62 2 58 7 63 5 66 6 64 0 57 0 69 0 56 9 50 0 72 0 quantile w 一般用法 quantile x probs seq 0 1 0 25 na rm FALSE 分散程度的度量 表示数据分散 或变异 程度的特征量有方差 标准差 极差 四分位极差 变异系数和标准误等 在R软件中 用var 和sd 计算方差 标准差 var x na rm FALSE sd x na rm FALSE 变异系数 平方和 对于变异系数 校正平方和 未校正平方和等指标 需要编写简单的程序 变异系数CV计算 cv 100 sd x mean x cv校正平方和CSS css sum x mean x 2 css未校正平方和USS uss sum x 2 uss 极差与标准误 样本极差 记为R 的计算 R max x min x 样本上 下四分位数之差称为四分位差 或半极差 记为R1 它也是度量样本分散性的重要数字特征 特别对于具有异常值的数据 它作为分散性具有稳健性 因此在稳健性数据分析中具有重要作用 半极差计算 R1 quantile x 0 75 quantile x 0 25 样本标准误 记为sm 定义为s sqrt n 样本标准误计算 sm sd x sqrt length x 分布形状的度量 偏度系数Kurtosis是刻划数据的对称性指标 关于均值对称的数据其偏度系数为0 右侧更分散的数据偏度系数为正 左侧更分散的数据偏度系数为负 当数据的总体分布为正态分布时 峰度系数Skewness近似为0 当峰度系数为正时 两侧极端数据较多 当峰度系数为负时 两侧极端数据较少 偏度系数Skewness 样本峰度系数sk计算程序n length x m mean x s sd x sk n n 1 n 2 sum x m 3 s 3计算公式 峰度系数Kurtosis计算 样本峰度系数ku计算程序n length xm mean x s sd x ku n n 1 n 1 n 2 n 3 sum x m 4 s 4 3 n 1 2 n 2 n 3 计算公式 相关分析 R软件采用用cov 函数计算协方差或协方差阵 用cor 函数计算相关矩阵 相关系数 函数cov 和cor 的使用格式为 cov x y NULL use all obs method c pearson kendall spearman cor x y NULL use all obs method c pearson kendall spearman 其中x是数值型向量 矩阵或数据框 y是空值 NULL 缺省值 向量 矩阵或数据框 但需要与x的维数相一致 与cov和cor有关的函数还有 cov wt 计算加权协方差 加权协方差矩阵 cor test 计算相关性检验 相关分析示例 例为了解某种橡胶的性能 今抽取10个样品 每个测量三项指标 硬度 变形和弹性 rubber txt 试计算样本均值 样本协方差阵和样本相关矩阵 并用Pearson相关性检验确认变量X1 X2 X3是否相关 rubber read table d rubber txt mean rubber cov rubber cor rubber cor test X1 X2 data rubber cor test X1 X3 data rubber cor test X2 X3 data rubber 回归分析 案例 根据经验 在人的身高相等的情况下 血压的收缩压Y与体重X1 千克 年龄X2 岁数 有关 现收集了13个男子的数据 见表 试建立Y关于X1 X2的线性回归方程 估计出Y b0 b1X1 b2X2F检验 H0 b1 b2 0 T检验 H0 bj 0j 0 1 2 求解程序 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 X1 X2 data blood 进行回归分析summary lm sol 汇总分析结果Y 62 96 2 136X1 0 4002X2 预测 X 80 40 时 相应Y的概率为0 95的预测区间 new data frame X1 c 80 75 X2 c 40 38 lm pred predict lm sol new interval prediction level 0 95 lm pred 回归分析的汇总结果 Call lm formula Y X1 X2 data blood Residuals Min1QMedian3QMax 4 0404 1 01830 46400 69084 3274Coefficients 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 1Residualstandarderror 2 854on10degreesoffreedomMultipleR squared 0 9461 AdjustedR squared 0 9354F statistic 87 84on2and10DF p value 4 531e 07预测结果如下 fitlwrupr1123 9699117 2889130 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位小数点的格式展示各变量的方差膨胀因子情况如下 X1X21 961 96可以看到所有参数估计的VIFj 1 1 Rj2 值都远远小于10 并且接近1 因此这里我们不用担心多重共线性的问题 二项选择模型 当我们考虑多个连续解释变量对某个取0 1值的响应变量的影响时 R中常用probit或logit回归来分析 probit 1 P Y 1 0 Xlogit logit P Y 1 log P Y 1 1 P Y 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 Ymat x family binomial data norell logit回归 glm sol glm Ymat x family binomial link probit data norell summary glm sol 预测 pre predict glm sol data frame x 3 5 p exp pre 1 exp pre pd seq 0 5 len 100 pre predict glm sol data frame x d p exp pre 1 exp pre norell y norell success norell nplot 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 Y X1 X2 X3 family binomial data life summary glm sol 定序回归模型 当我们考察多个连续解释变量对某个取有序离散值的响应变量的影响时 适用定序回归 R中常用MASS包中polr函数来实现 Probit定序回归 1 P Y k k kXLogit定序回归 logit P Y 1 k kXpolr formula data weights start subset na action contrasts NULL Hess FALSE model TRUE method c logistic probit cloglog cauchit 定序回归示例