矩阵的行和列的比较

Question

假设我们有p by n 矩阵。我想生成一个输出矩阵，w (p x p) 如w_ij 表示i_th 行数比j_th 大多少倍（显然最多可以是n）。

我的代码在这里，我正在寻找更快的方法。

p <- dim(dat)[1]
n <- dim(dat)[2]
w <- matrix(0,p,p)

for(i in 1:n){
     for(j in 1:(p-1)){
         for(k in (j+1):p){
             if(dat[j,i] > dat[k,i]){                     
                   w[j,k] <- w[j,k]+1
             }else{          
                   w[k,j] <- w[k,j]+1
                   }
          }  
      }
  }

---

一个小例子

如果输入数据是

dat <- matrix(1:9, 3)
dat
#     [,1] [,2] [,3]
#[1,]    1    4    7
#[2,]    2    5    8
#[3,]    3    6    9

预期的结果是

W <- matrix(c(0,3,3,0,0,3,0,0,0),3)
W
#     [,1] [,2] [,3]
#[1,]    0    0    0
#[2,]    3    0    0
#[3,]    3    3    0

Answer 1

这似乎可以快速提高速度，而无需太多额外的工作

newd <- t(dat)
for(i in 1:p) { 
  w[,i] <- colSums((newd - dat[i,]) > 0)
} 

快速比较：将代码包装在函数中

f1 <- function(dat){
p <- dim(dat)[1]
n <- dim(dat)[2]
w <- matrix(0,p,p)
for(i in 1:n){
     for(j in 1:(p-1)){
         for(k in (j+1):p){
             if(dat[j,i] > dat[k,i]){                     
                   w[j,k] <- w[j,k]+1
             }else{          
                   w[k,j] <- w[k,j]+1
                   }
          }  
      }
  }
w
}

f2 <- function(dat){
p <- dim(dat)[1]
w <- matrix(0,p,p)
newd <- t(dat)
for(i in 1:p) { 
  w[,i] <- colSums((newd - dat[i,]) > 0)
} ; w}

生成稍大的数据

set.seed(1)
dat <- matrix(rnorm(1e4), 100)

比较

all.equal(f1(dat), f2(dat))

基准测试

library(microbenchmark)
microbenchmark(f1(dat), f2(dat), times=10)
#    expr        min         lq       mean     median         uq        max neval  cld
# f1(dat) 1586.10589 1594.40701 1619.03102 1616.14899 1635.05695 1688.08589    10  b
# f2(dat)   22.56083   23.13493   23.98392   23.34228   24.39766   28.29201    10  a

当然，根据矩阵的大小，可能值得用 c++/Rcpp 编写循环以获得更大的速度增益