给定R中的2D数字"高度图"矩阵,我如何找到所有局部最大值？

Question

我有一个带有非负数值的R矩阵.矩阵实际上是2D高度图,我想找到该矩阵中的所有局部最大值.对于"平坦"峰值,其中相邻元素彼此相等(并且它们共同是局部最大值),只要在每个"平坦"区域内得到至少一个坐标,我不关心会发生什么.

有效地执行此操作是否有任何功能？显然,我可以手动编写循环来单独测试每个元素,但在R中执行此操作会非常慢.我需要为大约一百万个矩阵做这个,每个矩阵平均大约884个元素.

理想情况下,会有一个函数将矩阵作为输入并返回一个2列矩阵,其中第1列是行坐标,第2列是列坐标,矩阵中每个局部最大值有一行.

允许矩阵边缘的局部最大值.矩阵之外的区域可以被视为具有零高度.

要使用的可重现的示例矩阵:

set.seed(5)
msize <- 20 # Change this to whatever you like
x <- matrix(data=abs(rnorm(msize*2)), nrow=msize, ncol=msize)

Answer 1

包中的focal()函数raster是为这样的计算而设计的.以下代码返回所有局部最大值的坐标,包括边缘上的坐标和属于"高原"部分的坐标.

library(raster)

## Construct an example matrix
set.seed(444)
msize <- 10
x <- matrix(sample(seq_len(msize), msize^2, replace=TRUE), ncol=msize)

## Convert it to a raster object
r <- raster(x)
extent(r) <- extent(c(0, msize, 0, msize) + 0.5)

## Find the maximum value within the 9-cell neighborhood of each cell
f <- function(X) max(X, na.rm=TRUE)
ww <- matrix(1, nrow=3, ncol=3) ## Weight matrix for cells in moving window
localmax <- focal(r, fun=f, w=ww, pad=TRUE, padValue=NA)

## Does each cell have the maximum value in its neighborhood?
r2 <- r==localmax

## Get x-y coordinates of those cells that are local maxima
maxXY <- xyFromCell(r2, Which(r2==1, cells=TRUE))
head(maxXY)
#       x  y
# [1,]  8 10
# [2,] 10 10
# [3,]  3  9
# [4,]  4  9
# [5,]  1  8
# [6,]  6  8

# Visually inspect the data and the calculated local maxima
plot(r)   ## Plot of heights
windows() ## Open a second plotting device
plot(r2)  ## Plot showing local maxima