Numpy: Generando una sum 2D de gaussianos pdf como una matriz

Estoy tratando de generar una matriz numpy [600 x 600] que contenga la sum de 10 matrices de tipo gaussiano (cada una con un centro generado aleatoriamente).

He intentado usar un filtro gaussiano para generar los arreglos individuales de tipo gaussiano, luego los resumí, pero estoy seguro de que hay una forma vectorizada de abordar esto. Incluso con num_centers=10 es lento, y es posible que tenga que sumr hasta 20 gaussianos.

Aquí hay una pregunta similar, pero no parece tener una respuesta buena o concluyente y no estoy seguro de cómo aplicarlo a mi problema. Suma de gaussianos en Numpy rápido?

Esto es lo que he intentado.

 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() 

Puede eliminar el bucle y, en su lugar, crear una matriz con el valor 1 en cada centro y luego aplicar gaussian_filter una vez a esta matriz. Todos los pasos pueden ser vectorizados.

Aquí hay un ejemplo. Hice sigma más pequeño por lo que fue más fácil distinguir los centros, y aumenté el ancho a 800 (sin ninguna razón en particular :).

 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() 

Aquí está la ttwig:

trama

Puedes experimentar con el argumento de mode de gaussian_filter para ver qué modo funciona mejor para ti.

No estoy seguro de cómo manejar la creación de matrices gaussianas aleatorias de forma paralela, ya que eso es lo que más tiempo lleva en su código. (Utilicé timeit para determinar esto). Esto es de esperar, ya que gaussian_filter es una función de computación intensiva.

Sin embargo, vi un ligero aumento en el rendimiento al usar np.sum() en una variedad de gaussianos. Esto se debe a que llamar a np.sum() una vez es más eficiente que llamar a += desde dentro de un bucle.

Ejemplo

 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()