生物统计与 R:描述统计学

计量数据的位置

以一组舒张压数据为例子。

> data <- c(67,72,76,70,73,81,89,74,61,86,61,76,74,76,79,69,68,88,75,63,80,72,78,67,56,83,65,73,65,57,55,81,70,72,67,77,74)

均值 Arithmatic Mean

均值严格实际上指的是 average,常常说的均值用来代替的是 arithmatic mean, 即算术平均数。算数平均数是用所有数据的和除以数据的个数。 R 中有相应的函数mean,如果不嫌麻烦或者有特殊的目的,也可以自己去算。

> mean(data)
# [1] 72.16216
> sum(data) / length(data)
# [1] 72.16216

虽然使用R函数mean或者自己根据定义来算看起来结果没什么差别, 但是实际上使用 sum 来算和 mean 来算在效率上还是有差别的。 因为虽然 sum 和 mean 的底层都是用C写的函数来运算的, 不过mean函数包装了很多判断的部分,所以会降低一定效率。 然而,效率不是R主要考虑的问题, R用一定的效率去换取了更方便的交互使用体验是值得的。

> system.time(for(i in 1:10000) mean(data))
#    user  system elapsed
#   0.024   0.000   0.023
> system.time(for(i in 1:10000) .Internal(mean(data)))
#    user  system elapsed
#   0.004   0.000   0.002
> system.time(for(i in 1:10000) sum(data)/length(data))
#    user  system elapsed
#   0.004   0.000   0.004

均值的好处是可以使用所有的数字,并且最简单明了, 然而均值的一个缺点是非常受到极值的影响。 假设舒张压数据有一个数据记错了,记录为原来的10倍, 均值就会受到非常大的影响。

> faultdata <- data
> faultdata[7] <- data[7] * 10
> mean(faultdata)
# [1] 93.81081
> mean(data)
# [1] 72.16216

中位数 Median

中位数是位于所有样本值的排序后中间的值。 如果样本有2n+1个,中位数就是排序后的样本值的第n个, 如果样本有2n个,中位数就是第n个和第n+1个的均值,n为正整数。

因为中位数是用了排序之后的数据中中间的值(或者中间两个值的均值), 所以中位数不会像均值那样受到极大或者极小的极值影响。 然而,因为所有的数据中位数计算实际上只使用了一个或者两个数字, 所以中位数对数据分布的反映并不全面,数据利用率低。

因为计算中位数需要排序,所以使用计算机去计算中位数会方便很多。

> median(data)
# [1] 73

我们也可以自己计算中位数。

> length(data) %% 2
# [1] 1
> sort(data)[(length(data) + 1)/2]
# [1] 73

众数 Mode

众数是数据中出现频次最多的数字,因为对于正态分布来说, 接近其期望的数字出现概率最大,所以众数会接近于期望。 但是这并不能保证众数最接近于期望,对于双峰分布或多峰分布, 众数甚至可能并不唯一。

拿上面的数据为例,使用 R 语言的 table 函数可以计算每个数字的频次。

> table(data)
# data
# 55 56 57 61 63 65 67 68 69 70 72 73 74 75 76 77 78 79 80 81 83 86 88 89
#  1  1  1  2  1  2  3  1  1  2  3  2  3  1  3  1  1  1  1  2  1  1  1  1
> table(table(data))
#  1  2  3
# 15  5  4

可以看到对于例子中有四个值的频次最多,都是3次,所以有三个众数, 分别为:67,72,74,76。

> sort(table(data))
# data
# 55 56 57 63 68 69 75 77 78 79 80 83 86 88 89 61 65 70 73 81 67 72 74 76
#  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  2  2  2  2  2  3  3  3  3

计量数据的范围

如果仅仅使用均值等用于描述数据位置的数据, 我们并不能很好地了解数据的分布情况。 因为对于两组均值一样的数据,我们可以有完全不一样的分布类型。

> dat1 <- c(4, 5, 5, 2, 3, 4, 3, 4, 3, 4)
> mean(dat1)
# [1] 3.7
> dat2 <- c(2, 3, 3, 0, 1, 2, 1, 2, 11, 12)
> mean(dat2)
# [1] 3.7

所以我们也需要去描述数据分布的范围。

范围 Range

为了表示数据的范围,最简单的方法就是找出数据的最大值最小值, 我们就有信心,所有的数据都会落在最大值和最小值之间。

> c(min(dat1), max(dat1))
# [1] 2 5
> c(min(dat1), max(dat1))
# [1] 0 12

然而,最大值最小值非常受到极值的影响,两组数据虽然可能绝大多数相近, 但只要有一组数据有极大值,就会使得数据的分布范围不同。

分位数 Quantile

中位数可以用来显示数据的位置,中位数数实际上是二分之一分位数, 我们有信心是有一半的数据大于中位数,一半的数据小于中位数。 那么我们可以使用分位数来表示数据的分布。 例如使用四分之一分位数和四分之三分位数, 我们可以知道有二分之一的数据在四分之一分位数和四分之三分位数之间。 之所以不使用最大值最小值来反映数据的分布,也是为了避免极值的影响。

