这个问题是我读paper时引发的。首先介绍一下我提出问题的由来。
这篇文章是讲模型选择的。文章后面还有无数个人的评论。J. Sunil Rao and Robert Tibshirani的评论通过重复原文作者Shao的simulation study指出他们的方法Adaptive CCP比Shao的CVd要好。但是看完两张截然不同的simulation results tables之后,我对J. Sunil Rao and Robert Tibshirani的结果表示怀疑。(我自己的结果当然是和Shao的consistent)
他们集中于比较各种模型选择方法选中正确模型的概率p。下面我们先只看scenario1。R & T 用Shao的方法模拟得出Shao的CVd选中正确模型78次out of 100次。而Shao用自己的方法得出CVd选中正确模型934次out of 1000次。另外3个scenarios里,两组用CVd 选中正确模型的次数是
Shao: 947 965 948 (based on 1000 simulations)
R & T: 76 88 94 (based on 100 simulations)
于是我建立原假设:两组用的CVd是同一种方法。(也就是对于任意一个scenario i, 两组的pi应该相等)
备择假设:otherwise
用Fisher exact test + Bonferroni method可以得到如下结果
<br />
> total=c(100,1000)<br />
> shao=c(934,947,965,948)<br />
> RT=c(78,76,88,94)<br />
> a=c(.5,.1,.05)<br />
> pvalue=c()<br />
> for(i in seq(shao)){<br />
+ table <- matrix(c(RT[i],total[1]-RT[i],shao[i],total[2]-shao[i]), 2, 2)<br />
+ pvalue[i]=fisher.test(table)[[1]]<br />
+ }<br />
> pvalue<br />
[1] 2.962623e-06 7.926365e-09 5.802546e-04 6.426553e-01<br />
> pvalue*length(shao)<br />
[1] 1.185049e-05 3.170546e-08 2.321018e-03 2.570621e+00<br />
因此原假设在alpha level 小于等于3.170546e-08的情况下将被拒绝。
用中心极限定理作近似的检验如下:
<br />
> for(i in seq(shao)){<br />
+ p=(shao[i]+RT[i])/sum(total)<br />
+ x=(shao[i]/total[2]-RT[i]/total[1])/((1/total[1]+1/total[2])*p*(1-p))^.5<br />
+ pvalue[i]=2*pnorm(abs(x),lower.tail=F)<br />
+ }<br />
> pvalue<br />
[1] 6.220470e-08 2.787725e-12 6.141612e-05 7.328765e-01<br />
> pvalue*length(shao)<br />
[1] 2.488188e-07 1.115090e-11 2.456645e-04 2.931506e+00<br />
原假设在alpha level 小于等于1.115090e-11的情况下也将被拒绝。
因此,我怀疑R & T的模拟有问题。
参考文献:
[quote]http://www3.stat.sinica.edu.tw/statistica/password.asp?vol=7&num=2&art=1
[/quote]