Реализация дискретного гауссовского ядра в Python?

Question

Я ищу для реализации 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]) 

Выход выглядит следующим образом:

Из статьи Википедии, она должна выглядеть следующим образом:

Я предполагаю, что я делаю какую-то очень пустую ошибку. :/ Есть предположения? Спасибо!

EDIT: То, что я написал, эквивалентно scipy.special.ive(n, t). Я почти уверен, что это должна быть модифицированная функция Бесселя первого рода, а не вторая, но кто-то может подтвердить это?

Answer 1

Если вы хотите получить википедии сюжет, заменить

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

с

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

А именно, формула используется только имеет смысл в целых точках. Функция Бесселя как функция отрицательного порядка - это осциллирующая функция, http://dlmf.nist.gov/10.27#E2, поэтому сюжет отлично выглядит для меня.