Numpy：生成二维高斯总和 pdf 作为数组

Question

我正在尝试生成一个 [600 x 600] numpy 数组，其中包含 10 个类高斯数组的总和（每个数组都有一个随机生成的中心）。

我试过使用高斯滤波器生成单个类似高斯的数组，然后对它们求和，但我确信有一种矢量化方法可以解决这个问题。即使 num_centers=10 也很慢，我可能需要对多达 20 个高斯求和。

这里有一个类似的问题，但似乎没有一个好的或决定性的答案，我不确定如何将它应用到我的问题中。 Sum of Gaussians into fast Numpy?

这是我试过的方法。

import numpy as np
from scipy.ndimage import gaussian_filter
import matplotlib.pyplot as plt


num_centers = 10               # number of Gaussians to sum
sigma = 100                    # std. dev. of each Gaussian
result = np.zeros((600, 600))


for _ in range(num_centers):

    # Pick a random coordinate within the array as the center
    center = np.random.uniform(result.shape).astype(int)

    # Make array with 1 at the center and 0 everywhere else
    temp = np.zeros_like(result)
    temp[center[0], center[1]] = 1

    # Apply filter
    gaussian = gaussian_filter(temp, sigma)

    # Add to result
    result += gaussian


# Result should look like a contour map with several hills
plt.imshow(result * 1000)        # scale up to see the coloring
plt.show()

Answer 1

我不确定您如何以并行方式处理随机高斯数组的创建，因为这在您的代码中花费的时间最多。（我用 timeit 来确定这一点）。这是意料之中的，因为 gaussian_filter 是一个计算密集型函数。

但是，通过在高斯数组上使用 np.sum()，我确实看到了轻微的性能提升。这是因为调用 np.sum() 一次比在循环中调用 += 更有效。

例子

import numpy as np
from scipy.ndimage import gaussian_filter
import matplotlib.pyplot as plt


num_centers = 10               # number of Gaussians to sum
sigma = 100                    # std. dev. of each Gaussian
holder = np.zeros((num_centers, 600, 600))


for _ in range(num_centers):

    # Pick a random coordinate within the array as the center
    center = np.random.uniform(result.shape).astype(int)

    # Make array with 1 at the center and 0 everywhere else
    temp = np.zeros((600, 600))
    temp[center[0], center[1]] = 1

    # Apply filter
    holder[_] = gaussian_filter(temp, sigma)

result = np.sum(holder, axis=0)

# Result should look like a contour map with several hills
plt.imshow(result * 1000)        # scale up to see the coloring
plt.show()

Answer 2

您可以消除循环，而是在每个中心创建一个值为 1 的数组，然后对这个数组应用 gaussian_filter 一次。所有步骤都可以向量化。

这是一个例子。我把 sigma 变小了，这样更容易区分中心，我把宽度增加到 800（没有特别的原因:)。

import numpy as np
from scipy.ndimage import gaussian_filter
import matplotlib.pyplot as plt


num_centers = 10
sigma = 25
size = (600, 800)

impulses = np.zeros(size)

# rows and cols are the row and column indices of the centers
# of the gaussian peaks.
np.random.seed(123456)
rows, cols = np.unravel_index(np.random.choice(impulses.size, replace=False,
                                               size=num_centers),
                              impulses.shape)
impulses[rows, cols] = 1
# or use this if you want duplicates to sum:
# np.add.at(impulses, (rows, cols), 1)

# Filter impulses to create the result.
result = gaussian_filter(impulses, sigma, mode='nearest')

plt.imshow(result)
plt.show()

情节如下：

您可以尝试 gaussian_filter 的 mode 参数，看看哪种模式最适合您。

Numpy：生成二维高斯总和 pdf 作为数组

Numpy: Generating a 2D Sum of Gaussians pdf as an array

python

arrays

numpy

gaussian

vectorization

例子