Vol.1
方法介绍
了解什么是置换检验，以及我们为什么要使用置换检验这是我们学会使用置换检验的前提。在我们最常见的统计中，似乎不太涉及到置换检验（其实置换检验也很常见，杠精勿喷），大家常用的统计方法，例如T检验，方差分析等，这些都属于参数检验的范畴（参数检验：假设数据符合特定分布（如正态分布），适用于比例或均值比较）。而置换检验是一种非参数检验方法（非参数检验：不要求数据符合特定分布，适用于中位数或秩次比较，通常用于样本量小或数据不满足正态性假设的情况），用于比较两个或多个样本之间的差异。它通过对数据进行多次随机重新排列（或置换），从而生成一个分布来评估观察到的统计量的显著性。这种方法尤其适用于样本量较小或假设分布不明确的情况下。
置换检验能干什么呢？
其实一种研究中非常常见的情况就是我们在一批数据中分别拟合了男性和女性关于变量x和y的回归系数，但我们想要比较男性和女性群体中的回归系数是否存在显著性差异，然而遗憾的是，我们并没有另外的样本用于生成新的回归系数，然后进行T检验。因此这种思路其实是行不通的。那么研究中遇到这种情况最常用的两种方法就是Bootstrap test和Permutation test。关于Bootstrap test之前的推文有提及，这里不做赘述。此时，置换检验就是一种适合的方法。置换检验通过将男性和女性全部混合在一起，然后以一种不放回的重采样方式重新抽取与原始男性群体和女性群体相同数量的样本，将这一过程重复2000余次。每一次的重抽取就形成了一个对应的男性和女性新样本，然后就可以重新计算x与y的回归系数统计量了。上述随机置换步骤重复2000次就可以得到回归系数的经验分布，计算回归系数经验分布中大于原始回归系数的个数（计为n），那么P=n/随即置换次数。
Vol.2
如何实现
set.seed(1)
n <- 100
tr <- rbinom(100, 1, 0.5)
#这行代码使用rbinom()函数生成随机数据。具体来说：rbinom(n, size, prob)函数生成一个包含n个随机数的向量，这些随机数服从参数为size和prob的二项分布。在这行代码中，n设为100，size设为1，prob设为0.5。这意味着将生成100个伯努利试验（每次试验的结果为0或1）的结果，成功概率为0.5。所以，tr是一个长度为100的二进制向量，其中每个元素的值是0或1，均匀分布。
y <- 1 + tr + rnorm(n, 0, 3)
#1 + tr：首先，将tr中的每个元素加上1。由于tr的值只能是0或1，因此1 + tr的结果将是1（当tr为0时）或2（当tr为1时）。rnorm(n, 0, 3)：这部分代码生成一个长度为n的随机数向量，这些随机数服从正态分布，均值为0，标准差为3。这意味着生成的数值将在0附近分布，但可能会有较大的离散性。1 + tr + rnorm(n, 0, 3)：将1 + tr的结果与这个正态分布的随机数相加，形成y。如此，最终的y向量是由基于tr的结果和一些随机噪声组合而成的，因此y的取值将呈现出既受tr影响，又带有随机性的特征。综合起来，这几行代码生成了一个二元变量tr和一个受tr影响并加入了正态随机噪声的因变量y。这样可以用在各种模拟和统计建模之中。
diff(by(y, tr, mean))
#通过by(y, tr, mean)计算变量y分别在tr=0和tr=1的均值。然后，使用diff()函数来计算这两个均值之间的差值。 最终的结果是一个标量值，表示y在两组（tr=0 和 tr=1）之间的均值差。这通常用于检验两组数据的平均差异，并为后续的统计分析提供信息。
s <- sample(tr, length(tr), FALSE)
#这行代码用于创建一个随机的样本 s，其大小与 tr 一致。FALSE这个参数指定抽样方式为“无放回抽样”。s 将成为一个包含从 tr 中随机抽取的元素的向量，元素顺序会被打乱，并且每个元素只出现一次。由于 tr 只有0和1两个值，因此 s 也将包含0和1，但顺序不同。
diff(by(y, s, mean))
#这部分代码的作用是计算基于随机重抽样后的 s 的均值差异。
#在上述过程中，我们使用置换处理向量s而不是tr来计算差异，并找到一个非常小的差异。如果我们重复这个过程很多次，我们可以建立我们的近似置换分布（即平均差的采样分布）。我们将使用replicate函数来重复我们的置换过程。结果将是每个置换的差异向量（即我们的分布）：
dist <- replicate(2000, diff(by(y, sample(tr, length(tr), FALSE), mean)))
#这句代码使用replicate函数将diff(by(y, sample(tr, length(tr), FALSE), mean))这一过程重复了2000次，具体来说就是重复2000次计算基于随机重抽样后的 s 的均值差异，如此，就形成了一个关于均值差的分布了。
hist(dist, xlim = c(-3, 3), col = "black", breaks = 100)
abline(v = diff(by(y, tr, mean)), col = "red", lwd = 2)
#上面两句代码第一句是用于将dist用直方图表现出来，也就是相对于绘制出了dist的分布，其中横坐标是均值差的值，纵坐标是频数。而下面这句代码则是在绘制出分布图的基础上，在tr原始均值差的位置处画一条垂直于横坐标的红线。如下图所示：
#最后，其实从图中肉眼来看的花也会觉得tr=1和tr=0两组之间存在差异，但我们需要计算出P值。根据置换检验中P值的计算公式，计算dist中大于原始均值差的个数（计为n），那么P=n/随即置换次数。
sum(dist > diff(by(y, tr, mean)))/2000 # 单尾检验，P=0.009
sum(abs(dist) > abs(diff(by(y, tr, mean))))/2000 # 双尾检验，P=0.018
##############正式置换检验###################
#尽管上述代码可以实现置换检验的功能，但是我们没有必要每次都自行构建分布然后计算P值，幸运的是，我们可以站在前人的肩膀上去进行研究。我们可以通过Coin这个R包实现上述功能。
library(coin)
independence_test(y ~ tr, alternative = "greater") #单尾检验
#independence_test函数用于进行独立性检验，特别在两组数据中比较变量之间是否存在显著差异。y ~ tr：这是一个公式，表示通过tr来解释或预测y。在这个公式中，y是响应变量，而tr是分组变量或自变量。alternative = "greater"：这个参数指定了检验的性质。在这里，使用"greater"表示我们假设y在组tr=1中的均值大于tr=0中的均值，进行单尾检验。
independence_test(y ~ tr) # 双尾检验
参考文献
Good, P. (2013). Permutation tests: a practical guide to resampling methods for testing hypotheses. Springer Science & Business Media.
https://thomasleeper.com/Rcourse/Tutorials/
PSYCH统计实验室
通知公告
网络分析课程目前开放视频课啦！
单次课200元/讲（学生），250元/讲（非学生）
共有四讲内容：
①横断面网络分析简介与基础
②网络分析与因子分析
③交叉滞后网络分析
④时间序列网络分析
购买后开放视频权限14天，可多次申请。
并赠送所有课程相关资料（无PPT）
如果想申请购买，请联系M18812507626
更多资讯
关注我们
文稿：摘星
排版：Peruere
责编：Wink
审核：摘星
本文由“Psych统计自习室”课题组原创，欢迎转发至朋友圈。如需转载请联系后台，征得作者同意后方可转载。
