本文旨在解决上一篇所提问题，具体问题如下：应用中对于多组均数的比较，当方差不齐时用非参数检验是否恰当，请模拟回答。
一、检验方法介绍
（一）One-Way ANOVA（单因素方差分析）
- 定义与用途：
单因素方差分析（One-Way ANOVA）是一种用于比较三个或三个以上组别间均值差异的统计方法。它适用于一个自变量（也可称因子）下，有一个以上水平（或组别）的情况，以探究因素对因变量的影响是否显著。
- 基本原理：
单因素方差分析基于对变量的方差进行分解，将观测到的变异分解为组内变异（即组内各样本的变异）和组间变异（不同组别均值之间的差异）。通过计算适当的统计量（如F值）来比较这两种变异的大小，从而判断不同组别的均值是否存在显著性差异。
- 假设条件：
- 各组数据服从正态分布。
- 各组方差相等（方差齐性）。
- 观察值之间相互独立。
- 应用场景：
例如，为研究不同锻炼方法对青年亚健康人群血脂的影响，将人群随机分为多个组别，分别进行不同的锻炼方法，最后通过单因素方差分析比较各组血脂水平的差异。
（二）Welch近似检验
- 定义与用途：
韦尔奇检验（Welch's test），又称非等方差t检验，是一种用于比较两组独立样本均值差异的统计方法。它在两组样本方差明显不等时，比经典的t检验更稳健，能够提供更可靠的检验结果。
- 基本原理：
韦尔奇检验通过调整t检验的方差估计来适应样本方差不等的情况，从而得到更准确的统计推断。它使用Satterthwaite公式来计算有效的自由度，进而计算t值和p值。
- 假设条件：
韦尔奇检验不要求样本方差相等，但通常假设样本数据来自正态分布或近似正态分布。
- 应用场景：
当两组样本的方差明显不等时，使用韦尔奇检验可以更有效地判断两组样本均值是否存在显著差异。
（三）Kruskal-Wallis H检验
- 定义与用途：
克鲁斯卡尔-沃利斯检验（Kruskal-Wallis test），亦称“K-W检验”或“H检验”，是一种用于检验两个以上样本是否来自同一个概率分布的非参数方法。它不需要假设样本来自正态分布，适用于样本数据不满足ANOVA假设条件的情况。
- 基本原理：
克鲁斯卡尔-沃利斯检验是一种秩检验，通过将所有样本值合并后排序，并用秩代替原始值，然后计算各样本秩的和及相应的统计量（如H值），以判断各样本是否来自同一个概率分布。
- 假设条件：
- 被检验的几个样本必须是独立的或不相关的。
- 原假设是各样本服从的概率分布具有相同的中位数。
- 应用场景：
当样本数据不满足正态分布或方差齐性等ANOVA假设条件时，可以使用克鲁斯卡尔-沃利斯检验来判断多个样本是否来自同一个概率分布。例如，在研究不同药物对某种疾病治疗效果的影响时，如果样本数据不满足ANOVA的假设条件，可以使用K-W检验来进行分析。
二、R语言的实现
为了方便获得三组在3种标准差设置中三种多重比较方法在不同样本量下的一类错误以及检验效能， 我们构建了一个综合函数Hypothesis_testing_simulation()以至于能简便得到模拟结果中的一类错误以及检验效能等信息，为了直观比较不同比较方法下一类错误或检验效能，该函数加入了是否需要绘制一类错误或检验效能随样本量的变化图选项，以下为该函数参数相关信息：
- n: 模拟样本量设置
- iter: 模拟次数
- seed: 种子数
- mu: 三组的均值
- sd: 三组的方差
- plot.errpr: 是否绘制三种多重比较方法一类错误随样本量的变化图(默认FALSE)
- plot.power: 是否绘制三种多重比较方法检验效能随样本量的变化图(默认FALSE)
Hypothesis_testing_simulation <-function(n,iter,seed,mu,sd,
plot.error=FALSE,plot.power=FALSE){
value1 <- numeric()
value2 <- numeric()
value3 <- numeric()
for(i in 1:length(n)){
#模拟次数
p1 <- numeric()
p2 <- numeric()
p3 <- numeric()
#设置种子数
set.seed(seed)
for(j in 1:iter){
#创建模拟值
Group1 <- rnorm(n[i],mu[1],sd[1])
Group2 <- rnorm(n[i],mu[2],sd[2])
Group3 <- rnorm(n[i],mu[3],sd[3])
#创建为多重比较函数适应格式
data <- data.frame(value = c(Group1,Group2,Group3),
group = rep(1:3,each=n[i]))
#计算三种多重比较方法的p值
#One-way ANONA
p1[j] <- oneway.test(value ~ group,data,var.equal = TRUE)$p.value
#Welch近似
p2[j] <- oneway.test(value ~ group,data,var.equal = FALSE)$p.value
#Kruskal-Wallis H检验
p3[j] <- kruskal.test(value ~ group,data)$p.value
}
value1[i] <- mean(p1 <= 0.05)
value2[i] <- mean(p2 <= 0.05)
value3[i] <- mean(p3 <= 0.05)
}
#绘制各多重比较方法一类错误随样本量变化曲线图
if(plot.error == TRUE){
par(mar=c(5,5,2,2),xaxs="i",yaxs="i")
plot(n,value1,lwd=4,type="l",bty="l",lty=1,xaxt="n",yaxt="n",xlab="",
ylab="",col=1,ylim=c(min(c(value1,value2,value3)),max(c(value1,value2,
value3))))
axis(1,las=1,cex.axis=1.3,lwd=1.6)
axis(2,las=1,cex.axis=1.3,lwd=1.6)
lines(n,value2,col=2,lwd=4,lty=2)
lines(n,value3,col=3,lwd=4,lty=3)
title(xlab="样本量",font.lab=7,cex.lab=1.5,line=2.6)
title(main="不同多重比较方法一类错误随样本量变化")
legend("topright",c("One-wayANONA","Welch近似","Kruskal-Wallis H检验"),
lwd=4,col=c(1:3),lty=c(1:3),bty="n",xpd=TRUE,cex=1.1)
}
#绘制各多重比较方法检验效能随样本量变化曲线图
if(plot.power == TRUE){
par(mar=c(5,5,2,2),xaxs="i",yaxs="i")
plot(n,value1,lwd=4,type="l",bty="l",lty=1,xaxt="n",yaxt="n",xlab="",
ylab="",col=1,ylim=c(min(c(value1,value2,value3)),max(c(value1,value2,
value3))))
axis(1,las=1,cex.axis=1.3,lwd=1.6)
axis(2,las=1,cex.axis=1.3,lwd=1.6)
lines(n,value2,col=2,lwd=4,lty=2)
lines(n,value3,col=3,lwd=4,lty=3)
title(xlab="样本量",font.lab=7,cex.lab=1.5,line=2.6)
title(main="不同多重比较方法检验效能随样本量变化")
legend("topright",c("One-wayANONA","Welch近似","Kruskal-Wallis H检验"),
lwd=4,col=c(1:3),lty=c(1:3),bty="n",xpd=TRUE,cex=1.1)
}
return(cbind(n,value1,value2,value3))
}利用搭建好的函数进行我们所要做的模拟如下：
- 不同标准差下的一类错误以及变化图：
#标准差为1,1,1
Hypothesis_testing_simulation(n = c(20,50,100,500),iter = 10000,seed = 1234,
mu =c(0,0,0),sd = c(1,1,1),plot.error = TRUE)
#标准差为1,2,5
Hypothesis_testing_simulation(n = c(20,50,100,500),iter = 10000,seed = 1234,
mu =c(0,0,0),sd = c(1,2,5),plot.error = TRUE)
#标准差为1,5,10
Hypothesis_testing_simulation(n = c(20,50,100,500),iter = 10000,seed = 1234,
mu =c(0,0,0),sd = c(1,5,10),plot.error = TRUE)- 不同标准差下的检验效能以及变化图：
#标准差为1，1，1
Hypothesis_testing_simulation(n = c(20,50,100,500),iter = 10000,seed = 1234,
mu =c(0,0.5,0.5),sd = c(1,1,1),
plot.power = TRUE)
#标准差为1,2,5
Hypothesis_testing_simulation(n = c(20,50,100,500),iter = 10000,seed = 1234,
mu =c(0,0.5,0.5),sd = c(1,2,5),
plot.power = TRUE)
#标准差为1,5,10
Hypothesis_testing_simulation(n = c(20,50,100,500),iter = 10000,seed = 1234,
mu =c(0,0.5,0.5),sd = c(1,5,10),
plot.power = TRUE)三、结果显示
- 标准差为1，1，1
## n value1 value2 value3
## [1,] 20 0.0482 0.0472 0.0468
## [2,] 50 0.0505 0.0493 0.0487
## [3,] 100 0.0514 0.0523 0.0505
## [4,] 500 0.0499 0.0500 0.0497## n value1 value2 value3
## [1,] 20 0.3394 0.3319 0.3214
## [2,] 50 0.7234 0.7192 0.7032
## [3,] 100 0.9612 0.9606 0.9515
## [4,] 500 1.0000 1.0000 1.0000
- 标准差为1，2，5
## n value1 value2 value3
## [1,] 20 0.0725 0.0454 0.0673
## [2,] 50 0.0756 0.0477 0.0694
## [3,] 100 0.0715 0.0511 0.0705
## [4,] 500 0.0728 0.0504 0.0720## n value1 value2 value3
## [1,] 20 0.0925 0.1345 0.1119
## [2,] 50 0.1155 0.2914 0.2042
## [3,] 100 0.1636 0.5502 0.3729
## [4,] 500 0.8368 0.9985 0.9874
- 标准差为1,5,10
## n value1 value2 value3
## [1,] 20 0.0685 0.0459 0.0726
## [2,] 50 0.0747 0.0485 0.0755
## [3,] 100 0.0705 0.0518 0.0787
## [4,] 500 0.0698 0.0512 0.0802## n value1 value2 value3
## [1,] 20 0.0753 0.0646 0.0852
## [2,] 50 0.0819 0.0934 0.1034
## [3,] 100 0.0863 0.1443 0.1271
## [4,] 500 0.1878 0.5828 0.4086
个人见解，若有不对，请留言，鼠鼠会及时修正......
扫一扫
2016
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
