Hoe bereken je efficiënt een Gauss-kernelmatrix in numpy?

Dit is mijn huidige manier. Is er een manier waarop ik matrixbewerking kan gebruiken om dit te doen? X zijn de gegevenspunten.

Antwoord 1, autoriteit 100%

Wilt u de Gauss-kernel gebruiken voor b.v. beeld gladstrijken? Zo ja, dan is er een functie gaussian_filter()in scipy:

Bijgewerkt antwoord

Dit zou moeten werken – hoewel het nog steeds niet 100% nauwkeurig is, probeert het rekening te houden met de waarschijnlijkheidsmassa binnen elke cel van het raster. Ik denk dat het gebruik van de kansdichtheid in het midden van elke cel iets minder nauwkeurig is, vooral voor kleine kernels. Zie https://homepages.inf.ed.ac.uk/rbf /HIPR2/gsmooth.htmvoor een voorbeeld.

import numpy as np
import scipy.stats as st
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()

Testen op het voorbeeld in Afbeelding 3 van de link:

gkern(5, 2.5)*273

geeft

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

Het originele (geaccepteerde) antwoord hieronder is fout
De vierkantswortel is niet nodig en de definitie van het interval is onjuist.

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

Antwoord 2, autoriteit 90%

Ik heb zelf het geaccepteerde antwoord gebruikt voor mijn beeldverwerking, maar ik vind het (en de andere antwoorden) te afhankelijk van andere modules. Daarom is hier mijn compacte oplossing:

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)

Bewerken: Arange gewijzigd in linspace om even zijdelingse lengtes te verwerken

Antwoord 3, autoriteit 83%

Ik probeer het antwoord van FuzzyDuck’s antwoordhier te verbeteren. Ik denk dat deze benadering korter en gemakkelijker te begrijpen is. Hier gebruik ik signal.scipy.gaussianom de 2D Gauss-kernel te krijgen.

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

Plotten met matplotlib.pyplot:

import matplotlib.pyplot as plt
plt.imshow(gkern(21), interpolation='none')

Antwoord 4, autoriteit 62%

U kunt eenvoudig een eenvoudige 2D dirac-functiefilteren, het resultaat is dan het filter functie die werd gebruikt:

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)

Antwoord 5, autoriteit 21%

Een 2D Gauss-kernelmatrix kan worden berekend met numpy broadcasting,

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

Voor kleine kernelgroottes zou dit redelijk snel moeten zijn.

Opmerking:dit maakt het wijzigen van de sigma-parameter eenvoudiger ten opzichte van het geaccepteerde antwoord.

Antwoord 6, autoriteit 17%

Ik heb geprobeerd alleen numpy te gebruiken. Hier is de code

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)

U kunt het resultaat visualiseren met:

plt.imshow(get_gauss_kernel(5,1))

Antwoord 7, autoriteit 17%

Als u een computer vision-ingenieur bent en u een heatmap nodig heeft voor een bepaald punt als Gauss-distributie (vooral voor keypoint-detectie op afbeelding)

def gaussian_heatmap(center = (2, 2), image_size = (10, 10), sig = 1):
    """
    It produces single gaussian at expected center
    :param center:  the mean position (X, Y) - where high value expected
    :param image_size: The total image size (width, height)
    :param sig: The sigma value
    :return:
    """
    x_axis = np.linspace(0, image_size[0]-1, image_size[0]) - center[0]
    y_axis = np.linspace(0, image_size[1]-1, image_size[1]) - center[1]
    xx, yy = np.meshgrid(x_axis, y_axis)
    kernel = np.exp(-0.5 * (np.square(xx) + np.square(yy)) / np.square(sig))
    return kernel

het gebruik en uitvoer

kernel = gaussian_heatmap(center = (2, 2), image_size = (10, 10), sig = 1)
plt.imshow(kernel)
print("max at :", np.unravel_index(kernel.argmax(), kernel.shape))
print("kernel shape", kernel.shape)

Max op: (2, 2)

Kernelvorm (10, 10)

kernel = gaussian_heatmap(center = (25, 40), image_size = (100, 50), sig = 5)
plt.imshow(kernel)
print("max at :", np.unravel_index(kernel.argmax(), kernel.shape))
print("kernel shape", kernel.shape)

MAX op: (40, 25)

Kernel-vorm (50, 100)

Antwoord 8, Autoriteit 14%

linalg.normneemt een axisparameter. Met een beetje experimenteren vond ik dat ik de norm kon berekenen voor alle combinaties van rijen met

