tests/rowQuantiles.R

In HenrikBengtsson/matrixStats: Functions that Apply to Rows and Columns of Matrices (and to Vectors)

library("matrixStats")

## Create isFALSE() if running on an old version of R
if (!exists("isFALSE", mode="function")) {
  isFALSE <- function(x) is.logical(x) && length(x) == 1L && !is.na(x) && !x
}

rowQuantiles_R <- function(x, probs = c(0, 0.25, 0.50, 0.75, 1), na.rm = FALSE, drop = TRUE, type = 7L, ..., useNames = NA) {
  q <- apply(x, MARGIN = 1L, FUN = function(x, probs, na.rm) {
    if (!na.rm && any(is.na(x))) {
      na_value <- NA_real_
      if (type != 7L) storage.mode(na_value) <- storage.mode(x)
      rep(na_value, times = length(probs))
    } else {
      as.vector(quantile(x, probs = probs, na.rm = na.rm, type = type, names = FALSE, ...))
    }
  }, probs = probs, na.rm = na.rm)

  if (!is.null(dim(q))) q <- t(q)
  else dim(q) <- c(nrow(x), length(probs))

  colnames(q) <- matrixStats:::quantile_probs_names(probs)
  rownames(q) <- rownames(x)
  if (isFALSE(useNames)) dimnames(q) <- NULL

  if (drop) q <- drop(q)
  q
}


# - - - - - - - - - - - - - - - - - - - - - - - - - - - -
# Test with multiple quantiles
# - - - - - - - - - - - - - - - - - - - - - - - - - - - -
for (mode in c("logical", "integer", "double")) {
  cat("mode: ", mode, "\n", sep = "")
  x <- matrix(1:40 + 0.1, nrow = 8, ncol = 5)
  storage.mode(x) <- mode
  dimnames <- lapply(dim(x), FUN = function(n) letters[seq_len(n)])
  str(x)

  probs <- c(0, 0.5, 1)
  # Test with and without dimnames on x
  for (setDimnames in c(TRUE, FALSE)) {
    if (setDimnames) dimnames(x) <- dimnames
    else dimnames(x) <- NULL    
    # Check names attribute
    for (useNames in c(if (!matrixStats:::isUseNamesNADefunct()) NA, TRUE, FALSE)) {
      q0 <- rowQuantiles_R(x, probs = probs, useNames = useNames)
      print(q0)
      q1 <- rowQuantiles(x, probs = probs, useNames = useNames)
      print(q1)
      ## FIXME: Workaround for R (< 3.0.0)
      if (getRversion() < "3.0.0" && mode == "logical") storage.mode(q1) <- storage.mode(q0)
      stopifnot(all.equal(q1, q0))
      q2 <- colQuantiles(t(x), probs = probs, useNames = useNames)
      ## FIXME: Workaround for R (< 3.0.0)
      if (getRversion() < "3.0.0" && mode == "logical") storage.mode(q2) <- storage.mode(q0)
      stopifnot(all.equal(q2, q0))      
    }
  }
} # for (mode ...)


# - - - - - - - - - - - - - - - - - - - - - - - - - - - -
# Test with a single quantile
# - - - - - - - - - - - - - - - - - - - - - - - - - - - -
for (mode in c("logical", "integer", "double")) {
  cat("mode: ", mode, "\n", sep = "")
  x <- matrix(1:40, nrow = 8, ncol = 5)
  storage.mode(x) <- mode
  dimnames <- lapply(dim(x), FUN = function(n) letters[seq_len(n)])
  str(x)

  probs <- c(0.5)
  # Test with and without dimnames on x
  for (setDimnames in c(TRUE, FALSE)) {
    if (setDimnames) dimnames(x) <- dimnames
    else dimnames(x) <- NULL    
    # Check names attribute
    for (useNames in c(if (!matrixStats:::isUseNamesNADefunct()) NA, TRUE, FALSE)) {
      q0 <- rowQuantiles_R(x, probs = probs, useNames = useNames)
      print(q0)
      q1 <- rowQuantiles(x, probs = probs, useNames = useNames)
      ## FIXME: Workaround for R (< 3.0.0)
      if (getRversion() < "3.0.0" && mode == "logical") storage.mode(q1) <- storage.mode(q0)
      print(q1)
      stopifnot(all.equal(q1, q0))
      q2 <- colQuantiles(t(x), probs = probs, useNames = useNames)
      ## FIXME: Workaround for R (< 3.0.0)
      if (getRversion() < "3.0.0" && mode == "logical") storage.mode(q2) <- storage.mode(q0)
      stopifnot(all.equal(q2, q0))      
    }
  }
} # for (mode ...)


# - - - - - - - - - - - - - - - - - - - - - - - - - - - -
# Consistency checks
# - - - - - - - - - - - - - - - - - - - - - - - - - - - -
set.seed(1)

probs <- seq(from = 0, to = 1, by = 0.25)

