统计检验(基于R函数)

前言

非参数检验(卡方检验),参数检验(F检验,T检验,Z检验),方差分析(ANOVA)
非参数检验与参数检验区别:主要差异在于,非参数检验不需要假定总体分布形式,直接对数据的分布进行检验。由于不涉及总体分布的参数,故名「非参数」检验。比如,卡方检验。而参数检验一般需要正太性(符合正太分布)、方差齐次等假设,并已知总体均值,方差等值,或者从样本估计。

p value

一般在做假设检验的时候,我们可以通过如下步骤:
1、设定alpha
2、计算统计量t
3、根据alpha查统计量阈值来确定拒绝还是接受(这里是比较统计量t即可,t>t阈)

第二种方法,直接计算p value(对于不同假设分布有不同计算pvalue的公式)
p value的含义是:


image.png

在H0成立的情况下,Data产生的概率
非常好理解,p值为在H0假设下,最终会产生得到当前数据的概率。
如果p < alpha(比如0.05)则在H0的假设下,Data产生概率非常小(小于显著水平alpha),则应该拒绝H0
(注:显著性水平alpha与p值pval不是一回事,不要搞混。alpha为我们设定的显著性水平,pval为计算出来概率。alpha是我们主观去设定的值,pvalue为伴随数据客观得到的值)

第一类和第二类错误

第一类错误:H0为真,但是拒绝了H0(弃真)。概率为alpha(常用的显著性水平,即是达到5%概率犯第一类错误的水平【即很大概率不会错误地拒绝H0】)
一般我们计算时,p值就是真实情况下,犯第一类错误的概率。
第二类错误:H0为假,但是接受了H0(取伪)。当我们设定了alpha后,如果显著性水平没有达到alpha的值,那么我们会接受假设H0。在这种情况下,仍然有一定概率H0为假。这个概率为Beta。
1-β即为statistical power
通常,alpha设置的比较大,则beta也会比较大,如下图:
在统计量为蓝色竖线的时候,根据H0假设的分布,可以得到alpha,同时beta由H1的分布计算得来,但是通常我们不知道H1的真实分布,所以也无法精确计算beta
实践中,由于beta无法精确计算,但是我们还是想要降低第二类错误的概率,所以要做有如下的认知,beta的大小取决于两个因素:
1、H1假设与H0的距离,距离越大,beta越小。
2、alpha的值。alpha越大,beta越小。(更倾向于拒绝H0,一类错误概率增加,则二类错误率减少)。
通常在固定alpha,H1假设的情况下,要缩小beta需要增加样本量。

置信度和置信区间

置信水平(置信度)=1-显著性水平(alpha)
置信区间:U,V。指在一定的置信水平下(1-alpha),被观测参数的真实值会落在区间(U,V)之内。
PS:什么是统计量,统计量是样本的函数,且不依赖于任何未知的参数。比如样本均值就是个常用的统计量,mean(sample),它只依赖于所有样本的值。
PSS:一般来说,主要用于做区间估计。

自由度

统计学上,自由度是指当以样本的统计量来估计总体的参数时,样本中独立或能自由变化的数据的个数,称为该统计量的自由度。一般来说,自由度等于独立变量减掉其衍生量数。举例来说,变异数的定义是样本减平均值(一个由样本决定的衍生量),因此对N个随机样本而言,其自由度为N-1

