Project import generated by Copybara.

GitOrigin-RevId: 2f667fe27fbc2bc80379001010bf4fca9703b8b9

Author: scanisius

CHANGES.rst

@@ -21,13 +21,48 @@ of the Python package would have features that do not exist in version
 0.9 series
 ==========
 
-:Latest Python version:  0.9.3 (June 18, 2019)
-:Latest R version:       0.9.3 (June 18, 2019)
+:Latest Python version:  0.9.4 (July 28, 2021)
+:Latest R version:       0.9.4 (July 28, 2021)
 
 
 Python package
 --------------
 
+0.9.4 (July 28, 2021)
+~~~~~~~~~~~~~~~~~~~~~
+
+**Added**:
+
+- The new fdr_method argument of the pairwise_discover_test function
+  selects the false discovery rate estimation method used for multiple
+  testing correction. Its default value selects a Benjamini-Hochberg
+  procedure adapted for discrete test statistics. As an alternative,
+  the standard Benjamini-Hochberg procedure can be selected. While the
+  latter is much faster, it is also more conservative than the
+  discrete version. The discrete method was always used in previous
+  versions and is still the recommended choice.
+
+**Changed**:
+
+- Speed improvements in DiscoverMatrix and
+  pairwise_discover_test. DiscoverMatrix finishes instantaneously for
+  small to medium-sized data sets and takes only a few seconds for
+  very large data sets. pairwise_discover_test (with discrete
+  Benjamini-Hochberg) is about twice as fast as in previous versions.
+
+- Python 2.7 is no longer supported.
+
+**Fixed**:
+
+- If pandas >= 1.0 was installed, subsetting a DiscoverMatrix object
+  gave rise to "AttributeError: 'DataFrame' objects has no attribute
+  'ix'".
+
+- Fixed the underlying Fortran code to make it compile with GNU
+  Fortran >= 10. Due to these changes the minimum supported GNU
+  Fortran version is now version 5.
+
+
 0.9.3 (June 18, 2019)
 ~~~~~~~~~~~~~~~~~~~~~
 
@@ -68,6 +103,41 @@ First public release
 R package
 ---------
 
+0.9.4 (July 28, 2021)
+~~~~~~~~~~~~~~~~~~~~~
+
+**Added**:
+
+- The new fdr.method argument of the pairwise.discover.test function
+  selects the false discovery rate estimation method used for multiple
+  testing correction. Its default value selects a Benjamini-Hochberg
+  procedure adapted for discrete test statistics. As an alternative,
+  the standard Benjamini-Hochberg procedure can be selected. While the
+  latter is much faster, it is also more conservative than the
+  discrete version. The discrete method was always used in previous
+  versions and is still the recommended choice.
+
+**Changed**:
+
+- Speed improvements in discover.matrix and
+  pairwise.discover.test. discover.matrix finishes instantaneously for
+  small to medium-sized data sets and takes only a few seconds for
+  very large data sets. pairwise.discover.test (with discrete
+  Benjamini-Hochberg) is about twice as fast as in previous versions.
+
+- The matrixStats package is no longer a dependency.
+
+**Fixed**:
+
+- Re-enabled OpenMP, which distributes the computations in
+  pairwise.discover.test across multiple CPUs. OpenMP was disabled
+  since R 4.0.
+
+- Fixed the underlying Fortran code to make it compile with GNU
+  Fortran >= 10. Due to these changes the minimum supported GNU
+  Fortran version is now version 5.
+
+
 0.9.3 (June 18, 2019)
 ~~~~~~~~~~~~~~~~~~~~~
 

R/.Rbuildignore

@@ -0,0 +1 @@
+^\.bumpversion\.cfg$

R/.bumpversion.cfg

@@ -3,8 +3,8 @@ commit = True
 message = [R package] Bump version: {current_version} → {new_version}
 tag = True
 tag_name = r_v{new_version}
-current_version = 0.9.3
-parse = (?P<major>\d+)\.(?P<minor>\d+)\.(?P<patch>\d+)(\-(?P<release>[a-z]+))?
+current_version = 0.9.4
+parse = (?P<major>\d+)\.(?P<minor>\d+)\.(?P<patch>\d+)(\-(?P<release>release|9000))?
 serialize = 
 	{major}.{minor}.{patch}-{release}
 	{major}.{minor}.{patch}
@@ -14,6 +14,5 @@ serialize =
 [bumpversion:part:release]
 optional_value = release
 values = 
