在R中高效实现积分图/累加和表的方法

Question

我正在尝试构建一个积分图或者叫做面积和表，给定一个图像矩阵。对于那些不知道它是什么的人，可以参考维基百科的解释：

“积分图（也称为面积和表）是一种数据结构和算法，用于快速有效地生成网格子集中值的总和。”

换句话说，它用于在图像/矩阵中常数时间内求任何矩形区域的值之和。

我正在尝试在R中实现这个算法。然而，我的代码似乎运行时间过长。

以下是此链接中的伪代码。其中，in是输入的矩阵或图像，intImg是返回的结果。

对于i从0到w的循环，执行以下操作：
   sum←0
对于j从0到h的循环，执行以下操作：
      sum ← sum + in[i, j]
如果 i = 0，则
         intImg[i, j] ← sum
      否则，
         intImg[i, j] ← intImg[i − 1, j] + sum
      结束如果
   结束循环 
结束循环
这是我的实现：
w = ncol(im)
h = nrow(im)
intImg = c(NA)
length(intImg) = w*h
for(i in 1:w){ #x
  sum = 0;
  for(j in 1:h){ #y
    ind = ((j-1)*w)+ (i-1) + 1 #index
    sum = sum + im[ind]
    if(i == 1){
      intImg[ind] = sum
    }else{
      intImg[ind] = intImg[ind-1]+sum
    }
  }
}
intImg = matrix(intImg, h, w, byrow=T)
输入矩阵和输出矩阵的示例：

然而，在 480x640 矩阵上，这需要约4秒钟。在论文中，他们描述这些维度的运行时间为毫秒级别。

我在循环或索引方面做得不够高效吗？

我考虑用C++编写它并将其包装在R中，但我对C ++不是很熟悉。

谢谢

- Omar Wagih

---

一个指向 im 的链接可以使这个过程可重现。没有示例很难理解你在做什么，但你可能可以使用 colSums() 和 rowSums() 或其他一些矢量化函数。 - Chase

im 只是一个矩阵。我在问题中添加了输入和输出的示例。 - Omar Wagih

嗨，Joran，我已经尝试过了。没有太大的变化。有什么区别吗？ - Omar Wagih

@by0 - 它涉及到R如何为每次循环重新分配内存。 - Chase

变化不大，意思是仍然在4秒左右。确切地说是4.02秒。在这里，您可以在R控制台中复制粘贴以下内容，以获取一个非图像的示例矩阵：m = matrix(c(4,1,2,2,0,4,1,3,3,1,0,4,2,1,3,2), 4,4,byrow=T) - Omar Wagih

Answer 1

你可以尝试使用apply（如果你预分配内存，它不会比你的for循环更快）：

areaTable <- function(x) {
  return(apply(apply(x, 1, cumsum), 1, cumsum))
}

areaTable(m)
#      [,1] [,2] [,3] [,4]
# [1,]    4    5    7    9
# [2,]    4    9   12   17
# [3,]    7   13   16   25
# [4,]    9   16   22   33

- sgibb

---

虽然可能提升并不会很大，但这种方法会更快——因为你正在进行8个cumsum调用（加上形成循环代码和分配等操作），而OP的代码正在进行许多对+的调用来达到与你的cumsum调用相当的效果。由于它是解释性代码，所有这些函数调用都会累加。 - Gavin Simpson

@GavinSimpson 我其实很惊讶它的速度有多快，基于一个500x500的例子进行了快速计时。 - joran

是的，它快得多。我放弃了等待两个嵌套的for循环，双重应用不到一秒钟。 - Spacedman

@joran 我本应坚持我的初衷 - 在自我怀疑的时刻，我回去添加了警告之后才提交了评论。cumsum是非常高效的代码，所有这些+调用将随着输入维度的增加而增加。对OP的代码进行分析将会显示出来。 - Gavin Simpson

@GavinSimpson：你说得对，我没有考虑到解释性代码的这个事实。 - sgibb

@sgibb，请解释一下如何在C++中使其更快？ - Abc