Improve speed of rowwise computations

[frederikziebell]

Hi Henrik,

I think this is somewhat related to #200: I noticed that row-wise computations in matrixStats are slower than column-wise one and that this effect worsens with the size of the matrix. I think this could be improved in a similar fashion to the genefilter package, in which the row-wise t-test is about 10% slower than the column-wise for large matrices.

Here's some example code:

library("matrixStats")
library("genefilter")
library("bench")
library("tidyverse")

row_col_ratios <- function(n){
  
  set.seed(1)
  mat <- matrix(rnorm(n*n), ncol=n)
  
  res_means <- bench::mark(
    rowMeans2(mat),
    colMeans2(mat),
    iterations = 10,
    check = F
  )
  res_ttests <- bench::mark(
    rowttests(mat),
    colttests(mat),
    iterations = 10,
    check = F
  )
  
  
  data.frame(
    n,
    row_col_means = as.double(res_means$median[1])/as.double(res_means$median[2]),
    row_col_ttests = as.double(res_ttests$median[1])/as.double(res_ttests$median[2])
  )

}

c(1e2,2e2,5e2,1e3,2e3,5e3,1e4,2e4) %>% 
  map(row_col_ratios) %>% 
  bind_rows() %>% 
  pivot_longer(-n) %>% 
  ggplot(aes(n, value, color = name)) +
    geom_point() +
    geom_path() +
    labs(y = "runtime rowwise / runtime colwise")

[yaccos]

R stores matrices in a column-major order. This means that columns are stored continiously in memory, whereas rows are striped across the memory with a fixed offset. The column- and row-wise operations in matrixStats are currently implemented such that the result in one column or row is computed before proceeding to the next column or row. For columnwise operation, this is not an issue as the functions are looping over continious blocks of memory. However, for rowwise operations, data are retrieved from non-adjecent memory locations. This creates a memory bottleneck because the CPU has to wait for data to be retrieved from memory. For matrices small enough to fit in cashe the performance penalty is relatively small because the memory access from cache is quick. In the cases where the matrix is too big to fit in the cashe, the CPU must wait for data to be delivered from primary memory which is considerably slower, resulting in a larger performance penalty for rowwise functions.

[yaccos]

The issue could most likely be resolved by fine-tuning the code for rowwise operations to iterate through all rows at once and store the temporary results in an array. This would require a non-trivial effort and increase the complexity of the codebase (currently, the same macro template files are use both for column- and rowwise computations). I personally don`t think this will be implemented anytime soon due to the effort required.