给定一个二维数字“高度图"R中的矩阵，如何找到所有局部最大值? [英] Given a 2D numeric "height map" matrix in R, how can I find all local maxima?

问题描述

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

I have an R matrix with nonnegative numeric values. The matrix is effectively a 2D height map, and I want to find all the local maxima in this matrix. For "flat" peaks, where neighboring elements are equal to each other (and they are collectively a local maximum), I don't care what happens as long as I get at least one coordinate within each "flat" region.

是否有任何功能可以有效地执行此操作?显然，我可以手动编写循环来逐个遍历和测试每个元素，但是在 R 中这样做会很慢.我需要对大约一百万个矩阵执行此操作，每个矩阵平均大约有 884 个元素.

Are there any functions to do this efficiently? Obviously I could write the loop manually to go through and test every element individually, but doing that in R would be quite slow. I need to do this for around a million matrices, with an average of about 884 elements per matrix.

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

Ideally there would be a function that takes the matrix as input and returns a 2-column matrix with column 1 being row coordinates, column 2 being the column coordinates, and one row for each local maximum in the matrix.

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

Local maxima on the edges of the matrix are allowed. Areas outside of the matrix can be treated as having zero height.

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

Reproducible example matrix to use:

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

推荐答案

raster 包中的 focal() 函数专为此类计算而设计.以下代码返回所有局部极大值的坐标，包括边缘坐标和高原"部分的坐标.

The focal() function in the raster package is designed for calculations like this. The following code returns coordinates of all local maxima, including those on edges and those that are parts of "plateaus ".

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

这篇关于给定一个二维数字“高度图"R中的矩阵，如何找到所有局部最大值?的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！