如何在MATLAB中不使用循环比较矩阵元素与其相邻元素？

Question

我有一个MATLAB矩阵。我想要检查每个元素的4个相邻值（左、右、上、下）。如果当前元素小于任何相邻值，则将其设置为零；否则，保留其值。使用循环可以轻松完成此操作，但由于我有数千个这些矩阵，所以这样做非常昂贵。 你可能会认识到这是边缘检测后的非极大值抑制。

- Shan

Answer 1

如果您有图像处理工具箱，可以使用形态学膨胀来查找局部极值并抑制所有其他元素。

array = magic(6); %# make some data

msk = [0 1 0;1 0 1;0 1 0]; %# make a 4-neighbour mask

%# dilation will replace the center pixel with the 
%# maximum of its neighbors
maxNeighbour = imdilate(array,msk);

%# set pix to zero if less than neighbors
array(array<maxNeighbour) = 0;

array =
    35     0     0    26     0     0
     0    32     0     0     0    25
    31     0     0     0    27     0
     0     0     0     0     0     0
    30     0    34     0     0    16
     0    36     0     0    18     0

编辑使用与@gnovice相同的数据，并修复代码

- Jonas

---

谢谢你的回答，这是一个很酷的方法，但不是我要找的...请考虑以下循环：for (y = 1; y < height - 1; ++y) { for (x = 1; x < width - 1; ++x) { if (response[y][x] < response[y][x+1] || result[y][x] < result[y][x-1] || result[y][x] < result[y+1][x] || result[y][x] < result[y-1][x]) tmp[y][x] = 0.0; } } - Shan

@Shan：哎呀，是的。我修复了这个解决方案。 - Jonas

+1 对于一个比我的解决方案快一个数量级的解决方案，但这可能是因为当我运行它时，NLFILTER似乎会短暂地显示进度条（而不必要的图形肯定会减慢速度）。 - gnovice

Answer 2

其中一种方法是使用图像处理工具箱中的NLFILTER函数，该函数将给定的函数应用于矩阵的每个M×N块：

>> A = magic(6)  %# A sample matrix

A =

    35     1     6    26    19    24
     3    32     7    21    23    25
    31     9     2    22    27    20
     8    28    33    17    10    15
    30     5    34    12    14    16
     4    36    29    13    18    11

>> B = nlfilter(A,[3 3],@(b) b(5)*all(b(5) >= b([2 4 6 8])))

B =

    35     0     0    26     0     0
     0    32     0     0     0    25
    31     0     0     0    27     0
     0     0     0     0     0     0
    30     0    34     0     0    16
     0    36     0     0    18     0

上面的代码定义了一个匿名函数，它使用线性索引来获取3x3子矩阵b(5)的中心元素，并将其与其4个相邻元素b([2 4 6 8])进行比较。中心元素中的值乘以函数ALL返回的逻辑结果，当中心元素大于所有最近的邻居时，该结果为1，否则为0。

- gnovice

---

正是我所寻找的... 希望它的性能比循环更好 :) 谢谢gnovie - Shan

+1 如果您第一次就翻译正确了，尽管我怀疑imdilate可能仍然更快。 - Jonas

@Jonas：你可能是对的，因为在最近几个MATLAB版本中，IMDILATE的性能有了一些显著的改进。 - gnovice

我想从矩阵中删除所有至少有一个0邻居的1，所以我将这段代码更改为以下内容：B = nlfilter(A,[3 3],@(b) b(5)*(0<sum(b([2 4 6 8])) && sum(b([2 4 6 8]))<4));，它可以正常工作。非常感谢！ - bordart

Answer 3

如果您无法访问图像处理工具箱，则另一种实现方法是构建代表每个点的顶部、右侧、底部和左侧的四个矩阵，并搜索所有四个矩阵中对应元素非负的元素（即该元素超过其所有邻居）。 以下是拆分的思路： 生成一些测试数据： >> sizeA = 3; A = randi(255, sizeA) A = 254 131 94 135 10 124 105 191 84 使用零元素填充边框： >> A2 = zeros(sizeA+2) * -Inf; A2(2:end-1,2:end-1) = A A2 = 0 0 0 0 0 0 254 131 94 0 0 135 10 124 0 0 105 191 84 0 0 0 0 0 0 构造四个一阶差分矩阵： >> leftDiff = A2(2:end-1,2:end-1) - A2(2:end-1,1:end-2) leftDiff = 254 -123 -37 135 -125 114 105 86 -107 >> topDiff = A2(2:end-1,2:end-1) - A2(1:end-2,2:end-1) topDiff = 254 131 94 -119 -121 30 -30 181 -40 >> rightDiff = A2(2:end-1,2:end-1) - A2(2:end-1,3:end) rightDiff = 123 37 94 125 -114 124 -86 107 84 >> bottomDiff = A2(2:end-1,2:end-1) - A2(3:end,2:end-1) bottomDiff = 119 121 -30 30 -181 40 105 191 84 查找超过所有邻居的元素： indexKeep = find(leftDiff >= 0 & topDiff >= 0 & rightDiff >= 0 & bottomDiff >= 0) 创建结果矩阵：

>> B = zeros(sizeA);
B(indexKeep) = A(indexKeep)

B =

   254     0     0
     0     0   124
     0   191     0

将这些内容包装成一个函数并在1000个随机的100x100矩阵上进行测试后，该算法似乎非常快：

>> tic;
for ii = 1:1000
A = randi(255, 100);
B = test(A);
end; toc
经过时间为0.861121秒。

- b3.

---

非常好。当你可以使用高级工具箱时，很容易忘记这样简单的实现方式。 ;) 一个建议...你可以使用CIRCSHIFT来简化代码，像这样：index = A2>=circshift(A2,[0 1]) & A2>=circshift(A2,[0 -1]) & A2>=circshift(A2,[1 0]) & A2>=circshift(A2,[-1 0]); A2 = A2.*index; A2 = A2(2:end-1,2:end-1); 另外，我认为你想使用>=而不是>。 - gnovice

哦，另一个使算法更通用的建议是使用“-Inf”来填充矩阵，而不是0，以考虑矩阵中可能存在负值的情况（即您不希望填充出现为局部最大值）。 - gnovice

@gnovice：感谢您的建议。已经进行了编辑以包含它们。没有实现CIRCSHIFT只是为了保持更加繁琐的方法，但我同意使用CIRCSHIFT会是一个更好的解决方案。 - b3.