输入是一个方阵,大部分是0
,有些1
。目标是沿输入矩阵的对角线获取连续1
的(某种)累积总和。
#Input
ind = rbind(cbind(x = c(2, 3, 1, 2 , 3),
y = c(1, 2, 3, 4, 5)))
m1 = replace(matrix(0, 5, 5), ind, 1)
m1
# [,1] [,2] [,3] [,4] [,5]
#[1,] 0 0 1 0 0
#[2,] 1 0 0 1 0
#[3,] 0 1 0 0 1
#[4,] 0 0 0 0 0
#[5,] 0 0 0 0 0
#Desired Output
# [,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
我有一个for
循环可以完成这项工作,但是有更好的方法吗?
#Current Approach
m2 = m1
for (i in 2:nrow(m1)){
for (j in 2:nrow(m1)){
if (m1[i-1, j-1] == 1 & m1[i, j] == 1){
m2[i, j] = m2[i - 1, j - 1] + m2[i, j]
m2[i - 1, j - 1] = 0
}
}
}
m2
# [,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
从示例中可以看出,似乎每个对角线都是零,或者是一个 1 后跟零的序列。 我们假设情况总是如此。
首先形成一个函数cum
,该函数采用对角线x
并输出相同长度的零向量,但位置sum(x)
要设置为sum(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
如果对角线上的 1序列不需要从对角线的开头开始,但每个对角线上最多只有一个 1 序列仍然是正确的,那么请使用它代替上面的cum
:
cum <- function(x, s = sum(x)) replace(0 * x, s + which.max(x) - 1, s)
如果对角线上可以有多个序列,请使用它代替上面的cum
:
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