描述性统计分析汇总

简单描述性统计分析

通过summary()计算描述性统计量

对于基础安装,你可以使用summary()函数来获取描述性统计量。代码清单7-1展示了一个 示例。

#通过summary()计算描述性统计量
myvars <- c("mpg", "hp", "wt")
summary(mtcars[myvars])
##       mpg              hp              wt       
##  Min.   :10.40   Min.   : 52.0   Min.   :1.513  
##  1st Qu.:15.43   1st Qu.: 96.5   1st Qu.:2.581  
##  Median :19.20   Median :123.0   Median :3.325  
##  Mean   :20.09   Mean   :146.7   Mean   :3.217  
##  3rd Qu.:22.80   3rd Qu.:180.0   3rd Qu.:3.610  
##  Max.   :33.90   Max.   :335.0   Max.   :5.424

通过sapply()计算描述性统计量

使用apply()函数或sapply()函数计算所选择的任意描述性统计量。对于sapply()函数,其使用格式为: sapply(x, FUN, options) 其中的x是你的数据框(或矩阵),FUN为一个任意的函数。如果指定了options,它们将被传递 给FUN。你可以在这里插入的典型函数有mean()、sd()、var()、min()、max()、median()、length()、range()和quantile()。函数fivenum()可返回图基五数总括(Tukey’s five-number summary,即最小值、下四分位数、中位数、上四分位数和最大值)。令人惊讶的是,基础安装并没有提供偏度和峰度的计算函数,不过你可以自行添加。代码清 单7-2中的示例计算了若干描述性统计量,其中包括偏度和峰度。常用于定制化的描述性统计分析

#通过sapply()计算描述性统计量
mystats <- function(x, na.omit=FALSE){ 
 if (na.omit) 
 x <- x[!is.na(x)] 
 m <- mean(x) 
 n <- length(x) 
 s <- sd(x) 
 skew <- sum((x-m)^3/s^3)/n 
 kurt <- sum((x-m)^4/s^4)/n - 3 
 return(c(n=n, mean=m, stdev=s, skew=skew, kurtosis=kurt)) 
 } 
myvars <- c("mpg", "hp", "wt") 
sapply(mtcars[myvars], mystats)
##                mpg          hp          wt
## n        32.000000  32.0000000 32.00000000
## mean     20.090625 146.6875000  3.21725000
## stdev     6.026948  68.5628685  0.97845744
## skew      0.610655   0.7260237  0.42314646
## kurtosis -0.372766  -0.1355511 -0.02271075

通过Hmisc包中的describe()函数计算描述性统计量

Hmisc包中的describe()函数可返回变量和观测的数量、缺失值和唯一值的数目、平均值、 分位数,以及五个最大的值和五个最小的值。

#通过Hmisc包中的describe()函数计算描述性统计量
library(Hmisc) 
## Warning: package 'Hmisc' was built under R version 3.5.3

## Loading required package: lattice

## Loading required package: survival

## Loading required package: Formula

## Loading required package: ggplot2

## Warning: package 'ggplot2' was built under R version 3.5.3

## 
## Attaching package: 'Hmisc'

## The following objects are masked from 'package:base':
## 
##     format.pval, units
myvars <- c("mpg", "hp", "wt") 
describe(mtcars[myvars])
## mtcars[myvars] 
## 
##  3  Variables      32  Observations
## ---------------------------------------------------------------------------
## mpg 
##        n  missing distinct     Info     Mean      Gmd      .05      .10 
##       32        0       25    0.999    20.09    6.796    12.00    14.34 
##      .25      .50      .75      .90      .95 
##    15.43    19.20    22.80    30.09    31.30 
## 
## lowest : 10.4 13.3 14.3 14.7 15.0, highest: 26.0 27.3 30.4 32.4 33.9
## ---------------------------------------------------------------------------
## hp 
##        n  missing distinct     Info     Mean      Gmd      .05      .10 
##       32        0       22    0.997    146.7    77.04    63.65    66.00 
##      .25      .50      .75      .90      .95 
##    96.50   123.00   180.00   243.50   253.55 
## 
## lowest :  52  62  65  66  91, highest: 215 230 245 264 335
## ---------------------------------------------------------------------------
## wt 
##        n  missing distinct     Info     Mean      Gmd      .05      .10 
##       32        0       29    0.999    3.217    1.089    1.736    1.956 
##      .25      .50      .75      .90      .95 
##    2.581    3.325    3.610    4.048    5.293 
## 
## lowest : 1.513 1.615 1.835 1.935 2.140, highest: 3.845 4.070 5.250 5.345 5.424
## ---------------------------------------------------------------------------

通过pastecs包中的stat.desc()函数计算描述性统计量

