将二进制矩阵快速向量化为最后一个非零索引的向量

Question

假设在MATLAB中，我有一个矩阵A，其元素为0或1。

如何以更快、向量化的方式获取每列最后一个非零元素的索引向量？

我可以执行以下操作：

[B，I] = max(cumsum(A));

然后使用I，但是否有更快的方法？（我假设即使是对0和1求和，cumsum也会花费一些时间）。

编辑：我想我向量化的程度甚至比需要快 - Mr. Fooz的循环很棒，但MATLAB中的每个循环似乎都会让我在调试时间上付出很大的代价，即使它很快。

- Joe Soul-bringer

Answer 1

你应该关注的是速度，而不是完全矢量化。最近版本的Matlab在处理循环效率方面要聪明得多。如果有一种紧凑的矢量化表达方式，通常会更快，但是循环不应该像过去那样被害怕。

clc

A = rand(5000)>0.5;
A(1,find(sum(A,1)==0)) = 1; % make sure there is at least one match

% Slow because it is doing too much work
tic;[B,I1]=max(cumsum(A));toc

% Fast because FIND is fast and it runs the inner loop
tic;
I3=zeros(1,5000);
for i=1:5000
  I3(i) = find(A(:,i),1,'last');
end
toc;
assert(all(I1==I3));

% Even faster because the JIT in Matlab is smart enough now
tic;
I2=zeros(1,5000);
for i=1:5000
  I2(i) = 0;
  for j=5000:-1:1
    if A(j,i)
      I2(i) = j;
      break;
    end
  end
end
toc;
assert(all(I1==I2));

在 R2008a、Windows、x64 上，cumsum 版本需要 0.9 秒。而循环和查找版本只需 0.02 秒，双重循环版本仅需 0.001 秒。

编辑：哪个更快取决于实际数据。当将 0.5 更改为 0.999 时，双重循环需要 0.05 秒（因为平均命中时间更长）。cumsum 和循环&查找实现具有更一致的速度。

编辑 2：gnovice 的 flipud 解决方案很巧妙。不幸的是，在我的测试机上它需要 0.1 秒，所以它比 cumsum 更快，但比循环版本慢。

- Mr Fooz

---

哇，你的循环有什么特别之处使它更快，还是任何相似的循环都是最快的方式？ - Joe Soul-bringer

当我开始编写这些示例时，我预计双重循环最慢，而循环和查找最快。当内部循环必须运行到完成时，它会有点慢（请参见编辑2）。现在，Matlab对每个函数进行即时编译。这使得循环速度更快（但对于喜欢使用EVAL的人来说，会产生一些意想不到的后果）。一般来说，如果可以在不做额外工作的情况下完成向量化，则仍然最好使用向量化（flipud和cumsum解决方案是向量化的，但会做额外的工作）。 - Mr Fooz

在多核机器上，您甚至可以将PARFOR（并行for循环）纳入外部循环中，因为内部循环在每个列上独立运行。 - gnovice

我尝试了gnovice建议的parfor。第一次运行需要近1秒钟（可能是因为它需要加载parfor逻辑和/或创建工作线程池），但几次后，它的执行时间降到了0.04秒。所以至少在这个简单的循环（只有5000次迭代，在一台普通的双核处理器上）中，串行版本更快。 - Mr Fooz

嗯...我猜循环中所做的工作是如此微小，以至于并行化产生的额外开销使其相形见绌。 - gnovice

Answer 2

正如Fooz先生所示，使用较新版本的MATLAB可以使for循环变得相当快。但是，如果您真的想要紧凑的向量化代码，我建议尝试这个：

[B,I] = max(flipud(A));
I = size(A,1)-I+1;

这个方法比基于CUMSUM的答案要快，但仍然不如Fooz先生的循环选项快。

还有两件事需要考虑：

• 如果一列中没有1，你想得到什么结果？使用我上面给出的选项，我相信在这种情况下，您将得到一个索引size(A,1)（即A中的行数）。对于您的选项，在这种情况下，我相信您将得到一个1，而来自Fooz先生的嵌套for循环选项将给您一个0。

• 这些不同选项的相对速度可能会根据A的大小和您期望它具有的非零数的数量而有所不同。

- gnovice

---

聪明的想法。不幸的是，它比循环和查找慢了大约5倍。 - Mr Fooz

这差不多是我预期的结果：比CUMSUM快但仍然比循环慢...尽管这仍然取决于A的大小和填充分数（OP并没有真正定义）。 - gnovice