矩阵对角线上的累积和

从这个例子来看,似乎每个对角线都是零,否则就是一个由一和零组成的序列。我们假设情况总是这样。

首先形成一个函数 cum 是对角线 x sum(x) 设置为 总和(x) .

ave . row(m1)-col(m1) 在对角线上是常量,可用于分组。

cum <- function(x, s = sum(x)) replace(0 * x, s, s)
ave(m1, row(m1) - col(m1), FUN = cum)

##      [,1] [,2] [,3] [,4] [,5]
## [1,]    0    0    0    0    0
## [2,]    0    0    0    0    0
## [3,]    0    2    0    0    3
## [4,]    0    0    0    0    0
## [5,]    0    0    0    0    0

如果diaongal上的一个序列不需要从对角线的开头开始,但每个对角线上最多只有一个一个序列,则使用此序列代替 附有 上图:

cum <- function(x, s = sum(x)) replace(0 * x, s + which.max(x) - 1, s)

附有 上图:

library(data.table)
cum <- function(x) {
  ave(x, rleid(x), FUN = function(x, s = sum(x)) replace(0 * x, s, s))
}

Rcpp中的循环

library(Rcpp)

cppFunction('NumericMatrix diagcumsum( NumericMatrix m1 ) {

  int i = 0;
  int j = 0;
  int n_row = m1.nrow();

  NumericMatrix res = Rcpp::clone( m1 );

  for( i = 1; i < n_row; i++ ) {
    for( j = 1; j < n_row; j++ ) {
      if( m1( (i-1), (j-1) ) == 1 && m1( i, j ) == 1 ) {
        res(i, j) = res( (i-1), (j-1) ) + res(i, j);
        res( (i-1), (j-1) ) = 0;
      }
    }
  }
  return res;
}')

diagcumsum( m1 )

#      [,1] [,2] [,3] [,4] [,5]
# [1,]    0    0    0    0    0
# [2,]    0    0    0    0    0
# [3,]    0    2    0    0    3
# [4,]    0    0    0    0    0
# [5,]    0    0    0    0    0