-	dev
+	9000
 	release
-

R/DESCRIPTION

@@ -1,13 +1,11 @@
 Package: discover
-Version: 0.9.3
-Date: 2019-06-18
+Version: 0.9.4
+Date: 2021-07-27
 Title: DISCOVER mutual exclusivity and co-occurrence analysis
 Author: Sander Canisius
 Maintainer: Sander Canisius <s.canisius@nki.nl>
-Description: DISCOVER method for mutual exclusivity detection in genomic data.
+Description: DISCOVER method for mutual exclusivity and co-occurrence detection in genomic data.
 License: Apache License (== 2.0)
-Depends:
-    matrixStats
 RoxygenNote: 5.0.0
 VignetteBuilder: knitr
-Suggests: knitr
+Suggests: rmarkdown

R/R/analysis.R

@@ -1,15 +1,22 @@
 #' Perform many pairwise mutual exclusivity or co-occurrence tests.
 #'
 #' @param x A \code{discover.matrix} object.
-#' @param g An optional grouping vector for the rows of \code{events}. Pairs of rows within the same group are not tested.
-#' @param alternative Either \code{"less"} for mutual-exclusivity analysis, or \code{"greater"} for co-occurrence analysis.
-#' @param correct If \code{TRUE}, multiple testing correction is performed.
+#' @param g An optional grouping vector for the rows of \code{events}. Pairs of rows
+#'   within the same group are not tested.
+#' @param alternative Either \code{"less"} for mutual-exclusivity analysis, or
+#'   \code{"greater"} for co-occurrence analysis.
+#' @param fdr.method The false discovery rate procedure used for multiple testing correction.
+#'   If \code{"DBH"}, a Benjamini-Hochberg procedure adapted for discrete test statistics is
+#'   used. If \code{"BH"}, the standard Benjamini-Hochberg procedure is used. The latter is
+#'   much faster, but also more conservative than the discrete version.
 #' @return An object of type \code{pairwise.discover.out}.
 #'
 #' @useDynLib discover
 #' @export
