r - 在 R 中使用 Kolmogorov Smirnov 检验

我设计了3000个实验，所以在一个实验中有4组(治疗)，每组有50个人(受试者)。对于每个实验，我都做了一个标准的单向方差分析，并证明它们的 p.values 在原假设下是否具有 uni 概率函数，但是 ks.test 拒绝了这个假设，我不明白为什么？

subject<-50
treatment<-4
experiment<-list()
R<-3000
seed<-split(1:(R*subject),1:R)
for(i in 1:R){
  e<-c()
  for(j in 1:subject){
    set.seed(seed[[i]][j])
    e<-c(e,rmvnorm(mean=rep(0,treatment),sigma=diag(3,4),n=1,method="chol"))
   }
  experiment<-c(experiment,list(matrix(e,subject,treatment,byrow=T)))
 }

 p.values<-c()
for(e in experiment){
  d<-data.frame(response=c(e),treatment=factor(rep(1:treatment,each=subject)))
  p.values<-c(p.values,anova(lm(response~treatment,d))[1,"Pr(>F)"])
 }

 ks.test(p.values, punif,alternative = "two.sided")

最佳答案

我注释掉了您代码中更改随机种子的行，并获得了 0.34 的 P 值。那是一个未知的种子，所以为了可重复性，我做了 set.seed(1) 并再次运行它。这一次，我得到了 0.98 的 P 值。

至于为什么这会有所作为，我不是 PRNG 方面的专家，但任何体面的生成器都将确保连续抽奖在所有实际用途中在统计上都是独立的。对于更大的滞后，最好的将确保相同，例如，作为 R 的默认 PRNG 的 Mersenne Twister 保证它的滞后高达 623 (IIRC)。事实上，干预种子很可能会损害平局的统计特性。

您的代码也在以一种非常低效的方式做事。您正在为实验创建一个列表，并为每个实验添加一个项目。在每个实验中，您还创建一个矩阵，并为每个观察添加一行。然后你对 P 值做一些非常相似的事情。我看看能不能解决

这就是我替换您的代码的方式。严格来说，我可以通过避免公式、创建裸模型矩阵并直接调用 lm.fit 来使其更紧密。但这意味着必须手动编码 ANOVA 测试，而不是简单地调用 anova ，这比它的值(value)更麻烦。

set.seed(1) # or any other number you like

x <- factor(rep(seq_len(treatment), each=subject))
p.values <- sapply(seq_len(R), function(r) {
    y <- rnorm(subject * treatment, s=3)
    anova(lm(y ~ x))[1,"Pr(>F)"]
})
ks.test(p.values, punif,alternative = "two.sided")


        One-sample Kolmogorov-Smirnov test

data:  p.values
D = 0.0121, p-value = 0.772
alternative hypothesis: two-sided

关于r - 在 R 中使用 Kolmogorov Smirnov 检验，我们在Stack Overflow上找到一个类似的问题：https://stackoverflow.com/questions/17381214/