r - 对 Gamma 分布使用 anova() 给出看似随机的 p 值

Question

我试图确定两个 Gamm 分布之间是否存在显着差异。一种分布具有 (shape,scale)=(shapeRef,scaleRef)，而另一种分布具有 (shape,scale)=(shapeTarget,scaleTarget)。我尝试使用以下代码进行方差分析

n=10000
x=rgamma(n, shape=shapeRef, scale=scaleRef)
y=rgamma(n, shape=shapeTarget, scale=scaleTarget)
glmm1 <- gam(y~x,family=Gamma(link=log))
anova(glmm1)

生成的 p 值不断变化，可以在 <0.1 到 >0.9 之间。

我的处理方式是否错误？

编辑:我使用以下代码代替

f <- gl(2, n)
x=rgamma(n, shape=shapeRef, scale=scaleRef)
y=rgamma(n, shape=shapeTarget, scale=scaleTarget)
xy <- c(x, y)
anova(glm(xy ~ f, family = Gamma(link = log)),test="F")

但是，每次运行它时，我都会得到不同的 p 值。

Answer 1

如果您每次都选择不同的实现，那么您每次运行时确实会得到不同的 p 值。就像您的数据值是随机变量一样，您希望每次运行实验时数据值都会发生变化，p 值也是如此。如果原假设为真(您最初尝试的情况就是如此)，则 p 值将均匀分布在 0 和 1 之间。

生成模拟数据的函数:

simfun <- function(n=100,shapeRef=2,shapeTarget=2,
                    scaleRef=1,scaleTarget=2) {
    f <- gl(2, n)
    x=rgamma(n, shape=shapeRef, scale=scaleRef)
    y=rgamma(n, shape=shapeTarget, scale=scaleTarget)
    xy <- c(x, y)
    data.frame(xy,f)
}

运行 anova() 并提取 p 值的函数:

sumfun <- function(d) {
    aa <- anova(glm(xy ~ f, family = Gamma(link = log),data=d),test="F")
    aa["f","Pr(>F)"]
}

尝试 500 次:

set.seed(101)
r <- replicate(500,sumfun(simfun()))

p 值始终非常小(尺度参数的差异很容易区分)，但它们确实有所不同:

par(las=1,bty="l") ## cosmetic
hist(log10(r),col="gray",breaks=50)