> quantile(dat1)
#   0%  25%  50%  75% 100%
#    2    3    4    4    5
> quantile(dat2)
#    0%   25%   50%   75%  100%
#  0.00  1.25  2.00  3.00 12.00

方差和标准差 Variance and Standard Deviation

分位数和中位数一样对于数据的使用效率相对低。 可以使用具体的每一个数据相对于均值的差别的总和来表示数据整体相对于均值的变化程度。 有多种方法可以表示,但是(样本)方差和标准差是特性最好的方式。 方差是标准差的平方。样本方差除以的是n-1而不是n。

R 函数中 sd 函数计算标准差,var可以计算样本方差。

> sd(dat1)
# [1] 0.9486833
> var(dat1)
# [1] 0.9
> sd(dat1) ^ 2
# [1] 0.9
> sd(dat2)
# [1] 4.217688
> var(dat2)
# [1] 17.78889
> sd(dat2) ^ 2
# [1] 17.78889

如果需要计算样本相对于均值差的平方和,可以使用var函数辅助计算, 但应该记住样本方差计算中自由度为 n-1。

> mean(dat1)
# [1] 3.7
> sum((dat1 - mean(dat1))^2)
# [1] 8.1
> var(dat1)
# [1] 0.9
> var(dat1) * (length(dat1) - 1)
# [1] 8.1

变异系数 Coefficient of Variance

对于一个均值为 10 的分布和一个均值为 1000 的分布来说,100的标准差的含义是非常不同的。 所以需要矫正掉均值的影响。 变异系数,即 coefficient of variance,CV,的定义就是标准差除以均值:sd/mean。

> sd(dat1)/mean(dat1)
# [1] 0.2564009
> sd(dat2)/mean(dat2)
# [1] 1.139916

数据分组

对于众多的数据,直接很难去看明白,也可以通过把数据分成若干组的方式去看数据。 进行分组需要去注意: 1. 需要把数据分位k组,从最小的y1开始,到最大的yk+1结束 2. 每一组包括 y1 但不包括 y2,第二组包括y2但不包括y3,每个数据只能分在一组 3. 每组的范围一般应该相同 4. 计算每组中数据的个数 5. 总结每组数据频次

R 中可以使用 cut 函数去分组,使用 table 函数计算频次。

> cut(data, 5)
# [1] (61.8,68.6] (68.6,75.4] (75.4,82.2] (68.6,75.4] (68.6,75.4] (75.4,82.2]
# [7] (82.2,89]   (68.6,75.4] (55,61.8]   (82.2,89]   (55,61.8]   (75.4,82.2]
# [13] (68.6,75.4] (75.4,82.2] (75.4,82.2] (68.6,75.4] (61.8,68.6] (82.2,89]
# [19] (68.6,75.4] (61.8,68.6] (75.4,82.2] (68.6,75.4] (75.4,82.2] (61.8,68.6]
# [25] (55,61.8]   (82.2,89]   (61.8,68.6] (68.6,75.4] (61.8,68.6] (55,61.8]
# [31] (55,61.8]   (75.4,82.2] (68.6,75.4] (68.6,75.4] (61.8,68.6] (75.4,82.2]
# [37] (68.6,75.4]
# Levels: (55,61.8] (61.8,68.6] (68.6,75.4] (75.4,82.2] (82.2,89]

> table(cut(data, 5))
#   (55,61.8] (61.8,68.6] (68.6,75.4] (75.4,82.2]   (82.2,89]
#           5           7          12           9           4

图形描述

使用条形图(bar plot),茎叶图(stem-leaf plot), 盒型图(box plot)等可以描述数据。

条形图 bar plot

> barplot(table(cut(data, 5)))
> library(ggplot2)
> ggplot(data.frame(x=data), aes(x=x)) + geom_bar()

茎叶图 stem-leaf plot

> stem(data)
#   The decimal point is 1 digit(s) to the right of the |

#   5 | 567
#   6 | 113
#   6 | 5577789
#   7 | 0022233444
#   7 | 5666789
#   8 | 0113
#   8 | 689

对于茎叶图来说,数线左边的是高位,上面例子是十位, 右边每个数字代表一个数据点的个位数。 茎叶图和条形图有相似的功能,但是可以显示具体的数值。

盒型图 box plot

盒型图实际上是用图来表示分位数。 * 盒子最上面的横线表示四分之三分位数 * 盒子最下面的横线表示四分之一分位数 * 盒子中间的横线表示中位数 * 盒子的长度表示四分之三分位数和四分之一分位数的差 * 盒子最上的胡须最长为1.5倍的盒长 * 盒子最下的胡须最长为1.5倍的盒长 * 在四分之三分位数加1.5倍盒长之外的点为离群值outlier * 在四分之一分位数减1.5倍盒长之外的点为离群值outlier * 在四分之三分位数加3倍盒长之外的点为极端离群值 * 在四分之一分位数减3倍盒长之外的点为极端离群值

> boxplot(data)
By @Wolfson Liu in
Tags : #r, #statistics,