-pairwise.discover.test <- function (x, g=NULL, alternative=c("less", "greater"), correct=TRUE) {
+pairwise.discover.test <- function (x, g=NULL, alternative=c("less", "greater"), fdr.method=c("DBH", "BH")) {
   alternative <- match.arg(alternative)
+  fdr.method <- match.arg(fdr.method)
+  discrete.fdr <- fdr.method == "DBH"
 
   events <- x$events
   bg <- x$bg
@@ -18,9 +25,16 @@ pairwise.discover.test <- function (x, g=NULL, alternative=c("less", "greater"),
     result <- .Fortran("mutex",
                        as.integer(nrow(events)), as.integer(ncol(events)),
                        as.integer(events), as.double(bg), as.integer(alternative == "less"),
+                       as.integer(discrete.fdr),
                        p=double(nrow(events) * (nrow(events) - 1) / 2),
                        q=double(nrow(events) * (nrow(events) - 1) / 2),
                        pi0=double(1))
+
+    if (fdr.method == "BH") {
+      result$q <- p.adjust(result$p, "BH")
+      result$pi0 <- 1.0
+    }
+    
     p <- matrix(NA, nrow(events), nrow(events))
     p[lower.tri(p)] <- result$p
 
@@ -34,8 +48,14 @@ pairwise.discover.test <- function (x, g=NULL, alternative=c("less", "greater"),
                        as.integer(nrow(events)), as.integer(ncol(events)),
                        as.integer(events[i, ]), as.double(bg[i, ]), as.integer(alternative == "less"),
                        as.integer(length(block.sizes)), as.integer(block.sizes),
+                       as.integer(discrete.fdr),
                        p=double(nrow(events)**2), q=double(nrow(events)**2), pi0=double(1))
 
+    if (fdr.method == "BH") {
+      result$q <- p.adjust(result$p, "BH")
+      result$pi0 <- 1.0
+    }
+
     j <- order(i)
     p <- matrix(result$p, nrow=nrow(events))[j, j]
     p[is.nan(p)] <- NA
@@ -54,7 +74,8 @@ pairwise.discover.test <- function (x, g=NULL, alternative=c("less", "greater"),
     p.values=p,
     q.values=q,
     pi0=result$pi0,
-    alternative=alternative)
+    alternative=alternative,
+    fdr.method=fdr.method)
   class(result) <- "pairwise.discover.out"
   result
 }

R/R/background.R

@@ -8,6 +8,7 @@
 #' x <- matrix(as.integer(runif(100) < 0.3), 10, 10)
 #' bg <- estimateBackground(x)
 #'
+#' @noRd
 #' @useDynLib discover
 estimateBackground <- function (events, strata=NULL) {
   if (is.null(strata)) {
@@ -27,7 +28,7 @@ estimateBackground <- function (events, strata=NULL) {
                    mu=double(length(row.values) + length(col.values)), PACKAGE="discover")$mu
 
     num.rows <- length(row.values)
-    eA <- exp(mu[1:num.rows] / row.weights) %*% t(exp(mu[(num.rows+1):length(mu)] / col.weights))
+    eA <- exp(mu[1:num.rows]) %*% t(exp(mu[(num.rows+1):length(mu)]))
     bg <- 1.0 / (eA + 1)
 
     bg <- bg[row.inverse, col.inverse]

R/R/grouptests.R

@@ -13,7 +13,7 @@ groupwise.discover.test.coverage <- function (events, bg) {
 groupwise.discover.test.exclusivity <- function (events, bg) {
   log.p <- log(bg)
   log.not.p <- log1p(-bg)
-  p.exactly.one <- exp(matrixStats::colLogSumExps(rep.row(colSums(log.not.p), nrow(events)) + log.p - log.not.p))
+  p.exactly.one <- exp(colLogSumExps(rep.row(colSums(log.not.p), nrow(events)) + log.p - log.not.p))
   x <- sum(colSums(events) == 1)
   ppoisbinom(x, p.exactly.one, lower.tail=FALSE)
 }
@@ -22,7 +22,7 @@ groupwise.discover.test.exclusivity <- function (events, bg) {
 groupwise.discover.test.impurity <- function (events, bg) {
   log.p <- log(bg)
   log.not.p <- log1p(-bg)
-  p.exactly.one <- exp(matrixStats::colLogSumExps(rep.row(colSums(log.not.p), nrow(events)) + log.p - log.not.p))
+  p.exactly.one <- exp(colLogSumExps(rep.row(colSums(log.not.p), nrow(events)) + log.p - log.not.p))
   p.none <- exp(colSums(log.not.p))
   p.more.than.one <- 1 - p.none - p.exactly.one
   x <- sum(colSums(events) > 1)

R/R/logspace.R

@@ -0,0 +1,9 @@
+#' Compute the logarithm of the sum of exponentials for each column in a matrix.
+#'
+#' @param x a real-valued matrix.
+#' @return a vector with elements corresponding to the log sums of exponentials of the columns of \code{x}.
+#'
+#' @useDynLib discover
+colLogSumExps <- function (x) {
+  .Fortran("colLogSumExps", as.double(x), as.integer(nrow(x)), as.integer(ncol(x)), result=double(ncol(x)), NAOK=TRUE)$result
+}

R/R/print.R

@@ -4,6 +4,10 @@ print.pairwise.discover.out <- function (x, fdr.threshold=0.01) {
   cat("Pairwise DISCOVER", type, "test\n")
   cat("alternative hypothesis: observed overlap is", x$alternative, "than expected by chance\n")
 
+  fdr.method <- switch(x$fdr.method, DBH="discrete Benjamini-Hochberg", BH="Benjamini-Hochberg")
+  cat("FDR estimation method:", fdr.method)
+  cat("\n")
+
   cat("\n")
 
   cat("number of pairs tested:", sum(!is.na(x$p.values)))

R/src/Makevars

@@ -1,12 +1,11 @@
-PKG_FCFLAGS = -ffree-line-length-none $(SHLIB_OPENMP_FCFLAGS)
-PKG_LIBS = $(SHLIB_OPENMP_CFLAGS) $(BLAS_LIBS)
+PKG_FCFLAGS = -ffree-line-length-none $(SHLIB_OPENMP_FFLAGS)
+PKG_LIBS = $(SHLIB_OPENMP_FFLAGS) $(BLAS_LIBS)
 
 $(SHLIB) : rwrappers.o
 
 rwrappers.o : maxent.o fdr.o types.o poisbinom.o
 maxent.o : types.o
-fdr.o : poisbinom.o util.o types.o ieee_constants.o
+fdr.o : poisbinom.o util.o types.o
 poisbinom.o : logspace.o types.o
-logspace.o : ieee_constants.o types.o
+logspace.o : types.o
 util.o : uniinv.o types.o
-ieee_constants.o : types.o