美文网首页Data Analysis for the Life Sciences
DALS006-统计推断(Inference)05-置换检验(P

DALS006-统计推断(Inference)05-置换检验(P

作者: backup备份 | 来源:发表于2019-08-16 09:27 被阅读0次

    title: DALS006-统计推断(Inference)05-置换检验(Permutation Test)
    date: 2019-08-06 12:0:00
    type: "tags"
    tags:

    • 统计推断
      categories:
    • 生物统计

    置换检验(Permutation Test)

    现在考虑一种情况,在这种情况下,标准的数理统计近拟(standard mathematical statistical approximation)都不适用。但是,我们已经计算出了一个总结性的统计,例如均值的差值,但是,我们的这些数据不符合一些即定分布,也没有办法得到总体,也就是说无法进行一些数据模拟,例如Monte Carlo模拟,面对这种情况,我们可以采用转换检验(Permutation test)。

    转换检验过程

    现在来看一下转换检验,如下所示:

    library(dplyr)
    dir <- system.file(package = "dagdata")
    filename <- file.path(dir,"extdata/femaleMiceWeights.csv")
    dat <- read.csv(filename)
    control <- filter(dat, Diet == "chow") %>% dplyr::select(Bodyweight) %>% unlist
    treatment <- filter(dat, Diet == "hf") %>% dplyr::select(Bodyweight) %>% unlist
    obsdiff <- mean(treatment) - mean(control)
    

    在前面部分中,我们描述了参数模拟,这种方法有助于我们计算出观察到的差值是否有统计学的意义。而置换检验则会使用另外一种计算方法,在这种方法中,我们会随机置换干预组与对照组的分组标签。当我们置换了干预组与对照组的分级标签后,得到的新数据我们会假设它们服从零分布,现在我们通过置换这两组数据来计算一下零分布(置换1000次),如下所示:

    N <- 12
    avgdiff <- replicate(1000, {
      all <- sample(c(control, treatment))
      newcontrols <- all[1:N]
      newtreatments <- all[(N+1):(2*N)]
      return(mean(newtreatments)-mean(newcontrols))
    })
    hist(avgdiff)
    abline(v=obsdiff,col="red",lwd=2)
    
    image

    现在我们考虑一下,有多少个零均值大于观察值呢?这个比值就是零假设成立下的p值。我们向分子与分母都加1,用于校正一下p值(书中列出的参考文献,这里不解释),如下所示:

    # the proportion of permutations with larger difference
    (sum(abs(avgdiff) > abs(obsdiff))+1)/(length(avgdiff)+1)
    

    计算结果如下所示:

    > # the proportion of permutations with larger difference
    > (sum(abs(avgdiff) > abs(obsdiff))+1)/(length(avgdiff)+1)
    [1] 0.04795205
    

    现在我们使用小样本数据来重复一下这个实验,我们通过抽样的方式来创建一个小样本数据,如下所示:

    N <- 5
    control <- sample(control, N)
    treatment <- sample(treatment, N)
    obsdiff <- mean(treatment) - mean(control)
    
    avgdiff <- replicate(1000,{
      all <- sample(c(control, treatment))
      newcontrols <- all[1:N]
      newtreatments <- all[(N+1):(2*N)]
      return(mean(newtreatments) - mean(newcontrols))
    })
    hist(avgdiff)
    abline(v=obsdiff, col="red",lwd=2)
    
    image

    上面的直方图显示的是使用置换检验计算出来的小样本数据的均值之差,红线表示的是观察值的差异。

    现在我们发现,使用这种方法得到的观察差值并不显示,这里需要记住,没有理论能保证,置换的结果接近于实际的零分布。例如,如果不同的总体之间确实存在着差异,一些置换检验的结果可能就计算不出差异,也有一些能计算出差异。这说明,在置换检验中,零分布比实际的零分布有着一个更扁平的尾部。这也就说明了,为什么置换检验会得到一个比较保守的p值。出于这种原因,当我们的样本数目比较少时通常不做置换检验。

    练习

    P91

    相关文章

      网友评论

        本文标题:DALS006-统计推断(Inference)05-置换检验(P

        本文链接:https://www.haomeiwen.com/subject/nkaxsctx.html