def GaussianMatrix(X,sigma):
row,col=X.shape
GassMatrix=np.zeros(shape=(row,row))
X=np.asarray(X)
i=0
for v_i in X:
j=0
for v_j in X:
GassMatrix[i,j]=Gaussian(v_i.T,v_j.T,sigma)
j+=1
i+=1
return GassMatrix
def Gaussian(x,z,sigma):
return np.exp((-(np.linalg.norm(x-z)**2))/(2*sigma**2))
Это мой нынешний путь. Есть ли способ использовать матричную операцию для этого? X - это точки данных.
Вы хотите использовать ядро Гаусса, например, для сглаживания изображений? Если это так, в scipy есть функция gaussian_filter():
Обновленный ответ
Это должно сработать - хотя оно все еще не является точным на 100%, оно пытается учесть массу вероятности в каждой ячейке сетки. Я думаю, что использование плотности вероятности в средней точке каждой ячейки немного менее точно, особенно для небольших ядер. См. Пример https://homepages.inf.ed.ac.uk/rbf/HIPR2/gsmooth.htm.
def gkern(kernlen=21, nsig=3):
"""Returns a 2D Gaussian kernel."""
x = np.linspace(-nsig, nsig, kernlen+1)
kern1d = np.diff(st.norm.cdf(x))
kern2d = np.outer(kern1d, kern1d)
return kern2d/kern2d.sum()
Тестируем его на примере на рисунке 3 по ссылке:
gkern(5, 2.5)*273
дает
array([[ 1.0278445 , 4.10018648, 6.49510362, 4.10018648, 1.0278445 ],
[ 4.10018648, 16.35610171, 25.90969361, 16.35610171, 4.10018648],
[ 6.49510362, 25.90969361, 41.0435344 , 25.90969361, 6.49510362],
[ 4.10018648, 16.35610171, 25.90969361, 16.35610171, 4.10018648],
[ 1.0278445 , 4.10018648, 6.49510362, 4.10018648, 1.0278445 ]])
Исходный (принятый) ответ ниже принятый является неправильным Квадратный корень не нужен, и определение интервала является неправильным.
import numpy as np
import scipy.stats as st
def gkern(kernlen=21, nsig=3):
"""Returns a 2D Gaussian kernel array."""
interval = (2*nsig+1.)/(kernlen)
x = np.linspace(-nsig-interval/2., nsig+interval/2., kernlen+1)
kern1d = np.diff(st.norm.cdf(x))
kernel_raw = np.sqrt(np.outer(kern1d, kern1d))
kernel = kernel_raw/kernel_raw.sum()
return kernel
Вы можете просто gaussian-фильтровать простую функцию 2D dirac, результатом будет функция фильтра, которая используется:
import numpy as np
import scipy.ndimage.filters as fi
def gkern2(kernlen=21, nsig=3):
"""Returns a 2D Gaussian kernel array."""
# create nxn zeros
inp = np.zeros((kernlen, kernlen))
# set element at the middle to one, a dirac delta
inp[kernlen//2, kernlen//2] = 1
# gaussian-smooth the dirac, resulting in a gaussian filter mask
return fi.gaussian_filter(inp, nsig)
Я сам использовал принятый ответ для обработки изображений, но нахожу его (и другие ответы) слишком зависимым от других модулей. Поэтому вот мое компактное решение:
import numpy as np
def gkern(l=5, sig=1.):
"""\
creates gaussian kernel with side length l and a sigma of sig
"""
ax = np.linspace(-(l - 1) / 2., (l - 1) / 2., l)
xx, yy = np.meshgrid(ax, ax)
kernel = np.exp(-0.5 * (np.square(xx) + np.square(yy)) / np.square(sig))
return kernel / np.sum(kernel)
Редактировать: изменено значение arange на linspace для обработки равных сторон
Я пытаюсь улучшить ответ FuzzyDuck здесь. Я думаю, что этот подход короче и проще для понимания. Здесь я использую signal.scipy.gaussian для получения двумерного гауссовского ядра.
import numpy as np
from scipy import signal
def gkern(kernlen=21, std=3):
"""Returns a 2D Gaussian kernel array."""
gkern1d = signal.gaussian(kernlen, std=std).reshape(kernlen, 1)
gkern2d = np.outer(gkern1d, gkern1d)
return gkern2d
matplotlib.pyplot графика с использованием matplotlib.pyplot:
import matplotlib.pyplot as plt
plt.imshow(gkern(21), interpolation='none')
Двумерная гауссовская матрица ядра может быть вычислена с помощью широковещательной передачи numpy,
def gaussian_kernel(size=21, sigma=3):
"""Returns a 2D Gaussian kernel.
Parameters
----------
size : float, the kernel size (will be square)
sigma : float, the sigma Gaussian parameter
Returns
-------
out : array, shape = (size, size)
an array with the centered gaussian kernel
"""
x = np.linspace(- (size // 2), size // 2)
x /= np.sqrt(2)*sigma
x2 = x**2
kernel = np.exp(- x2[:, None] - x2[None, :])
return kernel / kernel.sum()
Для небольших размеров ядра это должно быть достаточно быстро.
Примечание: это упрощает изменение параметра сигмы относительно принятого ответа.
linalg.norm принимает параметр axis. При небольшом эксперименте я обнаружил, что могу рассчитать норму для всех комбинаций строк с
np.linalg.norm(x[None,:,:]-x[:,None,:],axis=2)
Он расширяет x в 3d-массив всех различий и принимает норму в последнем измерении.
Поэтому я могу применить это к вашему коду, добавив параметр axis к вашему Gaussian:
def Gaussian(x,z,sigma,axis=None):
return np.exp((-(np.linalg.norm(x-z, axis=axis)**2))/(2*sigma**2))
x=np.arange(12).reshape(3,4)
GaussianMatrix(x,1)
производит
array([[ 1.00000000e+00, 1.26641655e-14, 2.57220937e-56],
[ 1.26641655e-14, 1.00000000e+00, 1.26641655e-14],
[ 2.57220937e-56, 1.26641655e-14, 1.00000000e+00]])
Matching:
Gaussian(x[None,:,:],x[:,None,:],1,axis=2)
array([[ 1.00000000e+00, 1.26641655e-14, 2.57220937e-56],
[ 1.26641655e-14, 1.00000000e+00, 1.26641655e-14],
[ 2.57220937e-56, 1.26641655e-14, 1.00000000e+00]])
На основе ответа Тедди Хартанто. Вы можете просто вычислить свои собственные одномерные гауссовские функции, а затем использовать np.outer для вычисления двухмерного. Очень быстрый и эффективный способ.
С помощью приведенного ниже кода вы также можете использовать разные Sigmas для каждого измерения
import numpy as np
def generate_gaussian_mask(shape, sigma, sigma_y=None):
if sigma_y==None:
sigma_y=sigma
rows, cols = shape
def get_gaussian_fct(size, sigma):
fct_gaus_x = np.linspace(0,size,size)
fct_gaus_x = fct_gaus_x-size/2
fct_gaus_x = fct_gaus_x**2
fct_gaus_x = fct_gaus_x/(2*sigma**2)
fct_gaus_x = np.exp(-fct_gaus_x)
return fct_gaus_x
mask = np.outer(get_gaussian_fct(rows,sigma), get_gaussian_fct(cols,sigma_y))
return mask
Я пытался использовать только NumPy. Вот код
def get_gauss_kernel(size=3,sigma=1):
center=(int)(size/2)
kernel=np.zeros((size,size))
for i in range(size):
for j in range(size):
diff=np.sqrt((i-center)**2+(j-center)**2)
kernel[i,j]=np.exp(-(diff**2)/(2*sigma**2))
return kernel/np.sum(kernel)
Вы можете визуализировать результат, используя:
plt.imshow(get_gauss_kernel(5,1))
