R笔记：两独立样本的差异显著性检验

在生物医学研究中最常见的行为就是对两个组进行差异性比较，两个组的比较又可以分为两独立组间的比较和两相关组（配对组）的比较。两个独立组比较用什么方法？这需要看要分析的数据类型。

如果是两组定量资料间的比较，我们首先考虑的应该是两独立样本的t检验；如果是两组结局变量无序的分类资料间的比较，我们首先考虑的应该是卡方检验；如果两组定量资料不满足使用t检验的条件，或者结局变量为有序多分类的资料，我们可能要考虑使用非参数检验。

【1】两独立样本的t检验：正态、方差齐或不齐

t.test {stats}：Student's t-Test，Performs one and two sample t-tests on vectors of data.

• t.test(x, y = NULL,alternative = c("two.sided", "less", "greater"),mu = 0, paired = FALSE, var.equal = FALSE,conf.level = 0.95, ...)

• t.test(formula, data, subset, na.action, ...)

①数据载入

library(openxlsx)  #调用程序包openxlsx

ttest<-read.xlsx("D:/Temp/ints.xlsx",1)  #从名称为ints的excel文件中导入第1个sheet的数据到ttest数据框

②正态分布检验，本例使用shapiro检验，其他方法可参见<<正态分布的检验>>一文

library(stats)  #调用程序包stats。Stats默认载入，不进行专门的命令载入也可以

by(ttest$血糖下降值,ttest$分组,shapiro.test)  #对ttest数据框中的“血糖下降值”按“分组”变量分组，分别进行shapiro检验后输出。等同于shapiro.test(ttest[1:20,2]);shapiro.test(ttest[21:40,2])

qqnorm(ttest[1:20,2]);qqline(ttest[1:20,2], col = 3) #通过Q-Q_图判定试验组的正态性，并向Q-Q图中添加参考直线，直线颜色为3（蓝色）

qqnorm(ttest[21:40,2]);qqline(ttest[21:40,2], col = 4) #通过Q-Q_图判定对照组的正态性，并向Q-Q图中添加参考直线，直线颜色为4（绿色）

③方差齐性检验，本例采用leveneTest{car}进行检验，其他方法可参见<<方差齐性的检验>>一文。

library(car);library(carData)  #调用程序包car和carDate

leveneTest(血糖下降值~分组,data=ttest,center=mean)  #对数据框ttest中的“血糖下降值”按“分组“”变量分组，并检验两组的均值的方差齐性。采用的是方差分析的模型表达式形式

var.test(ttest[1:20,2],ttest[21:40,2]) #对数据框ttest中的第2列的1-20和第2列的21-40行进行方差齐性F检验（var.test {stats}）。采用的是数据直接比较形式

④两独立样本的t检验：两独立样本的t检验（正态、方差齐）

可以采用两组数据直接比较的输入形式或者模型公式的形式。

t.test(ttest[1:20,2],ttest[21:40,2],var.equal=T)

t.test(ttest$血糖下降值~ttest$分组,var.equal=T)

t.test(血糖下降值~分组,data=ttest,var.equal=T)

我们在载入数据后，先把要比较的数据赋值到一个新对象，后面的各种检验书写会变的简便易懂，如本例在数据载入后先进行变量赋值，后面的正态性、方差齐性和t检验输入上会简便一些。命令清单如下：

x<-ttest[1:40,1]

y<-ttest[1:40,2]

y1<-ttest[1:20,2]

y2<-ttest[21:40,2]

by(y,x,shapiro.test)  #shapiro正态性检验，等同于shapiro.test(y1);shapiro.test(y2)  

var.test(y1,y2)  #方差齐性F检验

library(car);library(carData) 

leveneTest(y~x,data=ttest)  #方差齐性lebene检验

t.test(y1,y2,var.equal=T)  #两独立样本的t检验

④两独立样本的t检验：两独立样本的t检验（正态、方差不齐）

方差不齐时，需要使用校正的t检验，t.test {stats}默认的是方差不同时的Welch双样本t检验，非齐性时只要不设定var.equal=TRUE就可以了。本例如果方差不齐，数据载入后命令为：

t.test(血糖下降值~分组,data=ttest)

结果如下：

【2.1】非参数检验：非正态分布和/或方差不齐

实际上当数据满足正态分布或方差齐性时，也是可以使用非参数检验来进行分析的，只不过t检验更容易发现存在的差异而已。两独立样本的非参数检验可以采用wilcox.test {stats}。

wilcox.test {stats}：Wilcoxon Rank Sum and Signed Rank Tests，Performs one- and two-sample Wilcoxon tests on vectors of data; the latter is also known as ‘Mann-Whitney’ test.【在R中Mann-Whitney U等同于Wilcoxon Rank Sum检验】

