在MATLAB中对三维矩阵的切片进行向量化访问

Question

我有一个大约这些尺寸的三维矩阵：

A = rand(20, 1000, 20);

第一和第三维始终具有相同的长度。我想将主对角线切片中的元素归零。这样可以实现我的意图。

for ii = 1:size(A, 1)
    A(ii, :, ii) = 0;
end

有没有矢量化或其他更快的方法来完成这个任务？这段代码大约会运行100,000次，每次的大小都大致相同，但不完全相同。

- buzjwa

---

同样的矩阵维度，它能运行10万次吗？ - Daniel

不完全是这样。它在每个矩阵上运行了约100次，而且有大约1000个这样大小略有不同的矩阵。如果它们总是相同的大小，我会使用单个逻辑索引矩阵。 - buzjwa

Answer 1

在对一个1000 20 20的矩阵进行操作时，您可以使用逻辑索引以针对多个尾维度进行操作，同时使用下标索引单独针对所有前一维度进行操作。这样做可以轻松地完成操作。但是，为了将此应用于您的矩阵，需要使用permute，这可能会导致速度较慢：

n=size(A,3)
A=permute(A,[2,1,3]);
A(:,diag(true(n,1)))=0;
A=permute(A,[2,1,3]);

如果您的代码可以永久交换A的维度并避免使用permute，那么这将导致最快的解决方案。

或者，您可以使用repmat将索引扩展到A的维度。

ix=repmat(reshape(diag(true(n,1)),n,1,n),[1,size(A,2),1])
A(ix)=0

对于相同大小的矩阵，您可以保留ix。由于我现在无法访问MATLAB，我不知道哪种解决方案更快。

- Daniel

---

我现在测试了这些方法。它们非常慢。我猜这是其中一个向量化不占优势的地方。 - buzjwa

Answer 2

你可以使用bsxfun来构建需要清零的元素的线性索引：

ind = bsxfun(@plus, (0:size(A,2)-1).'*size(A,1), 1:size(A,1)*size(A,2)+1:numel(A) );
A(ind) = 0;

- Luis Mendo

---

这比丹尼尔的建议更快，但仍然比for循环方法慢超过2倍。有趣的是，即使所有A矩阵的大小都是恒定的，并且只调用一次bsxfun，情况仍然如此。执行100,000次A(ind) = 0;仍然比执行100,000次for循环慢。我猜这是由于使用A(ii, :, ii)进行更快的内存访问，而不是使用A(ind)访问“随机”索引。如果这是真的，那么这是反对向量化此代码的情况。 - buzjwa

但是既然我要求一个矢量化的解决方案，我将接受这个答案作为最佳矢量化的解决方案。 - buzjwa