Nor*_*add 3 python numpy gaussian convolution
我有兴趣优化一个函数,该函数是两个函数的卷积。主要问题是我的结果函数完全是规模化的,我不明白 np.convolve 实际做了什么。
我编写了一个小脚本,应该对两个高斯函数进行卷积,但生成的高斯函数的大小比输入函数大得多:
from scipy.interpolate import interp1d
import matplotlib.pyplot as plt
import numpy as np
# /sf/ask/1266224291/
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()
Run Code Online (Sandbox Code Playgroud)
卷积后得到的高斯要高很多
超出我的预期(维基百科):
您必须重新规范化两个 x 刻度之间的 dx。
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
Run Code Online (Sandbox Code Playgroud)
!!注意:不要将重整化放在函数定义中。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()
Run Code Online (Sandbox Code Playgroud)