我正在尝试对拟合优度执行 Pearson 卡方检验。以下是拟合泊松分布的示例:
data <- rpois(200,50)
estimate <- mean(data)
freq.os<-table(data)
yfit <- dpois(as.integer(names(freq.os)), estimate)
chisq.test(x = freq.os, p = yfit)
# Error in chisq.test(x = freq.os, p = yfit) : probabilities must sum to 1.
当我评估 sum(yfit)
时,我得到 0.999839。
那么如何产生一组加起来为1的概率值呢?谢谢!
编辑
实际上我找到了解决办法:
chisq.test(freq.os, yfit)
但我很困惑 chisq.test()
是如何工作的,因为它告诉我 df = 429
。我想 df = n - k - 1
,在这种情况下应该是 35,其中 k = 1
代表 lambda,n = number
卡方和中的项。我哪里做错了?
最佳答案
上面的评论建议您在将 yfit
传递给 chisq.test
之前手动重新缩放它。实际上,您可以让 chisq.test()
为您完成此操作:
chisq.test(x = freq.os, p = yfit, rescale.p = TRUE)
关于您的编辑
chisq.test(freq.os, yfit)
不正确,因为它正在进行独立性测试。
chisq.test()
可以进行两种统计检验:
- 拟合优度检验,使用参数
x
和p
; - 独立性测试,使用参数
x
和y
。
请仔细阅读?chisq.test
。为了进行拟合优度测试,您必须像最初那样在函数中使用 p
参数。我上面的回答,使用 rescale.p = TRUE
将帮助您解决您看到的错误。
如何进行 Pearson 卡方检验
你说你不知道测试是怎么做的,那就看这部分。
您应该使用 set.seed()
,这样当其他人运行您的代码时,他们会得到与您相同的随机数。这是一个可重现的例子:
N <- 200 ## number of samples
set.seed(0) ## fix random seed for reproducibility
x <- rpois(N,50) ## generate sample
lambda <- mean(x) ## estimate mean
现在我们的目标是使用 Pearson 卡方检验来检验拟合优度。我们有
Null Hypothesis: x ~ Poisson(lambda)
O <- table(x) ## contingency table / observed frequency
n <- length(O) ## number of categories
# 31
x0 <- as.integer(names(O)) ## unique sample values
p <- dpois(x0, lambda); p <- p / sum(p) ## theoretical probability under Null Hypothesis
我们先来做卡方检验。
E <- p * N ## expected frequency under Null Hypothesis
R <- (O - E) / sqrt(E) ## standardised residuals
X <- sum(R * R) ## Chi-square statistic
# 36.13962
dof <- n - 1 ## degree of freedom
# 30
p.value <- 1 - pchisq(X, df = dof) ## p-value
# 0.2035416
p 值大于 0.05,因此我们不能拒绝原假设。概要情节:
z <- curve(dchisq(x, df = dof), from = 0, to = 100)
abline(v = X, lty = 3)
polygon(c(X, X, z$x[z$x > X]), c(0, dchisq(X, df = dof), z$y[z$x > X]), col = 2)
然后我们使用chisq.test()
chisq.test(O, p = p)
#Chi-squared test for given probabilities
#data: O
#X-squared = 36.14, df = 30, p-value = 0.2035
可以看到,结果是一样的。
关于r - 如何使用 chisq.test() 正确进行卡方检验?,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/38276295/