是否有R函数来获得矩阵相邻值的“加权”和(对于变化的半径)?

问题描述 投票:0回答:2

我有一个矩阵M(比如r行和c列),我想根据它的邻居得到每个矩阵元素的“加权”和,并创建一个新的矩阵M2。邻居这个词可以在1的半径范围内(通常称为元胞自动机理论中的摩尔邻域),或者半径可以是1,例如2,3等。

对于矩阵M中的特定单元格,请说出中间的某个位置。让我们说位置(i,j);然后第(i,j)个小区有“八个”邻居给出,

(i-1, j-1), (i-1, j), (i-1, j+1), (i, j-1), (i, j+1), (i+1, j-1), (i+1, j), (i+1, j+1).

我想创建一个矩阵M2,它计算第(i,j)个单元加上其八个邻居的“加权”和。基于细胞之间的欧几里德距离进行加权。所以,例如,

exp(-sqrt(2))*M[i-1,j-1] + exp(-1)*M[i-1,j] + exp(-sqrt(2))*M[i-1,j+1] + exp(-1)*M[i,j-1] + M[i,j] + exp(-1)*M[i,j+1] + exp(-sqrt(2))*M[i+1,j-1] + exp(-1)*M[i+1,j] + exp(-sqrt(2))*M[i+1,j+1]

对于所有细胞重复相同的想法(沿着边界的细胞需要被特别处理,因为它们不一定具有八个相邻细胞)。上面的想法适用于半径1,但我想要开发的代码需要对任何半径都是通用的。

r <- 4
c <- 4

n <- r*c

(M <- matrix(1:n, r, c))

addresses <- expand.grid(x = 1:r, y = 1:c)

#Got this code in the same forum

z <- rbind(c(-1,0,1,-1,1,-1,0,1),c(-1,-1,-1,0,0,1,1,1))

get.neighbors <- function(rw) {
  # Convert to absolute addresses 
  z2 <- t(z + unlist(rw))
  # Choose those with indices within mat 
  b.good <- rowSums(z2 > 0)==2  &  z2[,1] <= nrow(M)  &  z2[,2] <= ncol(M)
  M[z2[b.good,]]
}

apply(addresses,1 , get.neighbors) # Returns a list with neighbors

M

本质上,半径= 1的M2必须是当前单元加上邻居的“加权”和。当前当前单元格的权重始终为1。

M = [ 1  5   9  13
      2  6  10  14
      3  7  11  15
      4  8  12  16]

M2 = [ 5.033  13.803 .... ....
       ....   ....   .... ....
       ....   ....   .... ....
       ....   ....   .... ....]

如何在R中获取矩阵M2?如果半径超过1怎么样?我想在两个for循环中发生加权,所以我可以在代码中进一步使用[i,j]单元格的计算加权和来关闭两个for循环。

r for-loop matrix nearest-neighbor adjacency-matrix
2个回答
0
投票

我认为以下是你想要的加权和。我会以类似的方式找到邻居@ r2evans。

wtd_nbrs_sum <- function(input_matrix,
                         radius,
                         weight_matrix)
{
  temp_1 <- matrix(data = 0,
                   nrow = nrow(x = input_matrix),
                   ncol = radius)
  temp_2 <- matrix(data = 0,
                   nrow = radius,
                   ncol = ((2 * radius) + ncol(x = input_matrix)))
  input_matrix_modified <- rbind(temp_2,
                                 cbind(temp_1, input_matrix, temp_1),
                                 temp_2)
  output_matrix <- matrix(nrow = nrow(x = input_matrix),
                          ncol = ncol(x = input_matrix))
  for(i in seq_len(length.out = nrow(x = input_matrix)))
  {
    for(j in seq_len(length.out = nrow(x = input_matrix)))
    {
      row_min <- (radius + (i - radius))
      row_max <- (radius + (i + radius))
      column_min <- (radius + (j - radius))
      column_max <- (radius + (j + radius))
      neighbours <- input_matrix_modified[(row_min:row_max), (column_min:column_max)]
      weighted_sum <- sum(neighbours * weight_matrix)
      output_matrix[i, j] <- weighted_sum
    }
  }
  return(output_matrix)
}

r <- 4
c <- 4
n <- r*c
M <- matrix(data = 1:n,
            nrow = r,
            ncol = c)
R <- 1
wts <- matrix(data = c(exp(x = -sqrt(x = 2)), exp(x = -1), exp(x = -sqrt(x = 2)), exp(x = -1), 1, exp(x = -1), exp(x = -sqrt(x = 2)), exp(x = -1), exp(x = -sqrt(x = 2))),
              nrow = 3,
              ncol = 3)

wtd_nbrs_sum(input_matrix = M,
             radius = R,
             weight_matrix = wts)
#>           [,1]     [,2]     [,3]     [,4]
#> [1,]  5.033856 13.80347 24.16296 23.89239
#> [2,]  8.596195 20.66391 34.43985 32.84175
#> [3,] 11.186067 24.10789 37.88383 35.43163
#> [4,]  9.748491 19.86486 30.22435 28.60703

reprex package创建于2019-03-24(v0.2.1)

希望这可以帮助。


0
投票

编辑包括加权和。

可能有一个非常巧妙的技巧,但最直接(和可维护)的方式可能是一个简单的两个for循环实现。

M1 <- matrix(1:16, nr=4)
M1
#      [,1] [,2] [,3] [,4]
# [1,]    1    5    9   13
# [2,]    2    6   10   14
# [3,]    3    7   11   15
# [4,]    4    8   12   16

代码:

get_neighbors <- function(M, radius = 1) {
  M2 <- M
  M2[] <- 0
  nr <- nrow(M)
  nc <- ncol(M)
  eg <- expand.grid((-radius):radius, (-radius):radius)
  eg$wt <- exp(-sqrt(abs(eg[,1]) + abs(eg[,2])))
  for (R in seq_len(nr)) {
    for (C in seq_len(nc)) {
      ind <- cbind(R + eg[,1], C + eg[,2], eg[,3])
      ind <- ind[ 0 < ind[,1] & ind[,1] <= nr &
                    0 < ind[,2] & ind[,2] <= nc,, drop = FALSE ]
      M2[R,C] <- sum(M[ind[,1:2, drop=FALSE]] * ind[,3])
    }
  }
  M2
}

get_neighbors(M1, 1)
#           [,1]     [,2]     [,3]     [,4]
# [1,]  5.033856 13.80347 24.16296 23.89239
# [2,]  8.596195 20.66391 34.43985 32.84175
# [3,] 11.186067 24.10789 37.88383 35.43163
# [4,]  9.748491 19.86486 30.22435 28.60703

同样的事情,半径为2:

get_neighbors(M1, 2)
#          [,1]     [,2]     [,3]     [,4]
# [1,] 12.44761 25.64963 31.73247 32.70974
# [2,] 18.57765 35.96237 43.33862 43.51911
# [3,] 20.09836 37.80643 45.18268 45.03982
# [4,] 17.51314 31.88500 37.96784 37.77527

一个简单的测试,如果使用半径0,那么M1和M2应该是相同的(它们是)。

注意:这通常在基数R中表现不错,没有花哨使用apply或其表兄弟。由于这是一个非常简单的启发式算法,因此可以很容易地使用Rcpp实现更快。

© www.soinside.com 2019 - 2024. All rights reserved.