做统计学分析的其中一个重要目的就是寻找组间差异，在研究中，我们最常关注的问题莫过于处理组与对照组是否存在了显著不同。就两组间比较而言，t检验是常见的分析方法之一。
本文简介怎样在R中进行t检验，以实现两组间差异分析。
示例数据、R脚本等，已上传至百度盘（提取码6268）：
示例文件简要
文件“alpha.txt”为某16S测序所获得的alpha多样性指数（部分），其内容展示如下。
sample，样本名称；observed_species和shannon分别为两种类型的alpha多样性指数。
文件“group.txt”为各样本分组信息，其内容展示如下。
第一列（sample）为各样本名称；第二列（group）为各样本的分组信息。
接下来，我们期望在R中运行t检验，以查看不同分组间（两两分组之间）的各alpha多样性指数是否存在显著不同。
数据读取及正态性假设检验
首先将上述两个数据表读入R中，并合并在一起。
library(reshape2)#排列数据使用#读入文件alpha<-read.table('alpha.txt',sep='\t',header=TRUE,stringsAsFactors=FALSE,check.names=FALSE)group<-read.table('group.txt',sep='\t',header=TRUE,stringsAsFactors=FALSE,check.names=FALSE)alpha<-melt(merge(alpha,group,by='sample'),id=c('sample','group'))此时查看合并后的数据框“alpha”，其内容如下所示。
第一列（sample），样本名称；第二列（group），各样本所属的分组；第三列（variable），两种alpha多样性指数observed_species和shannon；第四列（value），各alpha多样性指数的数值。
（两组）数据提取
从中提取所期望比较的分组。例如在本示例中，我们期望查看group1和group2的observed_species指数是否存在显著差异。
#选择要比较的分组（此处查看group1与group2在observed_species指数上是否存在显著差异）richness_12<-subset(alpha,variable=='observed_species'&group%in%c('1','2'))richness_12$group<-factor(richness_12$group)提取分组数据后赋值为新数据框“richness_12”，并将分组列转化为因子类型。此时“richness_12”的内容如下所示。
只包含了group1和group2两个分组的observed_species指数数据，且分组列已转化为因子类型。
正态性假设检验
t检验的一个重要前提就是数据必须符合正态分布模型。因此在执行t检验之前必须验证数据分布的正态性。若数据不符合正态分布，则t检验将无法适用于该数据（此时可以考虑非参数的检验方法，例如wilcox秩和检验等；或者考虑对数据进行合理的转化，常见的如对数转化等，查看转化后的数据是否满足了t检验的前提假设）。验证数据是否符合正态分布的方法很多，以下展示两种常见方法。
正态QQ图
例如，可使用QQ图来观测数据分布。以下使用QQ图查看上述数据框“richness_12”中的数值分布，可使用car包中的qqPlot()来实现。
library(car)#QQ-plotqqPlot(lm(value~group,data=richness_12),simulate=TRUE,main='QQPlot',labels=FALSE)qqPlot()提供了精确的正态假设检验方法，它画出了在n-p-1个自由度的t分布下的学生化残差（studentized residual，也称学生化删除残差或折叠化残差）图形，其中n是样本大小，p是回归参数的数目（包括截距项）。
以上示例的qqPlot()作图结果如下。其中横坐标是标准的正态分布值，纵坐标是我们数据的值。如果两者基本相等，或者说所有的点都离直线很近，落在置信区间内（图中虚线部分，默认展示95%置信区间），即表明正态性假设符合得很好。由图可知，我们的数据符合正态分布模型，适合使用t检验进行分析。
QQ图这种类型的图示法方法简单，从图中可以直接判断，比较直观。但这种方法效率不是很高，不适用于很多样本的情形（当样本量很多时，我们不太可能通过QQ图逐一判断），它所提供的信息通常只作为正态性检验的重要补充。
Shapiro-Wilk检验
我们还可使用Shapiro-Wilk检验验证数据分布的正态性，它类似于线性回归的方法，是检验其于回归曲线的残差。R中提供了可用于执行Shapiro-Wilk检验的函数shapiro.test()，通过所得p值即可轻松地判断数据分布的正态性。该检验原假设（或称零假设）为数据集符合正态分布，通常情况下若通过Shapiro-Wilk检验得到p值小于0.05，则拒绝原假设，数据分布不符合正态性；反之接受原假设。
以下使用shapiro.test()继续查看上述数据框“richness_12”中的数值分布是否满足正态性。
#Shapiro-Wilk检验shapiro<-tapply(richness_12$value,richness_12$group,shapiro.test)shapiroshapiro$'1'$p.valueshapiro$'2'$p.value对于示例中的两个分组，group1（即“1”）和group2（即“2”）中的observed_species数值，分别给出了其检验p值。结果如下。可知两组的检验p值均大于0.05，即表明原假设成立（即数据分布符合正态性），我们的数据符合正态分布模型，适合使用t检验进行分析。
注：若两组数据中，某一组数据的Shapiro-Wilk检验p值小于了0.05（即便另一分组的Shapiro-Wilk检验p值均大于0.05），则对于整体数据而言，将不再适用t检验。
这种直接输出p值判断的方式，也比较直观。特别是当样本量较多时，不太适合通过上述提到的QQ图逐一判断时，可使用shapiro.test()去完成（例如通过循环的方式逐一提取两两数据并验证分布）。
以下展示一个批处理示例，使用上文的数据框“alpha”（包含了所有样本及所有alpha多样性指数），以循环的方式逐一挑选其中的两组数据后，并执行Shapiro-Wilk检验，最终查看其中各两两分组间的数据分布是否满足正态性。
#Shapiro-Wilk检验批量验证alpha$group<-factor(alpha$group)alpha_all<-levels(alpha$variable)group_all<-levels(alpha$group)result<-NULLfor(ainalpha_all){for(iin1:(length(group_all)-1)){for(jin(i+1):length(group_all)){dat<-subset(alpha,variable==a&group%in%c(group_all[i],group_all[j]))shapiro<-tapply(dat$value,dat$group,shapiro.test)p1<-shapiro[[group_all[i]]]$p.valuep2<-shapiro[[group_all[j]]]$p.valueifelse(p1<0.05|p2<0.05,pass<-'no',pass<-'yes')result<-rbind(result,c(a,group_all[i],p1,group_all[j],p2,pass))}}}result<-data.frame(result)names(result)<-c('alpha','group1','p.value1','group2','p.value2','pass')write.table(result,'shapiro_test.txt',sep='\t',row.names=FALSE,quote=FALSE)查看输出结果“shapiro_test.txt”，本示例如下所示。在各alpha多样性指数下，执行两两分组数据间的正态分布检验，并各自给出p值（分别为p.value1和p.value2）。若其中之一小于0.05，则标记为“no”，即数据分布未满足正态性；反之标记为正态性“yes”。后续即可根据数据分布，判断是否使用t检验（或者考虑使用非参数的检验方法，常见的如wilcox秩和检验等；或者考虑对数据进行合理的转化，常见的如对数转化等，查看转化后的数据是否满足了正态分布）。
t检验
当两组间的数据分布通过了正态假设检验后，即可执行t检验。可分为独立样本的t检验与非独立样本的t检验。
独立样本的t检验
假设样本间是相互独立的。例如，本示例数据来源于某16S测序数据集，假设不同group（1、2、3）的样本为土壤细菌群落样本，group1来自东北土壤，group2来自西北土壤，group1来自西南土壤，那么不同group的土壤群落之间可以认为是相互没有关联的，即各组相互独立。那么，可使用独立样本的t检验方法去处理。
在R中，t检验的调用方式如下示例。使用上文中的数据框“richness_12”，已知其中的数据满足了正态性分布，接下来想要获知group1和group2的observed_species指数间是否存在显著差异。
#独立样本的t检验t_test<-t.test(value~group,richness_12)t_testt_test$p.valuet.test()中默认两组间相互独立（默认参数paired = FALSE），执行独立样本的t检验。同时需注意，与其它多数统计软件不同的是，在R中的t检验默认假定方差不相等（默认参数var.equal = FALSE），并使用Welsh的修正自由度；可以通过参数“var.equal = TRUE”假定方差相等，并使用合并方差估计。此外，t.test()默认的备择假设是双侧的（默认参数alternative = 'two.sided'），即执行双侧检验；可分别使用“alternative = 'less'”或“alternative = 'greater'”执行单侧t检验。
上述示例的屏幕输出如下，根据p值可判断显著性。由于p值远小于0.05，即拒绝了原假设（原假设两组间没有差异），group1和group2的observed_species指数间存在显著不同。
非独立样本的t检验
假设样本间并非相互独立的。例如，假设测序数据样本来源于同一片区域的土壤，在土壤中添加了某化学物质后，分时期取样观察土壤细菌群落的时间变化趋势（以研究该化学物质对土壤细菌群落的影响），group1、2、3分别为第1、2、3个采样时期。此时，不同取样时期的土壤细菌群落之前并非毫无关联的，理论上来讲，组间数据并不独立。
这时可以考虑使用非独立样本的t检验方法去分析组间差异。当然，若此时仍然使用独立样本的t检验方法执行差异分析也是可行的，这种分析方法本身并没问题（仅仅是在统计算法上存在一些不同，相较之下可能非独立样本的t检验方法会更合适一些）。
继续对上文中的数据框“richness_12”中的数据执行非独立样本的t检验分析，如下示例。
#非独立样本的t检验t_test<-t.test(value~group,richness_12,paired=TRUE)t_testt_test$p.value此时在t.test()中设定参数“paired = TRUE”即可执行非独立样本的t检验。屏幕输出如下，根据p值（0.0002，远低于0.05）可知group1和group2的observed_species指数间存在显著不同。
可视化展示
若考虑作图将两组差异进行可视化展示。
例如，一个简单的箱线图示例，使用上文的数据框“richness_12”。
#boxplot()箱线图boxplot(value~group,data=richness_12,col=c('blue','orange'),ylab='Observed_species',xlab='Group',main='t-test:p-value<0.001')作图示例如下。
再例如，可考虑根据各样本的数值，分别计算各组中的均值以及标准差，展示为均值±标准差的柱形图样式。如下一个简单的示例，同样使用上文的数据框“richness_12”。
#ggplot2柱形图library(doBy)#使用其中的summaryBy()以方便按分组计算library(ggplot2)#ggplot2作图dat<-summaryBy(value~group,richness_12,FUN=c(mean,sd))ggplot(dat,aes(group,value.mean,fill=group))+geom_col(width=0.4,show.legend=FALSE)+geom_errorbar(aes(ymin=value.mean-value.sd,ymax=value.mean+value.sd),width=0.15,size=0.5)+theme(panel.grid=element_blank(),panel.background=element_rect(color='black',fill='transparent'),plot.title=element_text(hjust=0.5))+labs(x='Group',y='Observed_species',title='t-test:p-value<0.001')作图示例如下。
参考资料