如何选择一个合适的统计检验:(参考https://rpubs.com/Argaadya/bayesian_ab

image.png

连续型数据

基于正态分布的检验

均值检验(T检验)

主要用于小样本(样本容量小于30)的两个平均值差异程度的检验方法(主要适用与总体方差未知的情况:即要用小样本的方差预估总体方差)
T检验也是一种参数检验。(对方差齐性敏感,需要先检验方差齐性【F检验】)

  • 对于与总体配对的T检验,要求满足方差齐性,才能反应其均值的差异。【不齐性有特殊处理方法】
  • 适用于已知总体均值(单样本配对时需要总体均值,双样本不需要),样本均值,样本方差(即样本少总体方差未知),且大致来自于正太分布(一般除非明显的长尾多峰等分布以外,都大致可以检验)

基础延伸:
H0假设 :u与u0没有差异,u0其实就是一个comparison value,u=u0
H1假设(单边):样本A均值大于B(或者小于)u>u0
H1假设(双边):两个样本均值有差异u>u0 or u<u0
使用单边双边检验的区别:看你拒绝假设在两边,或者一边
均值检验中统计量

image.png

服从d=n-1的T分布。n为样本量
其泛化形式为
image.png

其中E(x0)为H0假设的对比值(comparison value),S(x)为对应变量x的标准差(一般总体参数未知)所以都是用样本标准差S,当我们对比的是均值时,均值的标准差用
image.png

  • 单总体样本:
    H0:样本均值与总体均值相同
    image.png

    ¯x为样本均值,u为总体均值,Sx为样本标准差,
    Sx/根号n 为样本均值的标准差
  • 双总体配对样本:(形式与上述一致)

    u0为0假设下差异值的均值。s为配对样本差值的标准差,s/根号n 为配对样本差值的均值的标准差
  • 双总体非配对样本:
    image.png

注:两个独立样本的差的variance是:


image.png

这个公式是一个比值。一个普通的比喻是,t值为信噪比。
配对t检验除了自由度有所差异,其实本身也差不多。
上述公式本质:分子都是均值,分母都是方差(均值的方差)。所以和z score本质也一样。
1)根据待检验的实验,计算出t值
2)根据t值查表获取p值,借此判断是否有统计学上的差异(不过通常做法是根据想要的p值【比如0.05】来查表获得临界的t值,再比较t值的大小来判断是否达到一定的显著性,p>临界t值,则效果显著)PS:p值代表原假设为真时,此事件发生的概率。【如果很小的话,则证明不太可能发生这样的情况,应该拒绝原假设,当然,拒绝也是有错误的可能性的,错误的概率也为p值,即有显著差异其误判的概率是p】

用于正态总体均值假设检验,单样本,双样本都可以。
结果意义:P值小于显著性水平时拒绝原假设,否则,接受原假设。具体的假设要看所选择的是双边假设还是单边假设(又分小于和大于)

单样本t检验

t.test(extra ~ 1, data = sleep)
#>
#>  One Sample t-test
#>
#> data:  extra
#> t = 3.413, df = 19, p-value = 0.002918
#> alternative hypothesis: true mean is not equal to 0
#> 95 percent confidence interval:
#>  0.5955845 2.4844155
#> sample estimates:
#> mean of x 
#>      1.54 

p=0.002918 < 0.05,拒绝原假设,接受备择假设,样本均值不等于0

双样本t检验

t.test(1:10, 10:20)
#> 
#>  Welch Two Sample t-test
#> 
#> data:  1:10 and 10:20
#> t = -7, df = 19, p-value = 2e-06
#> alternative hypothesis: true difference in means is not equal to 0
#> 95 percent confidence interval:
#>  -12.4  -6.6
#> sample estimates:
#> mean of x mean of y 
#>       5.5      15.0

p=2e-06 < 0.05,则拒绝原假设,接受备择假设。说明两样本存在均值差异(两样本的均值差不等于0)

配对 t 检验

## 使用公式进行配对t检验
## The sleep data are actually paired, so could have been in wide format:
sleep2 <- reshape(sleep, direction = "wide", 
                  idvar = "ID", timevar = "group")
t.test(Pair(extra.1, extra.2) ~ 1, data = sleep2)
#> 
#>  Paired t-test
#> 
#> data:  sleep2$extra.1 and sleep2$extra.2
#> t = -4.0621, df = 9, p-value = 0.002833
#> alternative hypothesis: true mean difference is not equal to 0
#> 95 percent confidence interval:
#>  -2.4598858 -0.7001142
#> sample estimates:
#> mean difference 
#>           -1.58 

## 不使用公式进行配对t检验
t.test(sleep2$extra.1,sleep2$extra.2,paired=T)
#> 
#>  Paired t-test
#> 
#> data:  sleep2$extra.1 and sleep2$extra.2
#> t = -4.0621, df = 9, p-value = 0.002833
#> alternative hypothesis: true mean difference is not equal to 0
#> 95 percent confidence interval:
#>  -2.4598858 -0.7001142
#> sample estimates:
#> mean difference 
#>           -1.58

方差检验(F检验-方差齐性检验)

检测两个及以上的样本总体方差差异是否显著的检验。(对正太性敏感,需要前提为正太分布)
基础延伸:
H0假设:两个样本的方差没有差异(齐次)σ1 = σ2
H1对立假设:两个样本的方差有差异σ1 > σ2
描述统计量F=σ1/σ2 服从 d1=n1-1, d2=n2-1 的F分布