pastecs包中有一个名为stat.desc()的函数,它可以计算种类繁多的描述性统计量。使用 格式为: stat.desc(x, basic=TRUE, desc=TRUE, norm=FALSE, p=0.95) 其中的x是一个数据框或时间序列。若basic=TRUE(默认值),则计算其中所有值、空值、缺失 值的数量,以及最小值、最大值、值域,还有总和。若desc=TRUE(同样也是默认值),则计算 中位数、平均数、平均数的标准误、平均数置信度为95%的置信区间、方差、标准差以及变异系 数。最后,若norm=TRUE(不是默认的),则返回正态分布统计量,包括偏度和峰度(以及它们 的统计显著程度)和Shapiro-Wilk正态检验结果。这里使用了p值来计算平均数的置信区间(默认 置信度为0.95)。

#通过pastecs包中的stat.desc()函数计算描述性统计量
library(pastecs) 
## Warning: package 'pastecs' was built under R version 3.5.3
myvars <- c("mpg", "hp", "wt") 
stat.desc(mtcars[myvars])
##                      mpg           hp          wt
## nbr.val       32.0000000   32.0000000  32.0000000
## nbr.null       0.0000000    0.0000000   0.0000000
## nbr.na         0.0000000    0.0000000   0.0000000
## min           10.4000000   52.0000000   1.5130000
## max           33.9000000  335.0000000   5.4240000
## range         23.5000000  283.0000000   3.9110000
## sum          642.9000000 4694.0000000 102.9520000
## median        19.2000000  123.0000000   3.3250000
## mean          20.0906250  146.6875000   3.2172500
## SE.mean        1.0654240   12.1203173   0.1729685
## CI.mean.0.95   2.1729465   24.7195501   0.3527715
## var           36.3241028 4700.8669355   0.9573790
## std.dev        6.0269481   68.5628685   0.9784574
## coef.var       0.2999881    0.4674077   0.3041285

通过psych包中的describe()计算描述性统计量

psych包也拥有一个名为describe()的函数,它可以计算非缺失值的数量、平均数、标准差、中位数、截尾均值、绝对中位差、最小值、最大值、值域、偏度、峰度和平均值的标准误。代码清单7-5中有一个示例。

library(psych) 
## Warning: package 'psych' was built under R version 3.5.3

## 
## Attaching package: 'psych'

## The following object is masked from 'package:Hmisc':
## 
##     describe

## The following objects are masked from 'package:ggplot2':
## 
##     %+%, alpha
myvars <- c("mpg", "hp", "wt") 
describe(mtcars[myvars])
##     vars  n   mean    sd median trimmed   mad   min    max  range skew
## mpg    1 32  20.09  6.03  19.20   19.70  5.41 10.40  33.90  23.50 0.61
## hp     2 32 146.69 68.56 123.00  141.19 77.10 52.00 335.00 283.00 0.73
## wt     3 32   3.22  0.98   3.33    3.15  0.77  1.51   5.42   3.91 0.42
##     kurtosis    se
## mpg    -0.37  1.07
## hp     -0.14 12.12
## wt     -0.02  0.17

注意
在前面的示例中,psych包和Hmisc包均提供了名为describe()的函数。R如何知道该 使用哪个呢?简言之,如代码清单7-5所示,最后载入的程序包优先。在这里,psych在 Hmisc之后被载入,然后显示了一条信息,提示Hmisc包中的describe()函数被psych 包中的同名函数所屏蔽(masked)。键入describe()后,R在搜索这个函数时将首先找 到psych包中的函数并执行它。如果你想改而使用Hmisc包中的版本,可以键入 Hmisc::describe(mt)。这个函数仍然在那里。你只是需要给予R更多信息以找到它。

分组计算描述性统计量

##使用aggregate()分组获取描述性统计量

在比较多组个体或观测时,关注的焦点经常是各组的描述性统计信息,而不是样本整体的描 述性统计信息。同样地,在R中完成这个任务有若干种方法。我们将以获取变速箱类型各水平的 描述性统计量开始。 在第5章中,我们讨论了整合数据的方法。你可以使用aggregate()函数(5.6.2节)来分组 获取描述性统计量,如代码清单7-6所示。

#代码清单7-6 使用aggregate()分组获取描述性统计量
myvars <- c("mpg", "hp", "wt") 
aggregate(mtcars[myvars], by=list(am=mtcars$am), mean)
##   am      mpg       hp       wt
## 1  0 17.14737 160.2632 3.768895
## 2  1 24.39231 126.8462 2.411000
aggregate(mtcars[myvars], by=list(am=mtcars$am), sd)
##   am      mpg       hp        wt
## 1  0 3.833966 53.90820 0.7774001
## 2  1 6.166504 84.06232 0.6169816

