作者toto01 (ur ...)
看板Statistics
標題[程式] 使用Monte Carlo估計power
時間Sat Jun 4 22:16:15 2011
小弟使用R去估計power,估計的方式是Monte Carlo,
mu <- 500 # The null hypothesis mean rqual 500
sigma <- 100 # Standard deviation
n <- 20 # The sample size of X is 20
nsim <- 10000 # The number of simulations
pv <- rep(NA, nsim) # empty vector to store p-values from t-test
for (i in 1:nsim){
set.seed(i)
x <- rnorm(n, mu, sigma)
pv[i] <- t.test(x, alternative = "two.sided",mu = 500)$p.value
}
mean(pv < 0.05) # The significat level is alpha = 0.05
hist(pv,prob=T)
最後算出power!
小弟想請教如何對這模擬10000次的t-test畫出empirical power curve
PS:模擬10000次的結果用長條圖呈現似乎不是常態分布
誠心求教
--
※ 發信站: 批踢踢實業坊(ptt.cc)
◆ From: 68.81.95.37
1F：→ bmka:如果你是simulate under the null,那麼這10000個p-value值 06/04 22:49
2F：→ bmka:畫histogram應該長得像uniform distribution,不會是normal 06/04 22:50
3F：→ bmka:(就算不是simulated under the null,也通常不會是normal) 06/04 22:51
4F：→ toto01:那如何對這模擬10000次的t-test畫出empirical power curve 06/04 22:56
5F：→ toto01:? 06/04 22:57
6F：→ toto01:沒錯!的確長得像uniform distribution 06/04 22:57
7F：→ bmka:你先搞清楚power 是什麼再說 06/05 10:36