Существует ли функция R для получения «взвешенной» суммы соседних значений матрицы (для различного радиуса)? - PullRequest
0 голосов
/ 24 марта 2019

У меня есть матрица M (скажем, r строк и c столбцов), и я хотел бы получить «взвешенную» сумму для каждого элемента матрицы на основе его соседей и создать новую матрицу M2. Слово «сосед» может быть в радиусе 1 (что часто называют соседством Мура в теории клеточных автоматов), или радиус может отличаться от 1, скажем, 2, 3 и т. Д.

Для конкретной ячейки в матрице M, скажем, где-то посередине. Скажем, положение (я, 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

По сути, M2 для радиуса = 1 должен быть «взвешенной» суммой текущей ячейки плюс соседей. Текущий текущий элемент всегда получает вес 1.

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

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

Как мне получить матрицу M2 в R? Как насчет, если радиус больше 1? Мне бы хотелось, чтобы взвешивание происходило внутри двух циклов for, чтобы я мог использовать вычисленную взвешенную сумму ячейки [i, j] далее в коде, закрывающем два цикла for.

Ответы [ 2 ]

0 голосов
/ 24 марта 2019

Я думаю, что следующее делает взвешенную сумму, которую вы хотите. Я буду искать соседей таким же образом, как и @ 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

Создано в 2019-03-24 с помощью пакета представительства (v0.2.1)

Надеюсь, это поможет.

0 голосов
/ 24 марта 2019

Отредактировано для включения взвешенных сумм.

Для этого может быть очень интересный трюк, но самый простой (и обслуживаемый) способ - это, вероятно, простой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, чтобы она была значительно быстрее.

...