数据框中每行的卡方检验

Question

我有一个数据框，其中包含同一过程的两个观察者独立计数。

obs.1 <- c(2,10,53,13,12,15,5)
obs.2 <- c(3,12,45,2,7,17,5)
df <- data.frame(obs.1,obs.2)

我希望在每一行上使用卡方检验（R中的chisq.test“MASS”）来查看obs.1到obs.2之间是否存在显着差异。 我想将结果（x-squared，p-value）添加到df中。 我觉得apply函数是实现这一目标的正确方法，但我没有成功。

- doncarlos

---

你尝试过 cbind(df, t(apply(df, 1, function(x) {ch <- chisq.test(x); c(unname(ch$statistic), ch$p.value)}))) 吗？ - akrun

@CathG 我使用 chisq，因为它在其他类似的例子中也被使用。Kappa 只适用于分类数据吗？ - doncarlos

如果您对使用哪种测试（在一般统计问题中）存在疑问，@doncarlos，那么http://stats.stackexchange.com/可能是一个更好的发布问题的地方。 - akrun

@akrun，在更深入的思考后，我会改变我的第一个kappa想法，转而使用Wilcoxon（或t检验，取决于点数），主要是因为kappa确实更适用于分类数据，观察者之间的微小差异可能会导致糟糕的kappa系数，而这可能并不是真正的显著差异。但我想，这确实取决于“数据的性质”。 - Cath

@CathG 感谢您的想法。只是为了澄清一下；两个观察者独立地观察同一个过程（物体通过），并计算他们所看到的内容。对于每一行，物体的颜色都不同，因此我想知道两个观察者/颜色组合之间是否存在任何统计差异。 - doncarlos

所以我肯定会选择成对测试（而绝对不是行卡方...）但正如@akrun所说，你可以在stats.exchange上问这个问题。 - Cath

Answer 1

以下是使用dplyr的另一种选项：

library(dplyr)

df %>%
  rowwise() %>% 
  mutate(
    test_stat = chisq.test(c(obs.1, obs.2))$statistic,
    p_val = chisq.test(c(obs.1, obs.2))$p.value
    )

- davechilders

Answer 2

您可以使用 "MARGIN =1" 的 apply 函数，然后进行 chisq.test。使用 $statistic 和 $p.value 提取值，并使用 cbind 将其添加到数据集中。

 df1 <- cbind(df, t(apply(df, 1, function(x) {
             ch <- chisq.test(x)
             c(unname(ch$statistic), ch$p.value)})))

 colnames(df1)[3:4] <- c('x-squared', 'p-value')

- akrun

---

这个有效。最初有一些行包含NA，但是通过(na.omit(data))解决了这个问题。 - doncarlos

Answer 3

有多种方法可以做到这一点。其中一种方法是使用 apply 循环遍历每一行 (MARGINE = 1)，然后提取您想要的输出的任何部分（我使用 lapply 遍历每个列表元素）。

xy <- data.frame(obs1 = c(3,12,45,2,7,17,5), obs2 = c(2,10,53,13,12,15,5))
result <- apply(X = xy, MARGIN = 1, FUN = chisq.test)

Warning message:
In FUN(newX[, i], ...) : Chi-squared approximation may be incorrect

# see where p-value is stored
str(chisq.test(xy[1, ]))

List of 9
 $ statistic: Named num 0.2
  ..- attr(*, "names")= chr "X-squared"
 $ parameter: Named num 1
  ..- attr(*, "names")= chr "df"
 $ p.value  : num 0.655 # thar she blows
 $ method   : chr "Chi-squared test for given probabilities"
 $ data.name: chr "xy[1, ]"
 $ observed : num [1:2] 3 2
 $ expected : num [1:2] 2.5 2.5
 $ residuals: num [1:2] 0.316 -0.316
 $ stdres   : num [1:2] 0.447 -0.447
 - attr(*, "class")= chr "htest"

Warning message:
In chisq.test(xy[1, ]) : Chi-squared approximation may be incorrect

unlist(lapply(result, "[", "p.value"), use.names = FALSE)

[1] 0.654720846 0.669815358 0.419020334 0.004508698 0.251349109 0.723673610 1.000000000

- Roman Luštrik