服从正态分布的两组样本总体方差比较

df <- data.frame(
  value = c(rnorm(10), rnorm(10, mean = 1)),
  group = c(rep("control", 10), rep("test", 10))
)
# 进行的是 F 检验
var.test(value ~ group, data = df)
#> 
#>  F test to compare two variances
#> 
#> data:  value by group
#> F = 0.4, num df = 9, denom df = 9, p-value = 0.2
#> alternative hypothesis: true ratio of variances is not equal to 1
#> 95 percent confidence interval:
#>  0.103 1.671
#> sample estimates:
#> ratio of variances 
#>              0.415

多组样本方差比较

使用 Bartlett 检验比较每个组(样本)数据的方差是否一致。

bartlett.test(count ~ spray, data = InsectSprays)
#> 
#>  Bartlett test of homogeneity of variances
#> 
#> data:  count by spray
#> Bartlett's K-squared = 25.96, df = 5, p-value = 9.085e-05
# k-squared检验统计量的近似卡方分布的自由度

多样本均值比较

对于两组以上数据间均值的比较,使用方差分析 ANOVA

aov(wt ~ factor(cyl), data = mtcars)
#> Call:
#>    aov(formula = wt ~ factor(cyl), data = mtcars)
#> 
#> Terms:
#>                 factor(cyl) Residuals
#> Sum of Squares         18.2      11.5
#> Deg. of Freedom           2        29
#> 
#> Residual standard error: 0.63
#> Estimated effects may be unbalanced

查看详细信息

model.tables(aov(wt ~ factor(cyl), data = mtcars))
#> Tables of effects
#> 
#>  factor(cyl) 
#>           4       6      8
#>     -0.9315 -0.1001  0.782
#> rep 11.0000  7.0000 14.000

通常先用 lm() 函数对数据建立线性模型,再用 anova() 函数提取方差分析的信息更方便。

md <- lm(wt ~ factor(cyl), data = mtcars)
anova(md)
#> Analysis of Variance Table
#> 
#> Response: wt
#>             Df Sum Sq Mean Sq F value  Pr(>F)    
#> factor(cyl)  2   18.2    9.09    22.9 1.1e-06 ***
#> Residuals   29   11.5    0.40                    
#> ---
#> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

ANOVA 分析假设各组样本数据的方差是相等的,如果知道(或怀疑)不相等,可以使用 oneway.test() 函数。

oneway.test(wt ~ cyl, data = mtcars)
#> 
#>  One-way analysis of means (not assuming equal variances)
#> 
#> data:  wt and cyl
#> F = 20, num df = 2, denom df = 19, p-value = 2e-05

多组样本的配对 t 检验

计算组间水平的两两比较,并进行多次检验校正

pairwise.t.test(mtcars$wt, mtcars$cyl)
#> 
#>  Pairwise comparisons using t tests with pooled SD 
#> 
#> data:  mtcars$wt and mtcars$cyl 
#> 
#>   4     6   
#> 6 0.01  -   
#> 8 6e-07 0.01
#> 
#> P value adjustment method: holm

p.adjust.methods   # 可以自定义 p 值校正方法
# c("holm", "hochberg", "hommel", "bonferroni", "BH", "BY",  "fdr", "none")

正态性检验

Shapiro-Wilk 检验:
检验数据是否符合正态分布,R函数:shapiro.test()。
结果含义:当p值小于某个显著性水平α(比如0.05)时,则认为样本不是来自正态分布的总体,否则则承认样本来自正态分布的总体。

shapiro.test(rnorm(30))
#> 
#>  Shapiro-Wilk normality test
#> 
#> data:  rnorm(30)
#> W = 1, p-value = 1

分布的对称性检验

用 Kolmogorov-Smirnov 检验查看一个向量是否来自对称的概率分布(不限于正态分布)
函数第 1 个参数指定待检验的数据,第 2 个参数指定对称分布的类型,可以是数值型向量、指定概率分布函数的字符串或一个分布函数。

ks.test(rnorm(10), "pnorm")
#> 
#>  One-sample Kolmogorov-Smirnov test
#> 
#> data:  rnorm(10)
#> D = 0.4, p-value = 0.09
#> alternative hypothesis: two-sided

ks.test(rpois(10, lambda = 1), "pnorm")
#> Warning in ks.test(rpois(10, lambda = 1), "pnorm"): ties should not be present
#> for the Kolmogorov-Smirnov test
#> 
#>  One-sample Kolmogorov-Smirnov test
#> 
#> data:  rpois(10, lambda = 1)
#> D = 0.5, p-value = 0.01
#> alternative hypothesis: two-sided

