【R分享|实战】参数分析~T检验与方差分析
T-test
#设置种子数，保重每次结果一致
set.seed (1234)
x <- rnorm(n=50,mean=122.8,sd=8.5) # 生成50个符合正态，均值为122.8，标准差为8.5的样本
hist(x,col="light blue") # 画直方图的基础函数单一样本，已知总体均值，所以采用单一样本t检验即可50例偷吃零食员工与认真工作的员工的体重是否有差异。
t.test(x,mu=140) # 进行t检验结果如下：
以0.05水平为标准，结果为p<0.05，具体统计学意义，表明偷吃与认真的员工体重是存在显著差异。
Example 2 -- 两独立样本t检验
前提条件：1）两个样本应该是相互独立的；2）样本来自的两个总体应该服从正态分布。
根据样本数据对两个样本来自的两个独立总体的均值是否有显著差异进行判断。
使用甲乙两台机床来加工同种零件，两种机床加工零件的尺寸服从正态分布，且方差相同，为检验甲乙两台机床加工的零件尺寸均值是否一致（p=0.05），从两种机床加工的零件中分别抽取若干零件测得的其尺寸如下：
#甲机床加工数据
x1<-c(20.5,19.8,19.7,20.4,20.1,20.0,19.0,19.9)
#乙机床加工数据
x2<-c(20.7,19.8,19.5,20.8,20.4,19.6,20.2)
#进行两样本t检验，指定方差相等
t.test(x1,x2,var.equal=T)结果：
从结果来看，p=0.4081>0.05，则不能拒绝原假设，即认为两台机床加工零件的尺寸均值没有显著性区别。同时结果也给出了两台机床加工零件尺寸样本数据的均值分别为：19.925和20.14286。
Example 3 -- 配对t检验
前提条件：1）两样本必须是配对的.这里有两层意思,即一是两样本的观察值数目相同；两样本的观察值的顺序不随意更改；2）样本来自的两个总体必须服从正态分布。
配对设计主要有4种情况：同一受试对象处理前后的数据，同一受试对象两个部位的数据，同一样品用两种方法检验的结果，配对的两个受试对象分别接受两种处理后的数据。
判断简便法和常规法测定尿铅含量的差别有无统计意义，对12份人尿同时用两种方法进行测定，所得结果如下表所示，请分析两种测定方法的测量结果是否不同？
#读取两组样本的数据
x <- c(2.41,2.90,2.75,2.23,3.67,4.49,5.16,5.45,2.06,1.64,1.06,0.77)
y <- c(2.80,3.04,1.88,3.43,3.81,4.00,4.44,5.41,1.24,1.83,1.45,0.92)
#配对样本t检验
t.test(x,y,paired=T) # paired=T 表示执行配对t检验配对t检验的结果为：p=0.874>0.05，不能拒绝原假设H0，显然两种方法测定尿铅含量的结果没有显著差异。
Anova analysis
在实验设计的时候，通常以单因素或者双因素为主。这里主要分享单因素和双因素方差分析（偶尔用到大于两个因素的实验，方法相同）。主要用到aov() 函数。
Example 1 -- 单因素方差分析
前提条件：1）样本独立性；2）样本总体必须服从正态分布；3）方差齐性。
我们选择R自带的数据集PlantGrowth，研究两个处理和一个对照组对植物产量的影响，每组10例共3记录，主要考察处理对提高植物产量有无影响。
数据构成：因变量weight，因子变量group，三个水平依次为ctrl、trt1、trt2。
常规的数据分析顺序，先进行三大检验（复习上期内容）。通常我用以下三个函数进行数据的检验，chisq.test-独立性检验、shapiro.test-正态性、bartlett.test-方差齐性。
#1 单因素方差分析
head(PlantGrowth)
# 独立性检验
mytable1 <- table(PlantGrowth$group,PlantGrowth$weight)
chisq.test(mytable1)
# 正态性检验
mydata <- PlantGrowth[,1]
shapiro.test(mydata[1:10]) # 对应group 每组水平下的检验
shapiro.test(mydata[11:20])
shapiro.test(mydata[21:30])
# 方差齐性
bartlett.test(weight~group,data = PlantGrowth)从三大检验的结果来看，p值均>0.05，也就是都接受了原假设，满足条件。
单因素方差分析
ANOVA <- aov(weight~group,data = PlantGrowth)# weight 因变量在~的左边，group 自变量/处理在~右边
summary(ANOVA)方差分析表的p=0.0159<0.05，拒绝原假设（H0：3组样本数据均值相等，无差异），说明3组样本均值有差异，不同处理方式对植物产量有显著影响。
特别注意：对于方差不齐的情况，我们可以尝试用一下oneway.test()函数。
oneway.test(weight~group,data = PlantGrowth,var.equal = F)# “var.equal = F” 表示两样本方差不同（缺省状态，默认状态），“var.equal = T”表示方差相同#由于前面检验了方差齐性，是相当的。因此，这里结果与上面aov的结果相同。
多重比较
这里介绍两种LSD.test() 和 Duncan() 函数
#1 LSD.test
library(agricolae) # 加载agricolae包
mutl2<-LSD.test(ANOVA, 'group',p.adj = 'bonferroni') #“p.adj = bonferroni”表示对P值进行bonferroni矫正。矫正方法还有：none, holm，hommel，hochberg，BH”, BY，fdr等#
print(mutl2$groups)
plot(mutl2)
#2 Duncan新复极差法
library(agricolae)
mult5<- duncan.test(ANOVA,"group")
print(mutl2$groups)
plot(mutl2)
这两种能够直接展示显著性差异的结果
Example 2 -- 双因素方差分析
前提条件：与单因素方差一致
# 双因素方差分析
head(ToothGrowth)
summary(ToothGrowth) # 查看数据集的总体情况
str(ToothGrowth) # 查看数据结构# 将dose变量转换成factor
ToothGrowth$dose <- factor(ToothGrowth$dose)
str(ToothGrowth)
# 利用table函数查看每个分组的样本量情况
table(ToothGrowth$supp,ToothGrowth$dose)library(ggplot2)
ggplot(ToothGrowth,aes(dose,len,color=supp))+
# data=Too... x=dose,y=len,对不同方式进行绘制颜色
geom_boxplot() # geom_boxplot 绘制默认的箱线图
进一步，通过双因素方差分析来判断不同因子对豚鼠牙齿生长的影响。
# 未包含dose*supp or dose:supp的交互项
anova.1 <- aov(len ~ dose + supp, data = ToothGrowth)
summary(anova.1)结果：
如果认为dose与supp之间可能存在交互作用；即，想要研究给药剂量（dose）与牙齿生长（len）之间的关系会不会随着给药方式的不同（supp）而改变，则添加交互项。
# 含交互项
anova.2 <- aov(len ~ dose + supp + dose:supp, data = ToothGrowth)
summary(anova.2)
anova.2 <- aov(len ~ dose + supp + dose*supp, data = ToothGrowth)
summary(anova.2)结果：
###########################
# 参数分析之t检验与方差分析
###########################
# t检验
#1 单样本t检验
#设置种子数，保证每次结果一致
set.seed (1234)
x <- rnorm(n=50,mean=122.8,sd=8.5) # 生成50个符合正态，均值为122.8，标准差为8.5的样本
hist(x,col="light blue") # 画直方图的基础函数
t.test(x,mu=140) # 进行t检验
#2 两独立样本t检验
#甲机床加工数据
x1<-c(20.5,19.8,19.7,20.4,20.1,20.0,19.0,19.9)
#乙机床加工数据
x2<-c(20.7,19.8,19.5,20.8,20.4,19.6,20.2)
#进行两样本t检验，指定方差相等
t.test(x1,x2,var.equal=T)
#3 配对t检验
#读取两组样本的数据
x <- c(2.41,2.90,2.75,2.23,3.67,4.49,5.16,5.45,2.06,1.64,1.06,0.77)
y <- c(2.80,3.04,1.88,3.43,3.81,4.00,4.44,5.41,1.24,1.83,1.45,0.92)
#配对样本t检验
t.test(x,y,paired=T) # paired=T 表示执行配对t检验
#方差分析
#1 单因素方差分析
head(PlantGrowth)
# 独立性检验
mytable1 <- table(PlantGrowth$group,PlantGrowth$weight)
chisq.test(mytable1)
# 正态性检验
mydata <- PlantGrowth[,1]
shapiro.test(mydata[1:10])
shapiro.test(mydata[11:20])
shapiro.test(mydata[21:30])
# 方差齐性
bartlett.test(weight~group,data = PlantGrowth)
# 单方差分析
ANOVA <- aov(weight~group,data = PlantGrowth)# weight 因变量在~的左边，group 自变量/处理在~右边
summary(ANOVA)
oneway.test(weight~group,data = PlantGrowth,var.equal = F)# “var.equal = F” 表示两样本方差不同（缺省状态，默认状态），“var.equal = T”表示方差相同#
# 多重比较
#1 LSD.test
library(agricolae) # 加载agricolae包
mutl2<-LSD.test(ANOVA, 'group',p.adj = 'bonferroni') #“p.adj = bonferroni”表示对P值进行bonferroni矫正。矫正方法还有：none, holm，hommel，hochberg，BH”, BY，fdr等#
print(mutl2$groups)
plot(mutl2)
#2 Duncan新复极差法
library(agricolae)
mult5<- duncan.test(ANOVA,"group")
print(mutl2$groups)
plot(mutl2)
# 双因素方差分析
head(ToothGrowth)
summary(ToothGrowth) # 查看数据集的总体情况
str(ToothGrowth) # 查看数据结构
# 将dose变量转换成factor
ToothGrowth$dose <- factor(ToothGrowth$dose)
str(ToothGrowth)
# 利用table函数查看每个分组的样本量情况
table(ToothGrowth$supp,ToothGrowth$dose)
library(ggplot2)
ggplot(ToothGrowth,aes(dose,len,color=supp))+
# data=Too... x=dose,y=len,对不同方式进行绘制颜色
geom_boxplot() # geom_boxplot 绘制默认的箱线图
# 未包含dose*supp or dose:supp的交互项
anova.1 <- aov(len ~ dose + supp, data = ToothGrowth)
summary(anova.1)
# 含交互项
anova.2 <- aov(len ~ dose + supp + dose:supp, data = ToothGrowth)
summary(anova.2)
anova.2 <- aov(len ~ dose + supp + dose*supp, data = ToothGrowth)
summary(anova.2)