将矩阵次对角线转换为列 r 的高效算法

Question

我有一个非方阵，需要对其次对角线进行一些计算。我发现最好的方法是将次对角线也变成 columns/rows 并使用像 cumprod 这样的函数。现在我使用 for 循环和 exdiag 定义如下：

exdiag <- function(mat, off=0) {mat[row(mat) == col(mat)+off]}

但是效率不高。你知道任何其他算法来实现这种结果吗？

一个小例子来展示我在做什么：

exdiag <- function(mat, off=0) {mat[row(mat) == col(mat)+off]}
mat <- matrix(1:72, nrow = 12, ncol = 6)
newmat <- matrix(nrow=11, ncol=6)
for (i in 1:11){
   newmat[i,] <- c(cumprod(exdiag(mat,i)),rep(0,max(6-12+i,0)))
}

此致，阿图尔

Answer 1

您可以修改 diag() 函数。

exdiag <- function(mat, off=0) {mat[row(mat) == col(mat)+off]}

exdiag2 <- function(matrix, off){diag(matrix[-1:-off,])}

速度测试：

mat = diag(10, 10000,10000)

off = 4

> system.time(exdiag(mat,4))

  user  system elapsed 
  7.083   2.973  10.054 

> system.time(exdiag2(mat,4))

   user  system elapsed 
  5.370   0.155   5.524 

> system.time(diag(mat))

   user  system elapsed 
  0.002   0.000   0.002

看起来矩阵的子集化需要很多时间，但它仍然比你的实现表现得更好。可能有很多其他的子集化方法，它们优于我的解决方案。 :)

Answer 2

最快但迄今为止最神秘的解决方案是从非方阵中获取所有可能的对角线，将矩阵视为向量并简单地为 [=37= 构造一个 id 向量]离子。最后你可以根据需要将它转换回矩阵。

以下函数执行此操作：

exdiag <- function(mat){

  NR <- nrow(mat)
  NC <- ncol(mat)
  smalldim <- min(NC,NR)
  if(NC > NR){
    id <- seq_len(NR) + 
      seq.int(0,NR-1)*NR +
      rep(seq.int(1,NC - 1), each = NR)*NR


  } else if(NC < NR){
    id <- seq_len(NC) + 
      seq.int(0,NC-1)*NR +
      rep(seq.int(1,NR - 1), each = NC)

  } else {
    return(diag(mat))
  }
  out <- matrix(mat[id],nrow = smalldim)
  id <- (ncol(out) + 1 - row(out)) - col(out) < 0
  out[id] <- NA
  return(out)
}

请记住，您必须考虑矩阵的形成方式。

在这两种情况下，我都遵循相同的逻辑：

首先构建一个序列，指示沿最小维度的位置
向该序列添加 0、1、2 ... 乘以行长度。

这将创建第一个对角线。这样做之后，您只需添加一个序列，将整个先前序列移动 1（向下或向右），直到到达矩阵的末尾。要右移，我需要将该序列乘以行数。

最后，您可以使用这些索引 select 来自 mat 的正确位置，并且 return 所有这些作为矩阵。由于此代码的矢量化性质，您必须检查最后的对角线是否正确。它们包含的元素比第一个少，因此您必须用 NA 替换不属于该次对角线的值。同样在这里，您可以简单地使用索引技巧。

您可以按如下方式使用：

> diag1 <- exdiag(amatrix)
> diag2 <- exdiag(t(amatrix))
> identical(diag1, diag2)
[1] TRUE

为了得出你的结果

amatrix <- matrix(1:72, ncol = 6)
diag1 <- exdiag(amatrix)
res <- apply(diag1,2,cumprod)
res[is.na(res)] <- 0
t(res)

将矩阵次对角线转换为列 r 的高效算法

Efficient algorithm to turn matrix subdiagonal to columns r

algorithm

r

transformation

matrix