检验两个向量是否服从同一分布

ks.test(rnorm(20), rnorm(30))
#> 
#>  Two-sample Kolmogorov-Smirnov test
#> 
#> data:  rnorm(20) and rnorm(30)
#> D = 0.1, p-value = 1
#> alternative hypothesis: two-sided

相关性检验

cor.test(mtcars$mpg, mtcars$wt)
#> 
#>  Pearson's product-moment correlation
#> 
#> data:  mtcars$mpg and mtcars$wt
#> t = -10, df = 30, p-value = 1e-10
#> alternative hypothesis: true correlation is not equal to 0
#> 95 percent confidence interval:
#>  -0.934 -0.744
#> sample estimates:
#>    cor 
#> -0.868

一共有 3 种方法("pearson", "kendall", "spearman")

cor.test(mtcars$mpg, mtcars$wt, method = "spearman", exact = F)
#> 
#>  Spearman's rank correlation rho
#> 
#> data:  mtcars$mpg and mtcars$wt
#> S = 10292, p-value = 1e-11
#> alternative hypothesis: true rho is not equal to 0
#> sample estimates:
#>    rho 
#> -0.886

不依赖分布的检验

均值检验

Wilcoxon 检验是 t 检验的非参数版本。

单样本秩和检验

x <- c(1.83,  0.50,  1.62,  2.48, 1.68, 1.88, 1.55, 3.06, 1.30)
y <- c(0.878, 0.647, 0.598, 2.05, 1.06, 1.29, 1.06, 3.14, 1.29)
wilcox.test(x, y, paired = TRUE, alternative = "greater")
#> 
#>  Wilcoxon signed rank exact test
#> 
#> data:  x and y
#> V = 40, p-value = 0.01953
#> alternative hypothesis: true location shift is greater than 0
wilcox.test(y - x, alternative = "less")    # The same.
#> 
#>  Wilcoxon signed rank exact test
#> 
#> data:  y - x
#> V = 5, p-value = 0.01953
#> alternative hypothesis: true location is less than 0

## 使用公式进行单样本检验
depression <- data.frame(first = x, second = y, change = y - x)
wilcox.test(change ~ 1, data = depression)
#> 
#>  Wilcoxon signed rank exact test
#> 
#> data:  change
#> V = 5, p-value = 0.03906
#> alternative hypothesis: true location is not equal to 0
wilcox.test(Pair(first, second) ~ 1, data = depression)
#> 
#>  Wilcoxon signed rank exact test
#> 
#> data:  Pair(first, second)
#> V = 40, p-value = 0.03906
#> alternative hypothesis: true location shift is not equal to 0

双样本秩和检验

x <- c(0.80, 0.83, 1.89, 1.04, 1.45, 1.38, 1.91, 1.64, 0.73, 1.46)
y <- c(1.15, 0.88, 0.90, 0.74, 1.21)
wilcox.test(x, y, alternative = "g")        # greater
#> 
#>  Wilcoxon rank sum exact test
#> 
#> data:  x and y
#> W = 35, p-value = 0.1272
#> alternative hypothesis: true location shift is greater than 0
wilcox.test(x, y, alternative = "greater",
+             exact = FALSE, correct = FALSE) # H&W large sample
#> 
#>  Wilcoxon rank sum test
#> 
#> data:  x and y
#> W = 35, p-value = 0.1103
#> alternative hypothesis: true location shift is greater than 0
wilcox.test(rnorm(10), rnorm(10, 2), conf.int = TRUE)
#> 
#>  Wilcoxon rank sum exact test
#> 
#> data:  rnorm(10) and rnorm(10, 2)
#> W = 2, p-value = 4.33e-05
#> alternative hypothesis: true location shift is not equal to 0
#> 95 percent confidence interval:
#>  -4.072355 -1.718803
#> sample estimates:
#> difference in location 
#>             -2.906765 

多均值比较

多均值比较使 Kruskal-Wallis 秩和检验。

x <- c(2.9, 3.0, 2.5, 2.6, 3.2) # normal subjects
y <- c(3.8, 2.7, 4.0, 2.4)      # with obstructive airway disease
z <- c(2.8, 3.4, 3.7, 2.2, 2.0) # with asbestosis
kruskal.test(list(x, y, z))
#> 
#>  Kruskal-Wallis rank sum test
#> 
#> data:  list(x, y, z)
#> Kruskal-Wallis chi-squared = 0.77143, df = 2, p-value = 0.68

