1

我有一个包含大约 20,000 个观察值的数据框。由此我创建了一个列联表,其中包含两个变量的频率。

有了这个,我想对独立性进行卡方检验,看看我的两个变量之间是否存在关系。通常这很容易,但许多单元格的预期值为 0,尽管原始数据框的大小很大。我想删除频率小于 5 的所有行。

我已经广泛搜索了堆栈交换,但我找不到解决这个特定问题的解决方案,我要么a)理解(我对 R 比较陌生),要么 b)使用列联表而不是原始数据框.

非常感谢任何帮助。

编辑:

感谢您的回复贾斯汀。

根据要求,我已经上传了数据框和列联表的摘录。我还上传了到目前为止我尝试过的少量代码,并有结果。

数据框

Department Super
AAP     1
ACS     4
ACE     1
AMA     1
APS     3
APS     2
APS     1
APS     1
ARC     5
ARC     7
ARC     1
BIB     6
BIB     6
BMS     2

所以有两列,第一列是三字母部门代码,第二列是一位整数 (1-7)。

列联表

table(department,super)

        1    2   3   4   5   6   7   8
ACS     32  10   7  24  50   7  24  14
AMA      0   4   2   6  10   3  11   1
...

所以有一个标准的带有频率的列联表。

到目前为止,我知道我可以创建一个逻辑测试来测试单元格内容是否小于 5:

depSupCrosstab <- depSupCrosstab[,2:8]>5

我不知道的是如何使用这行代码创建的矩阵来删除整行(如果它们有任何 FALSE 条目)。

希望有帮助。恐怕我是新手,但只有一种学习方法......

4

2 回答 2

1

恐怕你的问题更复杂。卡方检验的假设是每个单元格的预期频率大于 5。在您的示例中,您试图选择列联表的每个单元格的计数,即观察到的频率预期频率(在零假设下)是根据此处的基本示例所示的行和列总数计算得出的。

按照您的示例,假设的列联表可能如下所示:

ACS <- c(32, 10, 7, 24, 50, 7, 24, 14)
AMA <- c(0, 4, 2, 6, 10, 3, 11, 1)
ARC <- c(6, 10, 12, 3, 12, 23, 10, 2)

tab <- rbind(ACS, AMA, ARC)

如果您筛选观察到的计数等于或小于 5,您将删除 AMA 和 ARC:

apply(tab,1, function(x) any(x<=5))

  ACS   AMA   ARC 
FALSE  TRUE  TRUE 

这在概念上是错误的,因为如上所述,预期频率取决于整个数据。获取 exp。数:

chisq.test(tab, correct=F)$expected

         [,1]      [,2]      [,3]      [,4]      [,5]      [,6]      [,7]
ACS 22.558304 14.247350 12.466431 19.590106 42.742049 19.590106 26.713781
AMA  4.968198  3.137809  2.745583  4.314488  9.413428  4.314488  5.883392
ARC 10.473498  6.614841  5.787986  9.095406 19.844523  9.095406 12.402827
         [,8]
ACS 10.091873
AMA  2.222615
ARC  4.685512

Warning message:
In chisq.test(tab, correct = F): Chi-squared approximation may be incorrect

卡方检验发出警告消息,因为确实有一些带有 exp 的单元格。计数小于 5。但是,如果仅删除 AMA,则表的动态(行和列总计)会更改,并且所有 exp。计数高于 5:

chisq.test(tab[-2,], correct=F)$expected

        [,1]      [,2]     [,3]      [,4]     [,5]      [,6]     [,7]
ACS 25.95122 13.658537 12.97561 18.439024 42.34146 20.487805 23.21951
ARC 12.04878  6.341463  6.02439  8.560976 19.65854  9.512195 10.78049
         [,8]
ACS 10.926829
ARC  5.073171

因此,如果您同时删除 AMA 和 ARC,您将丢失重要信息。


您可以尝试运行Fisher 精确检验(请参阅下面的说明):

fisher.test(tab,simulate.p.value=TRUE,B=10000)

总结:

  1. 单个观察到的频率不能很好地指示预期频率。观察到的频率可能低于 5,但该小区的预期频率将高于 5。
  2. 在大型列联表中,最多可以有 20% 的 exp。频率低于 5,但结果是失去统计功效,因此测试可能无法检测到真正的效果。即使在那种情况下,exp。频率不应低于 1。
  3. 低 exp 的替代测试。频率是 Fisher 的精确检验。卡方检验统计量近似于卡方分布。如果样本量很大,这个近似值会变得更准确,因此需要 exp。频率 > 5。即使样本量很小,Fisher 精确检验也会计算卡方统计量的精确概率,但它可能计算量更大。不幸的是,对于大于 2x2 的列联表,您可能需要模拟 p 值,这有其自身的局限性(这里没有空间讨论它,但它是一个很好的研究主题)。选择大量重复进行模拟 (B),并对其进行调整以查看您的解决方案的稳健性。
于 2013-11-22T16:58:27.990 回答
0

我想我在相关问题中找到了答案。 apply在这种情况下是你的朋友,因为它可以遍历列或行。

要创建与您类似的数据框,然后仅选择所有列 > 5 的行,可以使用以下命令:

set.seed(1985)
tosub <- data.frame(matrix(round(runif(n = 80, min = 0, max = 100)), ncol = 8))
head(tosub,2)
x <- apply(tosub[,1:8] > 5, MARGIN = 1, all)
summary(x)
tosub[which(x),]

   X1 X2 X3 X4 X5 X6 X7 X8
1  66 30 72 59 26 69 76 47
2  27 42 26 95 66 14 67 18
4  42 28 93  7 35 35 95 23
5  38 89 69 91 98 91 60 69
9  89 31 91 72 28 31 58 58
10 53 87 27 89 95 37 98 20
于 2013-11-22T16:18:12.950 回答