如何使用 R 生成离散的二维高斯平滑核答案

【问题标题】：How to generate a discrete 2D Gaussian smoothing kernel using R如何使用 R 生成离散的二维高斯平滑核
【发布时间】：2015-02-15 02:16:29
【问题描述】：

我想知道 R 是否具有此功能，或者如何编写高斯平滑内核的 2D 离散版本。目标是使用它来平滑二维值矩阵。

【问题讨论】：

我很困惑你在这种情况下所说的离散是什么意思。从表面上看，将离散数据平滑为连续数据的内核没有任何意义（尽管对于一些简单的数据，即目的而言，这可能是合理的）。
使用矩阵（内核）来平滑更大的矩阵（数据）。我的数据以离散间隔采样，不代表连续函数，因此存储在矩阵中。
另一种思考方式是离散卷积。

标签： r smoothing

【解决方案1】：

我写了一个函数来解决我的问题。这是我称为kernelsmooth 的函数。它需要三个参数，一个要平滑的数据矩阵、一个内核和一个布尔标志 norm，如果为真，则通过其值的平均值对内核进行归一化（防止矩阵数据膨胀）。

kernelsmooth <- function(x, kern, norm=TRUE) {
  # how many rows/cols of zeroes are used to pad.
  width <- dim(kern)[1]
  pad <- floor(width / 2)

  # record the width and height the input data matrix
  x_w <- ncol(x)
  x_h <- nrow(x)

  # Are we normalizing the kernel?
  if (norm == TRUE) {
    k <- kern / sum(abs(kern))
  } else {
    k <- kern
  }

  # pad all around the matrix an equal width of zeros
  x_pad <- t(padzeros(data=x, nzeros=pad, side="both"))
  x_pad <- t(padzeros(data=x_pad, nzeros=pad, side="both"))

  # Pre-allocate the final (smoothed) data matrix
  s <- matrix(0, nrow = x_h, ncol = x_w)

  # Pre-allocate a temporary matrix for the iterative calculations
  temp <- matrix(0, width, width)

  # Loop through the data to apply the kernel.
  for (col in 1:x_w ) {
    for (row in 1:x_h ) {
      temp <- x_pad[row:(row + width - 1), col:(col + width - 1)]
      s[row,col] <-  sum(k * temp)
    }
  }

  # return the smoothed data
  return(s)
}

一些例子：

不要平滑矩阵（应用类似于单位矩阵的单位核）。

kernelsmooth(x = matrix(2,5,5), kern = matrix(c(0,0,0, 0,1,0, 0,0,0) , 3, 3))

输出：

     [,1] [,2] [,3] [,4] [,5] [,6]
[1,]    2    2    2    2    2    2
[2,]    2    2    2    2    2    2
[3,]    2    2    2    2    2    2
[4,]    2    2    2    2    2    2
[5,]    2    2    2    2    2    2
[6,]    2    2    2    2    2    2

平均平滑（每个元素通过其自身及其相邻元素的平均值进行平滑）：

kernelsmooth(x = matrix(4,5,5), kern = matrix(1,3,3))

输出：

         [,1]     [,2]     [,3]     [,4]     [,5]
[1,] 1.777778 2.666667 2.666667 2.666667 1.777778
[2,] 2.666667 4.000000 4.000000 4.000000 2.666667
[3,] 2.666667 4.000000 4.000000 4.000000 2.666667
[4,] 2.666667 4.000000 4.000000 4.000000 2.666667
[5,] 1.777778 2.666667 2.666667 2.666667 1.777778

高斯平滑（各向同性平滑）：

高斯

输出：

         [,1]     [,2]     [,3]     [,4]     [,5]     [,6]     [,7]     [,8]
[1,] 2.491349 3.321799 3.529412 3.529412 3.529412 3.529412 3.321799 2.491349
[2,] 3.321799 4.429066 4.705882 4.705882 4.705882 4.705882 4.429066 3.321799
[3,] 3.529412 4.705882 5.000000 5.000000 5.000000 5.000000 4.705882 3.529412
[4,] 3.529412 4.705882 5.000000 5.000000 5.000000 5.000000 4.705882 3.529412
[5,] 3.529412 4.705882 5.000000 5.000000 5.000000 5.000000 4.705882 3.529412
[6,] 3.529412 4.705882 5.000000 5.000000 5.000000 5.000000 4.705882 3.529412
[7,] 3.321799 4.429066 4.705882 4.705882 4.705882 4.705882 4.429066 3.321799
[8,] 2.491349 3.321799 3.529412 3.529412 3.529412 3.529412 3.321799 2.491349

【讨论】：

你没有包含函数padzeros的代码

【解决方案2】：

矩阵边缘的急剧下降，例如从 3.529412 到 0，是个坏消息。您需要将矩阵放大到边缘值尽可能接近零的大小。您可以通过计算每个矩阵的 2D 傅里叶变换并将其与在较大版本上完成的相同变换进行比较来检查此语句的有效性。应用于您的数据时，您还会看到更好的结果。

【讨论】：