## 相同的结果
x <- c(x, y, z)
g <- factor(rep(1:3, c(5, 4, 5)),
+             labels = c("Normal subjects",
+                        "Subjects with obstructive airway disease",
+                        "Subjects with asbestosis"))
kruskal.test(x, g)
#> 
#>  Kruskal-Wallis rank sum test
#> 
#> data:  x and g
#> Kruskal-Wallis chi-squared = 0.77143, df = 2, p-value = 0.68

## 使用公式
require(graphics)
boxplot(Ozone ~ Month, data = airquality)
kruskal.test(Ozone ~ Month, data = airquality)
#> 
#>  Kruskal-Wallis rank sum test
#> 
#> data:  Ozone by Month
#> Kruskal-Wallis chi-squared = 29.267, df = 4, p-value = 6.901e-06

方差检验

使用Fligner-Killeen(中位数)检验完成不同组别的方差比较

fligner.test(InsectSprays$count, InsectSprays$spray)
#> 
#>  Fligner-Killeen test of homogeneity of variances
#> 
#> data:  InsectSprays$count and InsectSprays$spray
#> Fligner-Killeen:med chi-squared = 14.483, df = 5, p-value = 0.01282
fligner.test(count ~ spray, data = InsectSprays)
#> 
#>  Fligner-Killeen test of homogeneity of variances
#> 
#> data:  count by spray
#> Fligner-Killeen:med chi-squared = 14.483, df = 5, p-value = 0.01282

scale参数差异

R 有一些检验可以用来确定scale参数的差异。分布的scale参数确定分布函数的尺度,如 t 分布的自由度。

Ansari-Bradley 检验(两样本)

## Serum iron determination using Hyland control sera
ramsay <- c(111, 107, 100, 99, 102, 106, 109, 108, 104, 99,
            101, 96, 97, 102, 107, 113, 116, 113, 110, 98)
jung.parekh <- c(107, 108, 106, 98, 105, 103, 110, 105, 104,
            100, 96, 108, 103, 104, 114, 114, 113, 108, 106, 99)
ansari.test(ramsay, jung.parekh)
#> 
#>  Ansari-Bradley test
#> 
#> data:  ramsay and jung.parekh
#> AB = 185.5, p-value = 0.1815
#> alternative hypothesis: true ratio of scales is not equal to 1
#> 
#> Warning message:
#> In ansari.test.default(ramsay, jung.parekh) : 无法精確計算带连结的p值

Mood 检验(两样本)

mood.test(rnorm(10), rnorm(10, 0, 2), conf.int = TRUE)
#> 
#>  Mood two-sample test of scale
#> 
#> data:  rnorm(10) and rnorm(10, 0, 2)
#> Z = -1.5595, p-value = 0.1189
#> alternative hypothesis: two.sided

离散数据

比例检验

未给定概率

使用 prop.test() 比较两组观测值发生的概率是否有差异。

heads <- rbinom(1, size = 100, prob = .5)   # 二项分布:随机生成一个概率为0.5,在100次试验中能成功发生的的次数 
prop.test(heads, 100)          # continuity correction TRUE by default
#> 
#>  1-sample proportions test with continuity correction
#> 
#> data:  heads out of 100, null probability 0.5
#> X-squared = 0.25, df = 1, p-value = 0.6171
#> alternative hypothesis: true p is not equal to 0.5
#> 95 percent confidence interval:
#>  0.3703535 0.5719775
#> sample estimates:
#>    p 
#> 0.47 
prop.test(heads, 100, correct = FALSE)
#> 
#>  1-sample proportions test without continuity correction
#> 
#> data:  heads out of 100, null probability 0.5
#> X-squared = 0.36, df = 1, p-value = 0.5485
#> alternative hypothesis: true p is not equal to 0.5
#> 95 percent confidence interval:
#>  0.3751082 0.5671114
#> sample estimates:
#>    p 
#> 0.47 

