【发布时间】：2015-03-25 13:19:53
【问题描述】：
早上好，
我正在尝试使用 R 非常快速地运行 100,000 次 Fisher 对模拟遗传数据的精确测试，最好在 30 秒内（因为我需要置换病例对照标签并迭代该过程 1,000 次，因此它会运行一夜） .
我尝试在已融化、整齐的数据上使用数据表，其中包含大约 200,000,000 行和四列（受试者 ID、疾病状态、位置和“值”[野生型等位基因的数量，一个 3 因素变量]） .该函数按位置分组，然后针对疾病值执行 Fisher 精确检验。
> head(casecontrol3)
ident disease position value
1: 1 0 36044 2
2: 2 0 36044 2
3: 3 0 36044 1
4: 4 0 36044 1
5: 5 0 36044 2
6: 6 0 36044 1
> setkey(casecontrol3,position)
> system.time(casecontrol4 <- casecontrol3[,list(p=fisher.test(value,
+ factor(disease))$p.value), by=position])
user system elapsed
215.430 11.878 229.148
> head(casecontrol4)
position p
1: 36044 6.263228e-40
2: 36495 1.155289e-68
3: 38411 7.842216e-19
4: 41083 1.272841e-69
5: 41866 2.264452e-09
6: 41894 9.833324e-10
但是，与在扁平、凌乱的病例对照表（100,000 行；列包含信息：疾病状态和野生型等位基因的数量）上使用简单的应用函数相比，它确实很慢，所以首先应用函数将每一行转换为 2x3 案例控制表，并使用 Fisher 精确检验的矩阵语法）。将数据从以前的（未融化的）形式转换为这种形式（未显示）大约需要 20 秒的运行时间。
> head(cctab)
control_aa control_aA control_AA case_aa case_aA case_AA
[1,] 291 501 208 521 432 47
[2,] 213 518 269 23 392 585
[3,] 170 499 331 215 628 157
[4,] 657 308 35 269 619 112
[5,] 439 463 98 348 597 55
[6,] 410 480 110 323 616 61
> myfisher <- function(row){
+ contab <- matrix(as.integer(row),nrow=2,byrow=TRUE)
+ pval <- fisher.test(contab)$p.value
+ return(pval)
+ }
> system.time(tab <- apply(cctab,1,"myfisher"))
user system elapsed
28.846 10.989 40.173
> head(tab)
[1] 6.263228e-40 1.155289e-68 7.842216e-19 1.272841e-69 2.264452e-09 9.833324e-10
如您所见，使用 apply 比使用 data.table 快得多，这真的让我感到惊讶。结果完全一样：
> identical(casecontrol4$p,tab)
[1] TRUE
有谁是使用 data.table 的专家，知道如何使用它来加快我的代码速度吗？还是数据太大了，我无法以融化的形式使用它（这排除了使用 data.table、dplyr 等）？请注意，我没有尝试过 dplyr，因为我听说 data.table 对于像这样的大数据集更快。
谢谢。
【问题讨论】：
-
除了速度问题之外，我会与统计学家交谈，以确定运行 100K 测试是否合适而无需进行某种形式的更正。它为我敲响了各种各样的警钟。
-
我不确定您的确切问题，但我在 R 中所做的是使用 R 分析器。它提取了一些堆栈样本，然后我在文本编辑器中手动查看这些样本。你不需要看很多样本。如果它花费了很大一部分时间，比如 50%，做一些可以避免的事情，那么这部分样本就会显示出来。
-
@thelatemail：稍后将进行多次测试更正（例如 Bonferroni）。
-
@AJP123：好的，你是说几乎所有的堆栈样本都包含fisher.test？然后我会做的是看看我是否需要 fisher.test 为我做的一切。它可能在一般性方面犯了错误 - 做的比我需要的多。然后我会看看是否有办法让它做到最低限度。最后，如果这还不够，我会考虑使用编译语言。每个人都希望 R 快，但它的主要目的不是快，而是易于使用。
-
@MikeDunlavey：是的 - 大多数堆栈样本都包含 fisher.test。我会和我的主管聊聊，如果我们需要它更快，我们可能会看看使用 Rcpp。
标签： r performance optimization time