在Python中實現離散高斯內核？

Question

我期待實施Lindeberg在其關於尺度空間理論的工作中定義的discrete Gaussian kernel。它被定義爲T（n，t）= exp（-t）* I_n（t），其中I_n是modified Bessel function of the first kind。

我想在Python中使用Numpy和Scipy實現這一點，但我遇到了一些麻煩。

def discrete_gaussian_kernel(t, n): 
    return math.exp(-t) * scipy.special.iv(n, t) 

我嘗試用繪圖：

import math 
import numpy as np 
import scipy 
from matplotlib import pyplot as plt 

def kernel(t, n): 
    return math.exp(-t) * scipy.special.iv(n, t) 
ns = np.linspace(-5, 5, 1000) 
y0 = discrete_gaussian_kernel(0.5, ns) 
y1 = discrete_gaussian_kernel(1, ns) 
y2 = discrete_gaussian_kernel(2, ns) 
y3 = discrete_gaussian_kernel(4, ns) 
plt.plot(ns, y0, ns, y1, ns, y2, ns, y3) 
plt.xlim([-4, 4]) 
plt.ylim([0, 0.7]) 

輸出看起來像：

從維基百科的文章，它應該看起來像：

我假設我犯了一些非常微不足道的錯誤。 ：/ 有什麼想法嗎？ 謝謝！

編輯： 我寫的相當於scipy.special.ive(n, t)。我很確定它應該是第一種改進的貝塞爾函數，而不是第二種，但是可以有人確認嗎？

Answer 1

如果你想獲得維基百科的情節，與

ns = np.arange(-5, 5+1) 

更換

ns = np.linspace(-5, 5, 1000) 

也就是說，您使用的公式只在整點纔有意義。 貝塞爾函數作爲負序的函數是一個振盪函數，http://dlmf.nist.gov/10.27#E2因此情節對我來說看起來很好。