第7章--基本统计分析

时间 2019-12-08

标签基本统计分析繁體版

原文原文链接

7.1 描述性统计分析

> vars <- c("mpg", "hp", "wt")
> head(mtcars[vars])
                   mpg  hp    wt
Mazda RX4         21.0 110 2.620
Mazda RX4 Wag     21.0 110 2.875
Datsun 710        22.8  93 2.320
Hornet 4 Drive    21.4 110 3.215
Hornet Sportabout 18.7 175 3.440
Valiant           18.1 105 3.460

以上述数据集为例，对于基础安装R的用户，可使用summary()函数来获取描述性统计量。这个函数提供了最小值、最大值、四分位数和数值型变量的均值，以及因子向量和逻辑型向量的频数统计。ide

> summary(mtcars[vars])
      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

分组计算描述性统计量函数

你可使用aggregate()函数来分组获取描述性统计量。ui

> aggregate(mtcars[vars],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[vars],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的使用。若是使用的是list(mtcars$am)，则am列将被标注为Group.1而不是am。你使用这个赋值指定了一个更有帮助的列标签。若是有多个分组变量，你可使用一下语句：spa

by=list(name1=groupvar1, name2=groupvar2,...)

而aggregate函数仅容许每次调用中使用平均数、标准差这样的单返回值函数。若是想要一次返回若干个统计量，可使用by()函数。设计

7.2 频数表和列联表

下面的示例中，假设A、B和C表明类别型变量。code

> library(vcd)
> library(grid)
> head(Arthritis)
  ID Treatment  Sex Age Improved
1 57   Treated Male  27     Some
2 46   Treated Male  29     None
3 77   Treated Male  30     None
4 17   Treated Male  32   Marked
5 36   Treated Male  46   Marked
6 23   Treated Male  58   Marked

1. 生成频数表

1.1 一维列联表blog

可使用table()函数生成简单的频数统计表。ip

> mytable <- with(Arthritis, table(Improved))
> mytable
Improved
  None   Some Marked 
    42     14     28

可使用prop.table()将这些频数转换为比例值，或使用prop.table()*100来转换为百分比：it

> prop.table(mytable)
Improved
     None      Some    Marked 
0.5000000 0.1666667 0.3333333

1.2 二维列联表io

对于二维表，table的使用格式为：

mytable <- table(A, B)

其中，A是行变量，B是列变量。

除此以外，xtabs()函数还可以使用公式风格的输入建立列联表。

mytable <- xtabs(~ A + B, data=mydata)

其中，mydata是一个矩阵或者数据框，要进行交叉分类的变量应出如今公式的右侧：

> mytable <- xtabs(~ Treatment+Improved, data=Arthritis)
> mytable
         Improved
Treatment None Some Marked
  Placebo   29    7      7
  Treated   13    7     21

可使用margin.table()和prop.table()函数分别生成边际频数和比例。

> margin.table(mytable,1)
Treatment
Placebo Treated 
     43      41 
> prop.table(mytable,1)
         Improved
Treatment      None      Some    Marked
  Placebo 0.6744186 0.1627907 0.1627907
  Treated 0.3170732 0.1707317 0.5121951

下标1指代table语句的第一个变量，列和列比例能够这样计算：

> margin.table(mytable,2)
Improved
  None   Some Marked 
    42     14     28 
> prop.table(mytable,2)
         Improved
Treatment      None      Some    Marked
  Placebo 0.6904762 0.5000000 0.2500000
  Treated 0.3095238 0.5000000 0.7500000

各单元格所占比例可用以下语句获取：

> prop.table(mytable)
         Improved
Treatment       None       Some     Marked
  Placebo 0.34523810 0.08333333 0.08333333
  Treated 0.15476190 0.08333333 0.25000000

可使用addmargins()函数为这些表格添加边际和。

3. 多维列联表

除了上面二维列联表里介绍的方法以外，还有ftable()函数能够以一种紧凑而吸引人的方式输出多维列联表。

> mytable <- xtabs(~Treatment+Sex+Improved, data=Arthritis)
> mytable
, , Improved = None

         Sex
Treatment Female Male
  Placebo     19   10
  Treated      6    7

, , Improved = Some

         Sex
Treatment Female Male
  Placebo      7    0
  Treated      5    2

, , Improved = Marked

         Sex
Treatment Female Male
  Placebo      6    1
  Treated     16    5

> ftable(mytable)
                 Improved None Some Marked
Treatment Sex                             
Placebo   Female            19    7      6
          Male              10    0      1
Treated   Female             6    5     16
          Male               7    2      5

> margin.table(mytable,1)
Treatment
Placebo Treated 
     43      41 
> margin.table(mytable,2)
Sex
Female   Male 
    59     25 
> margin.table(mytable,3)
Improved
  None   Some Marked 
    42     14     28 
> margin.table(mytable,c(1,3))
         Improved
Treatment None Some Marked
  Placebo   29    7      7
  Treated   13    7     21
> ftable(prop.table(mytable,c(1,2)))
                 Improved       None       Some     Marked
Treatment Sex                                             
Placebo   Female          0.59375000 0.21875000 0.18750000
          Male            0.90909091 0.00000000 0.09090909
Treated   Female          0.22222222 0.18518519 0.59259259
          Male            0.50000000 0.14285714 0.35714286
> ftable(addmargins(prop.table(mytable,c(1,2)),3))
                 Improved       None       Some     Marked        Sum
Treatment Sex                                                        
Placebo   Female          0.59375000 0.21875000 0.18750000 1.00000000
          Male            0.90909091 0.00000000 0.09090909 1.00000000
Treated   Female          0.22222222 0.18518519 0.59259259 1.00000000
          Male            0.50000000 0.14285714 0.35714286 1.00000000
>

2. 独立性检验

1. 卡方独立性检验

可使用chisq.test()函数对二维列表的行变量和列变量进行卡方独立性检验。

> mytable <- xtabs(~Treatment+Improved, data=Arthritis)
> chisq.test(mytable)

        Pearson's Chi-squared test

data:  mytable
X-squared = 13.055, df = 2, p-value = 0.001463

> mytable <- xtabs(~Treatment+Sex, data=Arthritis)
> chisq.test(mytable)

        Pearson's Chi-squared test with Yates' continuity correction

data:  mytable
X-squared = 0.38378, df = 1, p-value = 0.5356

Treatment和Improved有相关性，和Sex没有相关性。

2. Fisher精确检验

使用fisher.test()函数进行Fisher精确检验，Fisher精确检验的原假设是：边界固定的列联表中行和列是相互独立的。

> mytable <- xtabs(~Treatment+Improved, data=Arthritis)
> fisher.test(mytable)

        Fisher's Exact Test for Count Data

data:  mytable
p-value = 0.001393
alternative hypothesis: two.sided

3. Cochran-Mantel-Haenszel检验

mantelhaen.test()函数可用来进行该检验，其原假设是：两个名义变量在第三个变量的每一层中都是条件独立的。

> mytable <- xtabs(~Treatment+Improved+Sex, data=Arthritis)
> mantelhaen.test(mytable)

        Cochran-Mantel-Haenszel test

data:  mytable
Cochran-Mantel-Haenszel M^2 = 14.632, df = 2, p-value = 0.0006647

7.3 相关

1. 相关的类型

Pearson，Spearman和Kendall相关

Pearson积差关系数衡量了两个定量变量之间的线性相关程度。

Spearman等级相关系数则衡量分级定序变量之间的相关程度。

Kendall's Tau相关系数也是一种非参数的等级相关度量。

cor()函数能够计算这三种相关系数，而cov()函数可用来计算协方差。

cor(x, use= ,method= )

x指一个矩阵或者数据框；

use是指定数据缺失的处理方式；

method是指定相关系数的类型，默认为Pearson类型。

> states <- state.x77[,1:6]
> cov(states)
              Population      Income   Illiteracy     Life Exp      Murder      HS Grad
Population 19931683.7588 571229.7796  292.8679592 -407.8424612 5663.523714 -3551.509551
Income       571229.7796 377573.3061 -163.7020408  280.6631837 -521.894286  3076.768980
Illiteracy      292.8680   -163.7020    0.3715306   -0.4815122    1.581776    -3.235469
Life Exp       -407.8425    280.6632   -0.4815122    1.8020204   -3.869480     6.312685
Murder         5663.5237   -521.8943    1.5817755   -3.8694804   13.627465   -14.549616
HS Grad       -3551.5096   3076.7690   -3.2354694    6.3126849  -14.549616    65.237894
> cor(states)
            Population     Income Illiteracy    Life Exp     Murder     HS Grad
Population  1.00000000  0.2082276  0.1076224 -0.06805195  0.3436428 -0.09848975
Income      0.20822756  1.0000000 -0.4370752  0.34025534 -0.2300776  0.61993232
Illiteracy  0.10762237 -0.4370752  1.0000000 -0.58847793  0.7029752 -0.65718861
Life Exp   -0.06805195  0.3402553 -0.5884779  1.00000000 -0.7808458  0.58221620
Murder      0.34364275 -0.2300776  0.7029752 -0.78084575  1.0000000 -0.48797102
HS Grad    -0.09848975  0.6199323 -0.6571886  0.58221620 -0.4879710  1.00000000
> cor(states,method="spearman")
           Population     Income Illiteracy   Life Exp     Murder    HS Grad
Population  1.0000000  0.1246098  0.3130496 -0.1040171  0.3457401 -0.3833649
Income      0.1246098  1.0000000 -0.3145948  0.3241050 -0.2174623  0.5104809
Illiteracy  0.3130496 -0.3145948  1.0000000 -0.5553735  0.6723592 -0.6545396
Life Exp   -0.1040171  0.3241050 -0.5553735  1.0000000 -0.7802406  0.5239410
Murder      0.3457401 -0.2174623  0.6723592 -0.7802406  1.0000000 -0.4367330
HS Grad    -0.3833649  0.5104809 -0.6545396  0.5239410 -0.4367330  1.0000000

在默认状况下，获得的结果是一个方阵（全部变量之间两两计算相关）。你一样能够计算非方形的相关矩阵。

> x <- states[,c("Population", "Income", "Illiteracy", "HS Grad")]
> y <- states[,c("Life Exp", "Murder")]
> cor(x, y)
              Life Exp     Murder
Population -0.06805195  0.3436428
Income      0.34025534 -0.2300776
Illiteracy -0.58847793  0.7029752
HS Grad     0.58221620 -0.4879710

偏相关

偏相关是控制一个或多个定量变量时，另外两个定量变量之间的相关系数。可使用ggm包中的pcor()函数计算偏相关系数。

> pcor(c(1,5,2,3,6), cov(states))
[1] 0.3462724

7.3 相关性的显著性检验

经常使用的原假设为变量间不相关（即整体的相关系数为0）。你可使用cor.test()函数对单个的Pearson、Spearman和Kendall相关系数进行检验。简化后的格式为：

cor.test(x,y, alternative= , method= ,)

x和y为要检验相关性的变量；

alternative用来指定进行双侧检验或单侧检验；

method用来指定要计算的相关类型。

> cor.test(states[,3],states[,5])

        Pearson's product-moment correlation

data:  states[, 3] and states[, 5]
t = 6.8479, df = 48, p-value = 1.258e-08
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
 0.5279280 0.8207295
sample estimates:
      cor 
0.7029752

可是，cor.test()一次只能检验一种相关关系。psych包中提供的corr.test()函数能够一次作更多事情。

> library(psych)
> corr.test(states, use="complete")
Call:corr.test(x = states, use = "complete")
Correlation matrix 
           Population Income Illiteracy Life Exp Murder HS Grad
Population       1.00   0.21       0.11    -0.07   0.34   -0.10
Income           0.21   1.00      -0.44     0.34  -0.23    0.62
Illiteracy       0.11  -0.44       1.00    -0.59   0.70   -0.66
Life Exp        -0.07   0.34      -0.59     1.00  -0.78    0.58
Murder           0.34  -0.23       0.70    -0.78   1.00   -0.49
HS Grad         -0.10   0.62      -0.66     0.58  -0.49    1.00
Sample Size 
[1] 50
Probability values (Entries above the diagonal are adjusted for multiple tests.) 
           Population Income Illiteracy Life Exp Murder HS Grad
Population       0.00   0.59       1.00      1.0   0.10       1
Income           0.15   0.00       0.01      0.1   0.54       0
Illiteracy       0.46   0.00       0.00      0.0   0.00       0
Life Exp         0.64   0.02       0.00      0.0   0.00       0
Murder           0.01   0.11       0.00      0.0   0.00       0
HS Grad          0.50   0.00       0.00      0.0   0.00       0

 To see confidence intervals of the correlations, print with the short=FALSE option

参数use=的取值可为"pairwise"或"complete"（分别表示对缺失值执行成对删除或行删除）。

参数method=的取值是三种方法，默认为pearson。

7.4 t检验

1. 独立样本的t检验

一个针对两组的独立样本t检验能够用于检验两个整体的均值相等的假设。这里假设两组数据是独立的，而且是从正态整体中抽的的。

t.test(y ~ x, data)

> library(MASS)
> t.test(Prob ~ So, data=UScrime)

        Welch Two Sample t-test

data:  Prob by So
t = -3.8954, df = 24.925, p-value = 0.0006506
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 -0.03852569 -0.01187439
sample estimates:
mean in group 0 mean in group 1 
     0.03851265      0.06371269

2. 非独立样本的t检验

在两组观测之间相关时，你得到的是一个非独立组设计。非独立样本的t检验假定组间的差别呈正态分布。