cat("Consistency checks:\n")
n_sims <- if (Sys.getenv("_R_CHECK_USE_VALGRIND_") != "") 6L else 24L
for (kk in seq_len(n_sims)) {
  cat("Random test #", kk, "\n", sep = "")

  # Simulate data in a matrix of any shape
  dim <- sample(20:60, size = 2L)
  n <- prod(dim)
  x <- rnorm(n, sd = 100)
  dim(x) <- dim
  dimnames <- lapply(dim(x), FUN = function(n) rep(letters, length.out = n))

  # Add NAs?
  has_na <- ((kk %% 2) == 0L)
  if (has_na) {
    cat("Adding NAs\n")
    nna <- sample(n, size = 1)
    na_values <- c(NA_real_, NaN)
    t <- sample(na_values, size = nna, replace = TRUE)
    x[sample(length(x), size = nna)] <- t
  }

  # Logical, integer, or double?
  mode <- "numeric"
  if ((kk %% 6) %in% 1:2) {
    cat("Coercing to logical\n")
    mode <- "logical"
  } else if ((kk %% 6) %in% 3:4) {
    cat("Coercing to integers\n")
    mode <- "integer"
  }
  storage.mode(x) <- mode
  
  str(x)

  # rowQuantiles():
  for (type in 1:9) {
    cat(sprintf("type=%d, has_na=%s:\n", type, has_na))
    # Test with and without dimnames on x
    for (setDimnames in c(TRUE, FALSE)) {
      if (setDimnames) dimnames(x) <- dimnames
      else dimnames(x) <- NULL    
      # Check names attribute
      for (useNames in c(if (!matrixStats:::isUseNamesNADefunct()) NA, TRUE, FALSE)) {
        q0 <- rowQuantiles_R(x, probs = probs, na.rm = has_na, type = type, useNames = useNames)
        q1 <- rowQuantiles(x, probs = probs, na.rm = has_na, type = type, useNames = useNames)
        ## FIXME: Workaround for R (< 3.0.0)
        if (getRversion() < "3.0.0" && mode == "logical" && !has_na && type == 7L) storage.mode(q1) <- storage.mode(q0)
        stopifnot(all.equal(q1, q0))
        q2 <- colQuantiles(t(x), probs = probs, na.rm = has_na, type = type, useNames = useNames)
        ## FIXME: Workaround for R (< 3.0.0)
        if (getRversion() < "3.0.0" && mode == "logical" && !has_na && type == 7L) storage.mode(q2) <- storage.mode(q0)
        stopifnot(all.equal(q2, q0))
      }
    }
  }
} # for (kk ...)



for (mode in c("logical", "integer", "double")) {
  naValue <- NA_real_
  storage.mode(naValue) <- mode

  someValue <- 1
  storage.mode(someValue) <- mode

  for (type in 1:9) {
    # - - - - - - - - - - - - - - - - - - - - - - - - - - - -
    # All NA
    # - - - - - - - - - - - - - - - - - - - - - - - - - - - -
    x <- matrix(naValue, nrow = 3L, ncol = 4L)
    dimnames <- lapply(dim(x), FUN = function(n) letters[seq_len(n)])
    # Test with and without dimnames on x
    for (setDimnames in c(TRUE, FALSE)) {
      if (setDimnames) dimnames(x) <- dimnames
      else dimnames(x) <- NULL    
      # Check names attribute
      for (useNames in c(if (!matrixStats:::isUseNamesNADefunct()) NA, TRUE, FALSE)) {
        qr0 <- rowQuantiles_R(x, type = type, useNames = useNames)
        
        qr <- rowQuantiles(x, type = type, useNames = useNames)
        stopifnot(identical(qr, qr0))
        
        # x <- matrix(naValue, nrow = 4L, ncol = 3L)
        qc <- colQuantiles(t(x), type = type, useNames = useNames)
        
        stopifnot(identical(qc, qr))
      }
    }


    # - - - - - - - - - - - - - - - - - - - - - - - - - - - -
    # Empty matrices
    # - - - - - - - - - - - - - - - - - - - - - - - - - - - -
    probs <- c(0, 0.25, 0.75, 1)
    x <- matrix(naValue, nrow = 0L, ncol = 0L)
    dimnames(x) <- lapply(dim(x), FUN = function(n) letters[seq_len(n)])
    q <- rowQuantiles(x, probs = probs, type = type, useNames = TRUE)
    stopifnot(identical(dim(q), c(nrow(x), length(probs))))
    q <- colQuantiles(x, probs = probs, type = type, useNames = TRUE)
    stopifnot(identical(dim(q), c(ncol(x), length(probs))))
    
    x <- matrix(naValue, nrow = 2L, ncol = 0L)
    dimnames(x) <- lapply(dim(x), FUN = function(n) letters[seq_len(n)])
    q <- rowQuantiles(x, probs = probs, type = type, useNames = TRUE)
    stopifnot(identical(dim(q), c(nrow(x), length(probs))))
    
    x <- matrix(naValue, nrow = 0L, ncol = 2L)
    dimnames(x) <- lapply(dim(x), FUN = function(n) letters[seq_len(n)])
    q <- colQuantiles(x, probs = probs, type = type, useNames = TRUE)
    stopifnot(identical(dim(q), c(ncol(x), length(probs))))


    # - - - - - - - - - - - - - - - - - - - - - - - - - - - -
    # Single column matrices
    # - - - - - - - - - - - - - - - - - - - - - - - - - - - -
    probs <- c(0, 0.25, 0.75, 1)
    x <- matrix(someValue, nrow = 2L, ncol = 1L)
    dimnames(x) <- lapply(dim(x), FUN = function(n) letters[seq_len(n)])
    qr <- rowQuantiles(x, probs = probs, type = type, useNames = TRUE)
    print(qr)
    
    x <- matrix(someValue, nrow = 1L, ncol = 2L)
    dimnames(x) <- lapply(dim(x), FUN = function(n) letters[seq_len(n)])
    qc <- colQuantiles(x, probs = probs, type = type, useNames = TRUE)
    print(qc)
    
    stopifnot(identical(qc, qr))
  }
}