BioQC算法:加速Wilcoxon-Mann-Whitney检验

Gregor Sturm和jiao David Zhang

2022-11-01

在这个小插图中,我们解释了我们实现Wilcoxon-Mann-Whitney检验的底层算法细节。用于生成本文档的源代码可以在github存储库中找到BioQC

BioQC是一个R/Bioconductor包,用于从高通量基因表达谱数据中检测组织异质性。它实现了一种高效的Wilcoxon-Mann-Whitney测试,并提供了组织特异性基因签名,可以“开箱即用”。

算法的改进

Wilcoxon-Mann-Whitney (WMW)检验是一种非参数统计检验,用于检验两个总体的中位数是否相等。与t检验不同,它不需要正态分布的假设,这使得它对噪声更健壮。

与原生R实现相比,我们改进了Wilcoxon-Mann-Whitney检验的计算效率,基于以下三个修改:

  1. 近似wmw统计量(Zar, J. H.(1999)。Biostatistical分析。培生教育印度公司。页。174-177)。对于高通量基因表达数据,与准确版本的差异可以忽略不计。
  2. 核心算法是用C语言实现的,而不是用R作为编程语言。
  3. BioQC避免了徒劳而昂贵的分拣操作。

虽然(1)和(2)是直截了当的,但我们在下面详细说明(3)。

\ (W_ {a、b} \)为两个基因载体的近似WMW检验\ (a、b \),在那里\ \ ()感兴趣的基因集,通常包含少于几百个基因,以及\ (b \)是基因集以外的所有基因的集合(背景基因)通常包含\ (> 10000 \)基因。在生物质量控制的背景下,基因集被称为组织签名

给定一个\(m \乘n\)输入矩阵的基因表达数据\ \(米)基因和\ (n \)样品\(s_1 \dots, s_n\),\ (k \)基因集\(d_1, \dots, d_k\),每个样品都需要进行WMW-test\(s_i, I \in 1..n\)每个基因集\(d_j, j \in 1..k\).WMW-test的运行时间基本上是由对两个输入向量的排序操作决定的。使用原生Rwilcox.test,向量\ \ ()\ (b \)分别对每个基因进行排序。然而,在基因集分析的背景下,这是徒劳的,因为当在同一样本上测试不同的基因集时,(大)背景集相对于(小)基因集的变化不显著。

因此,我们通过扩展来近似WMW-test\ (b \)样本中的所有基因,在测试多个基因集时保持背景不变。像这样,\ (b \)每个样本只需要排序一次。单个基因集仍然需要排序,这不是一个主要问题,因为与背景基因集相比,它们很小。

BioQC通过避免对同一样品进行无用的分类操作,加快了Wilcoxon-Mann-Whitney测试的速度。

BioQC通过避免对同一样品进行无用的分类操作,加快了Wilcoxon-Mann-Whitney测试的速度。

时间基准

为了证明BioQC的卓越性能,我们同时应用BioQC和原生Rwilcox.test随机表达式矩阵并测量运行时间。

我们建立了77510个人类基因的1、5、10、50或100个样本的随机表达矩阵。基因是\(先验知识\)分布式后\ (\ mathcal {N} (0,1) \).原生的R和BioQC矩阵分别应用了Wilcoxon-Mann-Whitney检验的实现。

两种实现的数值结果,bioqcNumRes(从BioQC),rNumRes(从R),是等价的,如下条命令所示。

expect_equal(bioqcNumRes rNumRes)

BioQC实现速度超过500倍:BioQC计算100个样本中所有155个签名的富集分数大约需要1秒,而原生R实现大约需要20分钟:

**图2**:BioQC和R实施Wilcoxon-Mann-Whitney检验的时间基准结果。左面板:以秒为单位的运行时间(对数y轴)。右面板:两个实现占用时间的比率。所有结果都是在RedHat Linux服务器上的一个线程上实现的。

图2: BioQC时间基准结果,R实施Wilcoxon-Mann-Whitney检验。左面板:以秒为单位的运行时间(对数y轴)。右面板:两个实现占用时间的比率。所有结果都是在RedHat Linux服务器上的一个线程上实现的。

结论

我们已经证明了BioQC在两个数量级的时间内实现与本机实现相同的结果。这使得BioQC一种大规模高通量基因表达数据质量控制的高效工具。

R会话信息

sessionInfo()
## R版本4.2.1(2022-06-23)##平台:x86_64-pc-linux-gnu(64位)##运行在Ubuntu 20.04.5 LTS ## ##矩阵产品:默认## BLAS: /home/biocbuild/bbs-3.16-bioc/R/lib/libRblas。/home/biocbuild/bbs-3.16-bioc/R/lib/libRlapack。所以## ## locale: ## [1] LC_CTYPE=en_US。UTF-8 LC_NUMERIC= c# # [3] LC_TIME=en_GB LC_COLLATE= c# # [5] LC_MONETARY=en_US。utf - 8 LC_MESSAGES = en_US。UTF-8 ## [7] LC_PAPER=en_US。UTF-8 LC_NAME= c# # [9] LC_ADDRESS=C lc_phone = c# # [11] LC_MEASUREMENT=en_US。UTF-8 LC_IDENTIFICATION=C ## ##附加的基本包:## [1]stats4 stats graphics grDevices utils datasets methods ##[8]基础## ##其他附加包:[7] AnnotationDbi_1.60.0 IRanges_2.32.0 S4Vectors_0.36.0 ## [10] testthat_3.1.5 limma_3.54.0 RColorBrewer_1.1-3 ## [13] BioQC_1.26.0 Biobase_2.58.0 BiocGenerics_0.44.0 ## [16] knitr_1.40 ## ##通过命名空间加载(并且没有附加):# # [1] Rcpp_1.0.9 locfit_1.5 - 9.6 deldir_1.0-6 # # [4] png_0.1-7 Biostrings_2.66.0 gtools_3.9.3 # # [7] rprojroot_2.0.3 digest_0.6.30 R6_2.5.1 # # [10] GenomeInfoDb_1.34.0 RSQLite_2.2.18 evaluate_0.17 # # [13] httr_1.4.4 highr_0.9 zlibbioc_1.44.0 # # [16] rlang_1.0.6 jquerylib_0.1.4 blob_1.2.3 # # [19] rmarkdown_2.17 desc_1.4.2 stringr_1.4.1 # # [22] rcurl_1.98 - 1.9 bit_4.0.4 compiler_4.2.1 # # [25] xfun_0.34 pkgconfig_2.0.3 htmltools_0.5.3 # # [28] KEGGREST_1.38.0 GenomeInfoDbData_1.2.9 edgeR_3.40.0 # #[31] crayon_1.5.2 bitops_1.0-7 brio_1.1.3 ## [34] grid_4.2.1 jsonlite_1.8.3 gtable_0.3.1 ## [37] DBI_1.1.3 magrittr_2.0.3 KernSmooth_2.23-20 ## [40] cli_3.4.1 stringi_1.7.8 cachem_1.0.6 ## [43] XVector_0.38.0 bslib_0.4.0 vctrs_0.5.0 ## [46] tools_4.2.1 interp_1.1-3 bit64_4.0.5 ## [52] yaml_2.3.6 caTools_1.18.2 memoise_2.0.1 ## [55] sass_0.4.2