## H0原假设:从患者中抽取的四个人,在人群中吸烟者的真实比例相同。
## H1备择假设:抽取的这个比例在至少一个群体中是不同的。
smokers  <- c( 83, 90, 129, 70 )
patients <- c( 86, 93, 136, 82 )
prop.test(smokers, patients)
#> 
#>  4-sample test for equality of proportions without continuity correction
#> 
#> data:  smokers out of patients
#> X-squared = 12.6, df = 3, p-value = 0.005585
#> alternative hypothesis: two.sided
#> sample estimates:
#>    prop 1    prop 2    prop 3    prop 4 
#> 0.9651163 0.9677419 0.9485294 0.8536585

给定概率值

prop.test(heads, 100, p = 0.3, correct = FALSE) 
#> 
#>  1-sample proportions test without continuity correction
#> 
#> data:  heads out of 100, null probability 0.3
#> X-squared = 13.762, df = 1, p-value = 0.0002075
#> alternative hypothesis: true p is not equal to 0.3
#> 95 percent confidence interval:
#>  0.3751082 0.5671114
#> sample estimates:
#>    p 
#> 0.47 

二项式检验

在伯努利实验中对一个简单的原假设成功的概率进行精确检验。

## 在孟德尔遗传的假设下,两种特定基因型植物的杂交产生的后代中,1/4是“矮秆”,3/4是“高秆”。为了确定这个假设是否合理,一个杂交的结果是后代有243个矮秆植株和682个高秆植株。如果“高秆”被视为成功,原假设是p = 3/4,备择假设是p != 3/4。
binom.test(c(682, 243), p = 3/4)
#> 
#>  Exact binomial test
#> 
#> data:  c(682, 243)
#> number of successes = 682, number of trials = 925, p-value = 0.4
#> alternative hypothesis: true probability of success is not equal to 0.75
#> 95 percent confidence interval:
#>  0.708 0.765
#> sample estimates:
#> probability of success 
#>                  0.737
binom.test(682, 682 + 243, p = 3/4)  # The same
#> 
#>  Exact binomial test
#> 
#> data:  682 and 682 + 243
#> number of successes = 682, number of trials = 925, p-value = 0.4
#> alternative hypothesis: true probability of success is not equal to 0.75
#> 95 percent confidence interval:
#>  0.708 0.765
#> sample estimates:
#> probability of success 
#>                  0.737

与其他的检验函数不同,这里的 p 值表示试验成功率与假设值的最小差值

列联表检验

用来确定两个分类变量是否相关。

小的列联表

当列联表较小时,Fisher 精确检验有较好的检验结果。

## 一个关于喝茶的故事
## 有一位英国妇女,她声称能分辨出是先往杯子里倒牛奶还是先往杯子里倒茶。为了测试,她喝了8杯茶,其中4杯先加了牛奶。为了确定倒酒的顺序和女士的猜测结果是否有关联?可作出如下假设:
## 原假设是真实的倒酒顺序和女士的猜测之间没有关联,备择假设是有正关联(比值比大于1)。
TeaTasting <-
     matrix(c(3, 1, 1, 3),
            nrow = 2,
            dimnames = list(Guess = c("Milk", "Tea"),
                            Truth = c("Milk", "Tea")))
fisher.test(TeaTasting, alternative = "greater")
#> 
#>  Fisher's Exact Test for Count Data
#> 
#> data:  TeaTasting
#> p-value = 0.2429
#> alternative hypothesis: true odds ratio is greater than 1
#> 95 percent confidence interval:
#>  0.3135693       Inf
#> sample estimates:
#> odds ratio 
#>   6.408309 

大的列联表

当列联表较大时,Fisher 计算量很大,可以使用卡方检验替代。

chisq.test(TeaTasting)
#> Warning in chisq.test(TeaTasting): Chi-squared approximation may be incorrect
#> 
#>  Pearson's Chi-squared test with Yates' continuity correction
#> 
#> data:  TeaTasting
#> X-squared = 0.5, df = 1, p-value = 0.5

对于三变量的混合影响,使用 Cochran-Mantel-Haenszel 检验。

在三变量没有相互作用的情况下,进行Cochran-Mantel-Haenszel卡方检验,即每个层面中有两个变量是条件独立的。

## 青霉素和兔子:立即注射或1.5小时延迟注射青霉素对家兔抗-溶血性链球菌致死效果的研究。
Rabbits <-       # 生成5个2*2的三维数组
array(c(0, 0, 6, 5,
        3, 0, 3, 6,
        6, 2, 0, 4,
        5, 6, 1, 0,
        2, 5, 0, 0),
      dim = c(2, 2, 5),
      dimnames = list(
          Delay = c("None", "1.5h"),
          Response = c("Cured", "Died"),
          Penicillin.Level = c("1/8", "1/4", "1/2", "1", "4")))