注意
list(am=mtcars(am)的使用。如果使用的是list(mtcars)am),则am列将被标注为 Group.1而不是am。你使用这个赋值指定了一个更有帮助的列标签。如果有多个分组变量,可以 使用by=list(name1=groupvar1, name2=groupvar2, … , nameN=groupvarN)这样的 语句。

使用aggregate()分组获取描述性统计量

遗憾的是,aggregate()仅允许在每次调用中使用平均数、标准差这样的单返回值函数。 它无法一次返回若干个统计量。要完成这项任务,可以使用by()函数。格式为: by(data, INDICES, FUN) 其中data是一个数据框或矩阵,INDICES是一个因子或因子组成的列表,定义了分组,FUN是任 意函数。

#使用by()分组计算描述性统计量
dstats <- function(x)sapply(x, mystats) 
myvars <- c("mpg", "hp", "wt") 
by(mtcars[myvars], mtcars$am, dstats)
## mtcars$am: 0
##                  mpg           hp         wt
## n        19.00000000  19.00000000 19.0000000
## mean     17.14736842 160.26315789  3.7688947
## stdev     3.83396639  53.90819573  0.7774001
## skew      0.01395038  -0.01422519  0.9759294
## kurtosis -0.80317826  -1.20969733  0.1415676
## -------------------------------------------------------- 
## mtcars$am: 1
##                  mpg          hp         wt
## n        13.00000000  13.0000000 13.0000000
## mean     24.39230769 126.8461538  2.4110000
## stdev     6.16650381  84.0623243  0.6169816
## skew      0.05256118   1.3598859  0.2103128
## kurtosis -1.45535200   0.5634635 -1.1737358

这里,dstats()调用了代码清单7-2中的mystats()函数,将其应用于数据框的每一栏中。 再通过by()函数则可得到am中每一水平的概括统计量。

使用doBy包中的summaryBy()分组计算概述统计量

doBy包中summaryBy()函数的使用格式为: summaryBy(formula, data=dataframe, FUN=function), 其中的formula接受以下的格式:
var1 + var2 + var3 + ... + varN ~ groupvar1 + groupvar2 + ... + groupvarN, 在~左侧的变量是需要分析的数值型变量,而右侧的变量是类别型的分组变量。function可为 任何内建或用户自编的R函数。使用7.2.1节中创建的mystats()函数的一个示例如代码清单7-8 所示。

#使用doBy包中的summaryBy()分组计算概述统计量
library(doBy) 
## Warning: package 'doBy' was built under R version 3.5.3
summaryBy(mpg+hp+wt~am, data=mtcars, FUN=mystats)
##   am mpg.n mpg.mean mpg.stdev   mpg.skew mpg.kurtosis hp.n  hp.mean
## 1  0    19 17.14737  3.833966 0.01395038   -0.8031783   19 160.2632
## 2  1    13 24.39231  6.166504 0.05256118   -1.4553520   13 126.8462
##   hp.stdev     hp.skew hp.kurtosis wt.n  wt.mean  wt.stdev   wt.skew
## 1 53.90820 -0.01422519  -1.2096973   19 3.768895 0.7774001 0.9759294
## 2 84.06232  1.35988586   0.5634635   13 2.411000 0.6169816 0.2103128
##   wt.kurtosis
## 1   0.1415676
## 2  -1.1737358

使用psych包中的describeBy()分组计算概述统计量

psych包中的describeBy()函数可计算和describe()相同的描述性统计量,只是按照一 个或多个分组变量分层。

library(psych) 
myvars <- c("mpg", "hp", "wt") 
describeBy(mtcars[myvars], list(am=mtcars$am))
## 
##  Descriptive statistics by group 
## am: 0
##     vars  n   mean    sd median trimmed   mad   min    max  range  skew
## mpg    1 19  17.15  3.83  17.30   17.12  3.11 10.40  24.40  14.00  0.01
## hp     2 19 160.26 53.91 175.00  161.06 77.10 62.00 245.00 183.00 -0.01
## wt     3 19   3.77  0.78   3.52    3.75  0.45  2.46   5.42   2.96  0.98
##     kurtosis    se
## mpg    -0.80  0.88
## hp     -1.21 12.37
## wt      0.14  0.18
## -------------------------------------------------------- 
## am: 1
##     vars  n   mean    sd median trimmed   mad   min    max  range skew
## mpg    1 13  24.39  6.17  22.80   24.38  6.67 15.00  33.90  18.90 0.05
## hp     2 13 126.85 84.06 109.00  114.73 63.75 52.00 335.00 283.00 1.36
## wt     3 13   2.41  0.62   2.32    2.39  0.68  1.51   3.57   2.06 0.21
##     kurtosis    se
## mpg    -1.46  1.71
## hp      0.56 23.31
## wt     -1.17  0.17