• wilcox.test(x, y = NULL,alternative = c("two.sided", "less", "greater"),    mu = 0, paired = FALSE, exact = NULL, correct = TRUE,conf.int = FALSE, conf.level = 0.95, ...)

• wilcox.test(formula, data, subset, na.action, ...)

假如示例1的数据不满足正态分布或方差齐性，数据载入后采用非参数检验命令如下：

wilcox.test(ttest$血糖下降值~ttest$分组)

或

wilcox.test(血糖下降值~分组,data=ttest)

我们可以看到结果出现了一个警示信息（Warning message）：In wilcox.test.default(x = c(……) 无法精確計算带连结的p值，意思就是我们分析的数据中有重复数据（打结现象），相同的数据在排秩的时候会按照升幂排列后位置的平均值来计算，“结”多了，零分布的大样本公式就不准了，因此要做修正，修正方法也不难，但我不明白R中为什么不给出这个修正结果，而是只是给出个近似结果。wilcox.test {stats}的解释文件中表明”By default (if exact is not specified), an exact p-value is computed if the samples contain less than 50 finite values and there are no ties. Otherwise, a normal approximation is used”。如果你对出现这个信息恶心，可以做把命令改为：wilcox.test(血糖下降值~分组,data=ttest,exact =F)，这样这条恶心人的警示就没有了，但是没有了不代表P值就是修正过的，结果还是近似估计，只是阿Q精神而已！

【2.2】非参数检验：结局变量为等级资料的情况，包括分组变量无序、结局变量有序的行列表资料；双向有序属性不同行列表资料。

示例：药物A、B治疗某种疾病，基线无统计学差异，治疗3月后根据某指标判断无效、改善和显著改善的例数如下表。两药物的疗效是否有统计学差异。

该例使用卡方检验是不恰当的，我们在<<分类资料的差异比较常用统计方法选择及SPSS实现>>一文中做过演示，操作也很简单。该例结局变量为有序多分类，应该使用非参数检验，如果使用卡方会得出相反的统计学结论。比较遗憾的是在R语言中wilcox.test {stats}的参数中并没有频数加权的可用项，所以如果你的数据已经概括成表格化（如示例表格），这类结局变量为等级资料的非参数检验反而不简便。根据秩和检验的原料，笔者认为R可以对原始格式的数据进行Wilcoxon检验，这里说的“原始数据格式”是指按照每条记录占一行的原始格式。因此如果你要在R中实现结局变量为等级资料的Wilcoxon秩和检验，需要将表格化格式还原成原始数据格式，可能有好的方法来转换，但笔者尚未找到，笔者就用笨方法在excel表格中实现啦。可能还是更好的程序包可以直接对表格化的数据进行加权，然后进行非参分析，只是笔者尚未找到而已。

命令清单如下：

library(openxlsx)

nmpt<-read.xlsx("D:/Temp/ints.xlsx",4) #从名称为ints的excel文件中导入第4个sheet的数据到nmpt数据框

wilcox.test(effect~group,data=nmpt)

如果本例采用了卡方检验，结果如下：χ2=5.244，P=0.073>0.05，两组疗效没有统计学差异。

【3】卡方分析：四表格资料；双向无序的行列表资料；分组变量有序、结局变量无序的行列表资料

胞磷胆碱与神经节苷酯治疗脑血管疾病的有效率是否有差别？

①数据录入

#方法1

group<-c("胞磷胆碱","胞磷胆碱","神经节苷酯","神经节苷酯")

effect<-c("有效","无效","有效","无效")

fre<- c(46,6,18,8) 

chist <- data.frame(group, effect,fre) 

#方法2：手动输入

ds<-data.frame(group=character(0),effect=character(0),fre=numeric(0))  #创建一个不含任何实际数据数据框ds，该数据库有3个变量：group、effect和fre，数据类型分别是字符型、字符型和数值型

chist<-edit(ds)  #fix(ds)

在打开的表格中手动输入数据即可

library(openxlsx)

chist<-read.xlsx("D:/Temp/ints.xlsx",5)  #导入表格化数据

chist2<-read.xlsx("D:/Temp/ints.xlsx",6)  #导入原始格式数据，从名称为ints的excel文件中导入第6个sheet的数据到chist2数据框

Chi<-table(chist2$group,chist2$effect)

Chi<-xtabs(~group+effect,data=chist2)  #要分析的变量应在公式的右侧（即~符号的右方），第一个变量为行变量，第二个变量为列变量，以+作为分隔符

Chi<-xtabs(fre~group+effect,data=chist)  #如果导入的数据不是原始数据格式，而是带频数的概括化表格，频数变量写在公式的左侧。要分析的变量应在公式的右侧（即~符号的右方），第一个变量为行变量，第二个变量为列变量，以+作为分隔符

library(stats) 

chisq.test(Chi)

或者采用Fisher检验。

library(stats) 

fisher.test(Chi)

END