## 经典Mantel-Haenszel检验
mantelhaen.test(Rabbits)
#> 
#>  Mantel-Haenszel chi-squared test with continuity correction
#> 
#> data:  Rabbits
#> Mantel-Haenszel X-squared = 3.9286, df = 1, p-value = 0.04747
#> alternative hypothesis: true common odds ratio is not equal to 1
#> 95 percent confidence interval:
#>   1.026713 47.725133
#> sample estimates:
#> common odds ratio 
#>                 7
# P = 0.047,证据表明立即注射治愈率较高

## 精确条件检验
mantelhaen.test(Rabbits, exact = TRUE)
#> 
#>  Exact conditional test of independence in 2 x 2 x k tables
#> 
#> data:  Rabbits
#> S = 16, p-value = 0.03994
#> alternative hypothesis: true common odds ratio is not equal to 1
#> 95 percent confidence interval:
#>    1.077401 529.837399
#> sample estimates:
#> common odds ratio 
#>          10.36102 
# p = 0.040 选择的精确条件检验,直接注射治愈率较高。

 

用 McNemar 卡方检验-检验二维列联表的对称性

## 总统支持率:在相隔一个月的两项调查中,总统的执政表现得到了认可,该调查对1600名处于投票年龄的美国人进行了随机抽样调查。
Performance <-
matrix(c(794, 86, 150, 570),
       nrow = 2,
       dimnames = list("1st Survey" = c("Approve", "Disapprove"),
                       "2nd Survey" = c("Approve", "Disapprove")))
mcnemar.test(Performance)
#> 
#>  McNemar's Chi-squared test with continuity correction
#> 
#> data:  Performance
#> McNemar's chi-squared = 16.818, df = 1, p-value = 4.115e-05
# p-value = 4.115e-05 说明支持率有明显变化(实际上是下降)

列联表非参数检验

Friedman 秩和检验是一个非参数版本的双边 ANOVA 检验。(应用于未重复的块数据)

## 通过对三种方法(“四舍五入”、“窄角”和“广角”),来比较速度的差异,18名球员中,每一个人从距离本垒35英尺的一垒线上的一个点到距离二垒15英尺的一个点的两次跑垒的平均时间被记录下来。
RoundingTimes <-
matrix(c(5.40, 5.50, 5.55,
         5.85, 5.70, 5.75,
         5.20, 5.60, 5.50,
         5.55, 5.50, 5.40,
         5.90, 5.85, 5.70,
         5.45, 5.55, 5.60,
         5.40, 5.40, 5.35,
         5.45, 5.50, 5.35,
         5.25, 5.15, 5.00,
         5.85, 5.80, 5.70,
         5.25, 5.20, 5.10,
         5.65, 5.55, 5.45,
         5.60, 5.35, 5.45,
         5.05, 5.00, 4.95,
         5.50, 5.50, 5.40,
         5.45, 5.55, 5.50,
         5.55, 5.55, 5.35,
         5.45, 5.50, 5.55,
         5.50, 5.45, 5.25,
         5.65, 5.60, 5.40,
         5.70, 5.65, 5.55,
         6.30, 6.30, 6.25),
       nrow = 22,
       byrow = TRUE,
       dimnames = list(1 : 22,
                       c("Round Out", "Narrow Angle", "Wide Angle")))
friedman.test(RoundingTimes)
#> 
#>  Friedman rank sum test
#> 
#> data:  RoundingTimes
#> Friedman chi-squared = 11.143, df = 2, p-value = 0.003805
# p-value = 0.003805 强有力的证据证明这些方法在速度方面是等价的
最后编辑于
©著作权归作者所有,转载或内容合作请联系作者
禁止转载,如需转载请通过简信或评论联系作者。
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 211,884评论 6 492
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 90,347评论 3 385
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 157,435评论 0 348
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 56,509评论 1 284
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 65,611评论 6 386
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 49,837评论 1 290
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 38,987评论 3 408
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 37,730评论 0 267
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 44,194评论 1 303
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 36,525评论 2 327
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 38,664评论 1 340
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 34,334评论 4 330
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 39,944评论 3 313
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 30,764评论 0 21
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 31,997评论 1 266
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 46,389评论 2 360
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 43,554评论 2 349

推荐阅读更多精彩内容