原理参考 文章 ，主要思想我认为是求出所有分布的可能（中间的一般为零假设），出现这种分布的概率。

distribution= 参数可为exact（精确模式，即依据所有可能的排列组合，仅适用于两样本问题）、approxiamate(nresample=#)（蒙特卡洛抽样，#指需要重复的次数）、asymptotic（渐进分布抽样）

lmPerm包更擅长方差分析。示例实验设计仍为5组接受不同治疗方法的多组结果比较。

实验示例仍为关节炎的治疗（两种）与效果（无、部分、显著）间的关系

实验示例为研究文盲率与谋杀率是否相关

主要为 lmp() 、 aovp() 两个函数分别对应参数法的 lm() 线性回归、 aov() 方差分析。主要格式上的区别是添加了 perm= 参数。可以为Exact（精确模式）、Prob（不断从可能的序列中抽样，直至估计的标准差在估计的p值0.1之下）、SPR（使用贯序概率比检验来判断何时停止抽样）。值得注意的是当样本观测大于10，perm="Exact"自动默认转为"Prob"，因此精确检验只适用于小样本问题。

（1）简单线性回归

实验示例仍为以身高预测体重的设计

（2）多项式回归

高精度拟合身高体重回归关系

（3）多元回归

探究谋杀率与多因素的回归关系

（1）单因素方差分析

（2）单因素协方差分析

实验示例仍为药物对刚出生小鼠体重影响，协变量为怀孕时间

（3）双因素方差分析（交互效应）

实验示例：两种药物分别在不同剂量下对小鼠牙齿长度的影响。

核心思想是有放回的抽样多次（1000次）

（1）写一个能返回带研究统计量的函数；

（2）确定重复数，使用 boot() 函数处理；（一般重复1000次即可；此外有人认为初始样本大小为20-30即可得到足够好的结果）；

（3） boot.ci() 函数计算统计量置信区间。

实验示例：使用mtcar数据框，采用多元回归，根据车重和发动机排量来预测汽车的每加仑行驶的英里数。想获得95%的R平方值（预测变量对响应变量可解释的方差比）的置信区间

（1）首先写函数

（2）然后使用boot()函数

（3）最后boot.ci()函数求置信区间

实验示例：使用mtcar数据框，采用多元回归，根据车总和发动机排量来预测汽车的每加仑行驶的英里数。想获取一个统计量向量--三个回归系数（截距项、车总、发动机排量）95%的置信区间。

数据准备

许多实际情况中统计假设（假定观测数据抽样自正态分布或者其他性质较好的理论分布）并不一定满足，比如数据抽样于未知或混合分布、样本量过小、存在离群点、基于理论分布设计合适的统计检验过于复杂且数学上难以处理等情况，这时基于随机化和重抽样的统计方法就可派上用场。

置换检验的定义

置换检验(Permutation test)，也称随机化检验或重随机化检验，是Fisher于20世纪30年代提出的一种基于大量计算（computationally intensive），利用样本数据的全（或随机）排列，进行统计推断的方法，因其对总体分布自由，应用较为广泛，特别适用于总体分布未知的小样本资料，以及某些难以用常规方法分析资料的假设检验问题。

置换检验的原理

1、提出原假设，比如XX处理后结果没有变化

2、计算统计量，如两组的均值之差，记作t0

3、将所有样本放在一起，然后随机排序进行分组，再计算其统计量t1

4、重复第3步骤，直至所有排序可能性都齐全（比如有A组有n样本，B组有m样本，则总重复次数相当于从n+m中随机抽取n个的次数），得到一系列的统计量（t1-tn）

5、最后将这些统计量按照从小到大排序，构成抽样分布，再看t0是否落在分布的置信区间内（如95%置信区间），这时候可计算一个P值（如果抽样总体1000次统计量中大于t0的有10个，则估计的P值为10/1000=0.01），落在置信区间外则拒绝原假设

6、如果第3步骤是将所有可能性都计算了的话，则是精确检验；如果只取了计算了部分组合，则是近似结果，这时一般用蒙特卡罗模拟（Monte Carlo simulation）的方法进行置换检验

7、置换检验和参数检验都计算了统计量，但是前者是跟置换观测数据后获得的经验分布进行比较，后者则是跟理论分布进行比较。

请牢记：置换检验都是使用伪随机数来从所有可能的排列组合中进行抽样（当做近似检验时）。因此，每次检验的结果都有所不同。

coin包提供了一个进行置换检验的一般性框架。通过该包，你可以回答如下问题。

响应值与组的分配独立吗？

两个数值变量独立吗？

两个类别型变量独立吗？

表12-2列出来的每个函数都是如下形式：

function_name(formula, data, distribution=)

其中：

 formula描述的是要检验变量间的关系。示例可参见表12-2；

 data是一个数据框；

 distribution指定经验分布在零假设条件下的形式，可能值有exact，asymptotic和

approximate。

若distribution = "exact"，那么在零假设条件下，分布的计算是精确的（即依据所有可能的排列组合）。当然，也可以根据它的渐进分布（distribution = "asymptotic"）或蒙特卡洛重抽样（distribution = "approxiamate(B = #)"）来做近似计算，其中#指所需重复的次数。

distribution = "exact"当前仅可用于两样本问题。

传统t检验表明存在显著性差异（p <0.05），而精确检验却表明差异并不显著（p >0.072）。

第7章我用自己的数据进行了t检验，对比一下传统t检验和置换检验，结果如下：

两种检验方式下结果都是显著的

Wilcoxon-Mann-Whitney U检验

coin包规定所有的类别型变量都必须以因子形式编码。

wilcox.test()默认计算的也是精确分布。

K样本检验的置换检验

通过chisq_test()或cmh_test()函数，我们可用置换检验判断两类别型变量的独立性。 当数据可根据第三个类别型变量进行分层时，需要使用后一个函数。若变量都是有序型，可使用lbl_test()函数来检验是否存在线性趋势。

卡方独立性检验

卡方独立性检验的置换检验

你可能会有疑问，为什么需要把变量Improved从一个有序因子变成一个分类因子？（好问题！）这是因为，如果你用有序因子，coin()将会生成一个线性与线性趋势检验，而不是卡方检验。

结果解读：两种检验下p值都是小于0.05，说明Treatment和Improved之间相互不独立

自己数据的演示

结果解读：p值均为1，表明nitrogen和variety相互独立。

spearman_test()函数提供了两数值变量的独立性置换检验。

当处于不同组的观测已经被分配得当，或者使用了重复测量时，样本相关检验便可派上用场。

对于两配对组的置换检验，可使用wilcoxsign_test()函数；多于两组时，使用friedman_test()函数。

自己数据演示

lmPerm包可做线性模型的置换检验。比如lmp()和aovp()函数即lm()和aov()函数的修改版，能够进行置换检验，而非正态理论检验。

lmp()和aovp()函数的参数与lm()和aov()函数类似，只额外添加了perm =参数。

perm =选项的可选值有"Exact"、"Prob"或"SPR"。Exact根据所有可能的排列组合生成精确检验。Prob从所有可能的排列中不断抽样，直至估计的标准差在估计的p值0.1之下，判停准则由可选的Ca参数控制。SPR使用贯序概率比检验来判断何时停止抽样。注意，若观测数大于10，perm = "Exact"将自动默认转为perm = "Prob"，因为精确检验只适用于小样本问题。

简单线性回归的置换检验

R语言实战的例子：

多项式回归的置换检验

R语言实战的例子：

自己数据集的例子：

R语言实战的例子：

自己数据集的例子：

当两种方法所得结果不一致时，你需要更加谨慎地审视数据，这很可能是因为违反了正态性假设或者存在离群点。

R语言实战的例子：

自己数据集的例子：

R语言实战的例子：

自己数据集的例子：

R语言实战的例子：

自己数据集的例子：

值得注意的是，当将aovp()应用到方差分析设计中时，它默认使用唯一平方和法（SAS也称为类型III平方和）。每种效应都会依据其他效应做相应调整。R中默认的参数化方差分析设计使用的是序贯平方和（SAS是类型I平方和）。每种效应依据模型中先出现的效应做相应调整。对于平衡设计，两种方法结果相同，但是对于每个单元格观测数不同的不平衡设计，两种方法结果则不同。不平衡性越大，结果分歧越大。若在aovp()函数中设定seqs = TRUE，可以生成你想要的序贯平方和。

你可能已经注意到，基于正态理论的检验与上面置换检验的结果非常接近。在这些问题中数据表现非常好，两种方法结果的一致性也验证了正态理论方法适用于上述示例。

当然，置换检验真正发挥功用的地方是处理非正态数据（如分布偏倚很大）、存在离群点、样本很小或无法做参数检验等情况。不过，如果初始样本对感兴趣的总体情况代表性很差，即使是置换检验也无法提高推断效果。

置换检验主要用于生成检验零假设的p值，它有助于回答“效应是否存在”这样的问题。不过，置换方法对于获取置信区间和估计测量精度是比较困难的。幸运的是，这正是自助法大显神通的地方。

所谓自助法，即从初始样本重复随机替换抽样，生成一个或一系列待检验统计量的经验分布。 无需假设一个特定的理论分布，便可生成统计量的置信区间，并能检验统计假设。

倘若你假设均值的样本分布不是正态分布，该怎么办呢？可使用自助法。

(1)从样本中随机选择10个观测，抽样后再放回。有些观测可能会被选择多次，有些可能一直都不会被选中。

(2)计算并记录样本均值。

(3)重复1和2一千次。

(4)将1000个样本均值从小到大排序。

(5)找出样本均值2.5%和97.5%的分位点。此时即初始位置和最末位置的第25个数，它们就限定了95%的置信区间。

boot包扩展了自助法和重抽样的相关用途。你可以对一个统计量(如中位数)或一个统计量向量(如一列回归系数)使用自助法。

一般来说，自助法有三个主要步骤。

(1)写一个能返回待研究统计量值的函数。如果只有单个统计量（如中位数），函数应该返回一个数值；如果有一列统计量（如一列回归系数），函数应该返回一个向量。

(2)为生成R中自助法所需的有效统计量重复数，使用boot()函数对上面所写的函数进行处理。

(3)使用boot.ci()函数获取第(2)步生成的统计量的置信区间。

主要的自助法函数是boot()，它的格式为：

bootobject <- boot(data=, statistic=, R=, ...)

参数见下表：

boot()函数调用统计量函数R次，每次都从整数1:nrow(data)中生成一列有放回的随机指标，这些指标被统计量函数用来选择样本。统计量将根据所选样本进行计算，结果存储在bootobject中。

你可以用bootobject t0和bootobject t来获取这些元素。

一旦生成了自助样本，可通过print()和plot()来检查结果。如果结果看起来还算合理， 使用boot.ci()函数获取统计量的置信区间。格式如下:

boot.ci(bootobject, conf=, type= )

type参数设定了获取置信区间的方法。perc方法（分位数）展示的是样本均值，bca将根据偏差对区间做简单调整。

回归的R平方值

1000次自助抽样

输出结果

结果可视化

95%的置信区间获取

回归系数向量函数

自助抽样1000次

获得车重和发动机排量95%的置信区间

置换检验和自助法并不是万能的，它们无法将烂数据转化为好数据。当初始样本对于总体情况的代表性不佳，或者样本量过小而无法准确地反映总体情况，这些方法也是爱莫能助。

参考资料：

R语言中，和排序相关的常用函数有： order() ， sort() ， rank() ，一般是对向量进行操作，也可以对数据框的列进行操作。

1. order(..., decreasing = FALSE)

“...” 中可以是一个向量（数值型，字符型，逻辑型， 因子型 均可），也可以是多个向量（ 长度必须相同 ）

它返回的是排列后（默认是升序）的元素在该向量中所处的位置，即 索引 ，所以返回的不是原来向量的那些数值，而是排序后那些数值所对应的位置。它在常在数据框中运用，可以根据某一列和某几列来调整数据框。

1.1 单个向量

1.2 多个向量（数据框的多列）

上面最左边的一列数值（不是var1这一列）就是 order(df$var1) 生成的，列var1的顺序就是order后的，注意和前面df 的列var1进行比较，还有一点要牢记的就是这种在数据框里的调整，是 整行变动(都按列var1来) 。故而本质上，对于数据框而言， order函数出来的是原始数据框中的行号，行号顺序一变，意味着行号代表的整行跟着变 。

接着再按列var2排序（注意是在列var1已经排好的基础上再按列var2重新排序，即此时列var1里的 非重复值的行顺序不会再变了，只有列var1里有重复值（数值相同的）的行才会变换 ，而且是按照列var2来变换，。如还要按其它列再排，以此类推）：

最后再按列var3排序，此时只有列var2中 有重复值的行 （当然这时var1肯定也是重复的）才会变换------这里是含有9的那两行，并且是按照这些重复值对应到列var3的那些行的数值（数值1和2）来变的

总结提升，order()函数中，如果第一个向量（或者说是数据框里你想要根据它来调整的那一列）里没有重复值，那么按照后面的所有向量（不管有多少个）的重排都不会进行（或者说重排的结果不会变），order(vec1,vec2,vec3,...) 返回的行号及其顺序和 order(vec1) 的是一样。

同理，也可以处理对应的行（比如矩阵或是数据框的行）

1.3 factor（因子型向量）

2. sort(x, decreasing = FALSE)

x 只能是 一个 向量（数值型，字符型，逻辑型， 因子型 均可）

返回的是排序后（默认是升序）的那个数值向量（ 还是那些数值，只不过是排序过了的 ）

3. rank(x, na.last = TRUE)

求秩的函数，x 只能是 一个**向量（数值型，字符型，逻辑型）， 该向量一般不会有重复值 ，返回的是该向量中对应元素的“排名”，即元素顺序它不会改变，只是告诉你每个元素在整个向量中的名次（如果要排序（默认是升序）的话）。

如果向量有重复值，出现的结果会有些不太好处理

匹配两个向量，返回的是第一个向量 x中的各个元素在第二个向量 y中所匹配的元素的位置值（索引，下标值），即 返回的是第二个向量的下标值组成的向量 。

注意事项：

1. 返回的下标值向量的长度与第一个向量相等，即 length(x) == length(match(x,y)) 为 TRUE。

2. 第一个向量可以是只有一个元素的向量。

3. 两个向量的长度不一定要相等。

4. 返回的是 x中每个元素在y中的位置，可用来提取y中的元素，没有匹配上的会返回 NA。

拓展用法：

y [match(x, y)] ： x中 和y 相同的元素都是哪些，前提是x中的元素在y中都有（即x是y的子集），否则会返回 NA。（ 最好是用x[ x %in y] ）

y[-match(x, y)] ： 找出 y里面有 x里面没有的数值，前提是x中的元素在y中都有（即x是y的子集），否则会报错。

当两个向量类型不一样时, match函数会进行类型转换，然后再匹配。

x %in% y

判断x中的元素是否都在y中，返回的是个长度和 x 一样的逻辑值向量，存在的话返回TRUE，否则返回FALSE。

sum(x %in% y) ： 统计 x 中有多少个元素在y中，或者说 x 和 y 有几个相同元素

x[x %in% y] ：x中 和y 相同的元素都是哪些

---