如何在numpy中有效地计算高斯核矩阵? [英] How to calculate a Gaussian kernel matrix efficiently in numpy?

问题描述

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是数据点.

This is my current way. Is there any way I can use matrix operation to do this? X is the data points.

推荐答案

是否要将高斯内核用于例如图像平滑?如果是这样，则有一个 gaussian_filter() :

Do you want to use the Gaussian kernel for e.g. image smoothing? If so, there's a function gaussian_filter() in scipy:

更新后的答案

这应该起作用-尽管它仍然不是100％准确，但它尝试考虑网格每个像元内的概率质量.我认为在每个像元的中点使用概率密度的准确性稍差，尤其是对于小内核.参见 https://homepages.inf.ed.ac.uk/rbf /HIPR2/gsmooth.htm 为例.

This should work - while it's still not 100% accurate, it attempts to account for the probability mass within each cell of the grid. I think that using the probability density at the midpoint of each cell is slightly less accurate, especially for small kernels. See https://homepages.inf.ed.ac.uk/rbf/HIPR2/gsmooth.htm for an example.

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

通过链接在图3的示例中对其进行测试:

Testing it on the example in Figure 3 from the link:

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

下面接受的原始(接受)答案是错误的 平方根是不必要的，并且间隔的定义不正确.

The original (accepted) answer below accepted is wrong The square root is unnecessary, and the definition of the interval is incorrect.

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

这篇关于如何在numpy中有效地计算高斯核矩阵?的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！