np.linalg.norm(x[None,:,:]-x[:,None,:],axis=2)

Het breidt xuit in een 3D-reeks van alle verschillen en neemt de norm op de laatste dimensie.

Dus ik kan dit toepassen op uw code door de axisParameter toe te voegen aan uw 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)

produceert

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

Overeenkomst:

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

Antwoord 9, autoriteit 10%

Voortbouwend op het antwoord van Teddy Hartanto. Je kunt gewoon je eigen eendimensionale Gauss-functies berekenen en vervolgens np.outergebruiken om de tweedimensionale te berekenen. Zeer snelle en efficiënte manier.

Met onderstaande code kun je ook voor elke dimensie verschillende Sigma’s gebruiken

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

Antwoord 10, autoriteit 3%

Het door FuzzyDuck geaccepteerde antwoord aanpassen aan de resultaten van deze website: http:// dev.theomader.com/gaussian-kernel-calculator/Ik presenteer u nu deze definitie:

import numpy as np
import scipy.stats as st
def gkern(kernlen=21, sig=3):
    """Returns a 2D Gaussian kernel."""
    x = np.linspace(-(kernlen/2)/sig, (kernlen/2)/sig, kernlen+1)
    kern1d = np.diff(st.norm.cdf(x))
    kern2d = np.outer(kern1d, kern1d)
    return kern2d/kern2d.sum()
print(gkern(kernlen=5, sig=1))

uitvoer:

[[0.003765   0.015019   0.02379159 0.015019   0.003765  ]
 [0.015019   0.05991246 0.0949073  0.05991246 0.015019  ]
 [0.02379159 0.0949073  0.15034262 0.0949073  0.02379159]
 [0.015019   0.05991246 0.0949073  0.05991246 0.015019  ]
 [0.003765   0.015019   0.02379159 0.015019   0.003765  ]]

Antwoord 11, autoriteit 3%

Een goede manier om dat te doen is door de functie gaussian_filter te gebruiken om de kernel te herstellen.
Bijvoorbeeld:

indicatrice = np.zeros((5,5))
indicatrice[2,2] = 1
gaussian_kernel = gaussian_filter(indicatrice, sigma=1)
gaussian_kernel/=gaussian_kernel[2,2]

Dit geeft

array[[0.02144593, 0.08887207, 0.14644428, 0.08887207, 0.02144593],
       [0.08887207, 0.36828649, 0.60686612, 0.36828649, 0.08887207],
       [0.14644428, 0.60686612, 1.        , 0.60686612, 0.14644428],
       [0.08887207, 0.36828649, 0.60686612, 0.36828649, 0.08887207],
       [0.02144593, 0.08887207, 0.14644428, 0.08887207, 0.02144593]]

Antwoord 12

Omdat ik niet vond wat ik zocht, heb ik mijn eigen oneliner gecodeerd. U kunt het dienovereenkomstig aanpassen (volgens de afmetingen en de standaarddeviatie).

Hier is bijvoorbeeld de one-liner-functie voor een 3×5 patch.

from scipy import signal
def gaussian2D(patchHeight, patchWidth, stdHeight=1, stdWidth=1):
    gaussianWindow = signal.gaussian(patchHeight, stdHeight).reshape(-1, 1)@signal.gaussian(patchWidth, stdWidth).reshape(1, -1)
    return gaussianWindow
print(gaussian2D(3, 5))

Je krijgt een uitvoer als volgt:

[[0.082085   0.36787944 0.60653066 0.36787944 0.082085  ]
[0.13533528  0.60653066 1.         0.60653066 0.13533528]
[0.082085   0.36787944 0.60653066 0.36787944 0.082085  ]]

Je kunt meer lezen over scipy’s Gaussiaanse hier.

Antwoord 13

Weer een implementatie.

Dit is genormaliseerd zodat voor sigma > 1en voldoende grote win_size, de totale som van de kernelelementen is gelijk aan 1.

def gaussian_kernel(win_size, sigma):
    t = np.arange(win_size)
    x, y = np.meshgrid(t, t)
    o = (win_size - 1) / 2
    r = np.sqrt((x - o)**2 + (y - o)**2)
    scale = 1 / (sigma**2 * 2 * np.pi)
    return scale * np.exp(-0.5 * (r / sigma)**2)

Een 5×5 kernel genereren:

gaussian_kernel(win_size=5, sigma=1)