set.seed(139)
k <- 1000
m1 <- 10
m2 <- 10
m3 <- 10
power <- c(rep(0,5))
pvalue <- c(rep(0,k))
for(i in 1:k){
res1 <- -3 + rnorm(m1,0,1)
res2 <- -2 + rnorm(m2,0,1)
res3 <- -1.5 + rnorm(m3,0,1)
res <- c(res1,res2,res3)
fac <- c(rep("A",m1),rep("B",m2),rep("C",m3))
mod <- aov(res~fac)
mod.summary <- summary(mod)
pvalue[i] <- mod.summary[[1]][,5][[1]]
}
power[1] <- length(pvalue[which(pvalue < 0.05)])/k
set.seed(139)
pvalue <- c(rep(0,k))
for(i in 1:k){
res1 <- -3 + rnorm(m1,0,2)
res2 <- -2 + rnorm(m2,0,1)
res3 <- -1.5 + rnorm(m3,0,1)
res <- c(res1,res2,res3)
fac <- c(rep("A",m1),rep("B",m2),rep("C",m3))
mod <- aov(res~fac)
mod.summary <- summary(mod)
pvalue[i] <- mod.summary[[1]][,5][[1]]
}
power[2] <- length(pvalue[which(pvalue < 0.05)])/k
set.seed(139)
pvalue <- c(rep(0,k))
for(i in 1:k){
res1 <- -3 + rnorm(m1,0,1)
res2 <- -2 + rnorm(m2,0,1)
res3 <- -1.5 + rnorm(m3,0,2)
res <- c(res1,res2,res3)
fac <- c(rep("A",m1),rep("B",m2),rep("C",m3))
mod <- aov(res~fac)
mod.summary <- summary(mod)
pvalue[i] <- mod.summary[[1]][,5][[1]]
}
power[3] <- length(pvalue[which(pvalue < 0.05)])/k
set.seed(139)
pvalue <- c(rep(0,k))
for(i in 1:k){
res1 <- -3 + rnorm(m1,0,2)
res2 <- -2 + rnorm(m2,0,1)
res3 <- -1.5 + rnorm(m3,0,2)
res <- c(res1,res2,res3)
fac <- c(rep("A",m1),rep("B",m2),rep("C",m3))
mod <- aov(res~fac)
mod.summary <- summary(mod)
pvalue[i] <- mod.summary[[1]][,5][[1]]
}
power[4] <- length(pvalue[which(pvalue < 0.05)])/k
set.seed(139)
pvalue <- c(rep(0,k))
for(i in 1:k){
res1 <- -3 + rnorm(m1,0,2)
res2 <- -2 + rnorm(m2,0,2)
res3 <- -1.5 + rnorm(m3,0,2)
res <- c(res1,res2,res3)
fac <- c(rep("A",m1),rep("B",m2),rep("C",m3))
mod <- aov(res~fac)
mod.summary <- summary(mod)
pvalue[i] <- mod.summary[[1]][,5][[1]]
}
power[5] <- length(pvalue[which(pvalue < 0.05)])/k
data.frame(Power=power, SD.no.temp=c(1,1,1,1,2), SD.pos.treat=c(1,2,1,2,2), SD.neg.treat=c(1,1,2,2,2))