0

我需要帮助在 R 中编写一个循环，它可以让我对不同的参数进行一系列功率计算。Alpha 将始终 = 0.05，N 是将采用各种值的样本大小，因此：

N= 1000, 2000, 3000, 4000, 5000, 6000, 7000, 8000, 9000, 1000.

对于每个分析的 N，H2 也将采用不同的值，它们是：0.005、0.01、0. 02

例如，对于 N=1000，我想在 alpha=0.05 时计算 H2=0.005、H2=0.01 和 H2=0.02 的功率。然后我想对上面提到的所有NI值做这个分析。

这是我在下面仅运行一次的代码： 

N = 1000
alpha = 0.05
H2 = 0.005

threshold = qchisq(alpha, df = 1, lower.tail = FALSE)
power = pchisq(threshold, df = 1, lower.tail = FALSE, ncp = N * H2)

有人可以帮我把它变成一个循环，这样我就可以在结构化表格中一次得到这些结果吗？谢谢。

4

2 回答 2

1 
N <- c(1000, 2000, 3000, 4000, 5000, 6000, 7000, 8000, 9000, 10000)
H2 <- c(0.005, 0.01, 0.02)
alpha <- 0.05

threshold <- qchisq(alpha, df = 1, lower.tail = FALSE)

paraComb <- expand.grid(N, H2)
ncp <- with(paraComb, Var1 * Var2)

setNames(cbind(paraComb,
               sapply(ncp, 
                      function(ncp) pchisq(threshold, 
                                           df = 1, lower.tail = FALSE, 
                                           ncp = ncp))
              ),
         c("N", "H2", "power"))

#        N    H2     power
# 1   1000 0.005 0.6087795
# 2   2000 0.005 0.8853791
# 3   3000 0.005 0.9721272
# 4   4000 0.005 0.9940005
# 5   5000 0.005 0.9988173
# <snip>
于 2013-08-28T20:06:44.577 回答
0 
N <- c(1000, 2000, 3000, 4000, 5000, 6000, 7000, 8000, 9000, 10000)
N<-as.list(rep(N,each=3)) # 3 is the length of original H2
H2 <- c(0.005, 0.01, 0.02)
H2<-as.list(rep(H2,10)) # 10 is the length of original N

myresult<-Map(function(x,y) cbind(x,y,power=pchisq(qchisq(0.05, df = 1, lower.tail = FALSE), df = 1, lower.tail = FALSE, ncp = x * y)),N,H2)
myout<-do.call(rbind,myresult)
colnames(myout)[1:2]<-c("N","H2")
> head(myout,10) 
         N    H2  power
 [1,] 1000 0.005 0.6088
 [2,] 1000 0.010 0.8854
 [3,] 1000 0.020 0.9940
 [4,] 2000 0.005 0.8854
 [5,] 2000 0.010 0.9940
 [6,] 2000 0.020 1.0000
 [7,] 3000 0.005 0.9721
 [8,] 3000 0.010 0.9998
 [9,] 3000 0.020 1.0000
[10,] 4000 0.005 0.9940
于 2013-08-28T20:43:53.627 回答