Python: размер результирующей функции свертки двух гауссианов с np.convolve

Question

Мне интересно оптимизировать функцию, которая является сверткой двух функций.Основная проблема заключается в том, что моя результирующая функция полностью масштабируется, и я не понимаю, что на самом деле делает np.convolve.

Я написал небольшой скрипт, который должен объединить два гаусса, но получающийся в результате гауссиан намного больше по размеручем функции ввода:

from scipy.interpolate import interp1d
import matplotlib.pyplot as plt
import numpy as np


# https://stackoverflow.com/questions/18088918/combining-two-gaussians-into-another-guassian
def gauss(x, p): # p[0]==mean, p[1]==stdev, p[2]==heightg, p[3]==baseline                   
    a = p[2]
    mu = p[0]
    sig = p[1]
    #base = p[3]
    return a * np.exp(-1.0 * ((x - mu)**2.0) / (2.0 * sig**2.0)) #+ base

p0 = [0, 0.3, 1] # Inital guess is a normal distribution
p02 = [0, 0.2, 0.5]


xp = np.linspace(-4, 4, 2000)
convolved = np.convolve(gauss(xp, p0),gauss(xp, p02), mode="same")


fig = plt.figure()



plt.subplot(2, 1, 1)

plt.plot(xp, gauss(xp, p0), lw=3, alpha=2.5)
plt.plot(xp, gauss(xp, p02), lw=3, alpha=2.5)
plt.xlim([-2, 2])

plt.subplot(2, 1, 2)
plt.plot(xp, gauss(xp, p0), lw=3, alpha=2.5)
plt.plot(xp, gauss(xp, p02), lw=3, alpha=2.5)
plt.plot(xp, convolved, lw=3, alpha=2.5,label="too damn high?")
plt.legend()
plt.xlim([-2, 2])

plt.tight_layout()
plt.show()

Результирующий гауссиан после свертки намного выше

, чем я ожидал (Википедия):

Answer 1

Вы должны перенормировать dx между двумя тиками x.

Numpy заменяет суммирование на интеграцию, но поскольку функция принимает только значения Y, ей не нужно учитывать элемент громкости при интеграции.ось, которую вам нужно включить вручную.
Мне тоже приходилось сталкиваться с этой проблемой, и это неприятно, когда вы начинаете делать что-то с dx = 1, и внезапно вы получаете неправильный результат из-за разныхРаспределение по оси X.

xp = np.linspace(-4, 4, 2000)
dx = xp[1] - xp[0]
convolved = np.convolve(gauss(xp, p0),gauss(xp, p02), mode="same") * dx

!! ПРИМЕЧАНИЕ: не помещайте перенормировку в определение функции.dx должен учитываться только один раз, поскольку интеграл входит в сумму.Если вы поместите его внутрь функции, он будет фактически посчитан дважды, потому что с помощью него генерируется бот-гауссиан.

PS: попытаться понять это лучшевы можете попытаться сгенерировать данные оси x с разным интервалом и без перенормировки вы увидите, что высота вашей свертки будет отличаться (чем меньше расстояние, тем больше высота)

fig = plt.figure()
ax = fig.add_subplot(111)
for spacing in (100,500,1000,2000):
    spacing += 1
    xp = np.linspace(-4, 4, spacing)
    dx = xp[1] - xp[0]
    convolved = np.convolve(gauss(xp, p01),gauss(xp, p02), mode="same") * dx
    ax.plot(xp, convolved, lw=3, alpha=2.5,label="spacing = {:g}".format(8/spacing))

ax.set_title("Convolution with different x spacing. With renormalization")
fig.legend()
plt.show()