Permutation test(置换检验)以及在R中的应用

Permutation test 可以称作是置换检验,Fisher于20世纪30年代提出的一种基于大量计算(computationally intensive),利用样本数据的全(或随机)排列,进行统计推断的方法,因其对总体分布自由,应用较为广泛,特别适用于总体分布未知的小样本资料,以及某些难以用常规方法分析资料的假设检验问题。在具体使用上它和Bootstrap Methods类似,通过对样本进行顺序上的置换,重新计算统计检验量,构造经验分布,然后在此基础上求出P-value进行推断。

总体上其实就是重新采样增加小样本的整体样本量,然后看其的概率分布来预测假设是否成立。这个基于t-test的,类似于与t-test,但是我觉的在样本量非常少的情况下,用置换检验可以更好的说明问题。这个p-value主要是sobs值在置换经验中均值的分布概率,最终也是利用p-value的值来判断假设是否成立的,看p值的大小,p值小于0.05时,是说明拒绝H0,大于0.05.则是说明服从0假设。


下面说一下,置换检验在R中应用。

1.先举R中带的例子。


排列测试在实验研究中特别相关,我们常常对治疗组之间无差异的拒绝零假设感兴趣。在这些情况下,置换检验很好的地代表了我们的推理过程,因为我们的零假设是两个治疗组在结果上没有差异(即,结果是独立于治疗分配而观察到的)。当我们在测试期间置换结果值时,因此我们看到我们可能具有的所有可能的替代治疗分配排列以及我们观察到的数据的平均差异(相对于我们可以看到的结果是独立的所有差异的治疗任务的位置)。虽然排列测试要求我们看到数据的所有可能排列(可能变得非常大),通过简单地进行大量的重采样,我们可以轻松地进行“近似置换测试”。在期望中,该过程应该近似于排列分布。假设我们的研究中有20个单位那么他的排列的数量是:

> factorial(20)

[1] 2.432902e+18

这个时候,数字完全超出了我们合理计算的数字,但是我们可以从该排列分布中随机抽样以获得近似排列分布,只需运行大量重新采样即可。让我们看一下使用一些组成数据的例子:

> set.seed(1)

> n <- 100tr <- rbinom(100, 1, 0.5)

> y <- 1 + tr + rnorm(n, 0, 3)

> diff(by(y, tr, mean))   #得到差异的结果 1.341

> s <- sample(tr, length(tr), FALSE)

> diff(by(y, s, mean))  #-0.2612

当我们重新取样而无需替换再次重新计算差异:为-0.2612

在这里,我们使用置换处理向量s计算差异并找到非常小的差异而不是用tr。如果我们重复这个过程很多次,我们可以建立我们的近似置换分布(即均值差的采样分布)。我们将使用replicatedo重复我们的排列过程。结果将是每个排列(即我们的分布)的差异向量:

dist <- replicate(2000,diff(by(y, sample(tr, length(tr), FALSE), mean)))

hist(dist, xlim = c(-3, 3), col = "black", breaks = 100)####绘制直方图

abline(v = diff(by(y, tr, mean)), col = "blue", lwd = 2)###画差异值的直线位置

然后我们可以用hist函数看一下这个分布然后画一个直方图来看它的差异


通过看这个结果,看到我们的0假设是被拒绝的,也就是说治疗组之间是有差异的。

然后我们可以量化这个结果,也就是生成一个p-value值,通过p-value来更直观的观测结果:

> sum(dist > diff(by(y, tr, mean)))/2000 # one-tailed test 单边检测## 0.009

> sum(abs(dist) > abs(diff(by(y, tr, mean))))/2000 # two-tailed test 双边检测##0.018

单边检验和双边检验的结果是不同的,但是这个结果差异不可能是由于改变治疗与结果无关的变量集中可观测变量相关的。


2:在R中我们不可能每次都要构造自己的置换检验的分布集,R那我们可以这用coin包中的independence_test函数,但是一个问题是实际问题考虑的可能只是单侧的置换检验的结果:下面是这个函数针对于上述例子的应用。

library(coin)

independence_test(y ~ tr, alternative = "greater") # one-tailed

independence_test(y ~ tr) # two-tailed    # two-tailed,默认就是

结果显示### ## Asymptotic General Independence Test

## ## data: y by tr

## Z = 2.315, p-value = 0.01029

## alternative hypothesis: greater

对于这个函数的应用,它的具体的参数及应用如下



independence_test(asat ~ group, data = asat,

                  ## exact null distribution

                  distribution = "exact",

                  ## one-sided test

                  alternative = "greater",

                  ## apply normal scores to asat$asat

                  ytrafo = function(data)

                      trafo(data, numeric_trafo = normal_trafo),

                  ## indicator matrix of 1st level of asat$group

                  xtrafo = function(data)

                      trafo(data, factor_trafo = function(x)

                          matrix(x == levels(x)[1], ncol = 1)))



除此之外,还可以用oneway_test做,效果类似。

但是oneway_test计算的默认是双侧的P值,这时的计算最好是用Deducer包中的perm.t.test函数

这个函数的可以很简单的设置是单边的还是双边的,简单的例子使用如下

library(Deducer)

x<-c(20,34,67,53,12,13,55,89)

y<-c(23,45,12,56,23,67,22,66)

perm.t.test(x,y,alternative = "greater",midp = TRUE,B = 1000) #进行1000次置换检验

下面简单介绍下perm.t.test的 参数

如下所示

perm.t.test(x,y,statistic=c("t","mean"),alternative=c("two.sided", "less", "greater"), midp=TRUE, B=10000)

Arguments

x  第一个参数向量(数字类型)

y   第二个参数向量(数字类型)

statistic 统计用的标准,t或者均值

alternative 进行统计排列的方式,主要有三种,如上面所示

midp 确定p-value是否应用

B    进行随机置换取样的次数

这就是置换检验在R中的简单介绍

©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 211,743评论 6 492
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 90,296评论 3 385
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 157,285评论 0 348
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 56,485评论 1 283
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 65,581评论 6 386
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 49,821评论 1 290
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 38,960评论 3 408
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 37,719评论 0 266
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 44,186评论 1 303
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 36,516评论 2 327
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 38,650评论 1 340
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 34,329评论 4 330
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 39,936评论 3 313
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 30,757评论 0 21
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 31,991评论 1 266
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 46,370评论 2 360
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 43,527评论 2 349

推荐阅读更多精彩内容