两分量高斯拟合无法正常工作答案

【问题标题】：two component gaussian fit not working correctly两分量高斯拟合无法正常工作
【发布时间】：2021-06-28 13:44:34
【问题描述】：

我正在尝试使用以下代码进行两分量高斯拟合：

def double_gaussian(velo_peak,a1, mu1, sigma1, a2, mu2, sigma2):
                      
            
            res_two = a1 * np.exp(-(velo_peak - mu1)**2/(2 * sigma1**2))  \
                      + a2 * np.exp(-(velo_peak - mu2)**2/(2 * sigma2**2))

            return res_two
#Guess parameters:
guess = [5, 115.2, 0.7, 4, 115.7, 0.7]
popt,pcov = curve_fit(double_gaussian, velo_peak, spec_peak, guess)
plt.plot(velo_peak*1e-9, double_gaussian(velo_peak, *popt), 'r-', label='fit')

我用两个分量定义了方程，然后输入了猜测值并尝试绘制，但我得到以下错误： warnings.warn('Covariance of the parameters could not be estimated', 因此，我得到了一条适合的直线。我会说我发现代码合乎逻辑，但当然，有问题。提供的任何帮助将不胜感激。

【问题讨论】：

你为什么要这么做velo_peak*1e-9？还有，velo/spec_peak里面是什么，能提供样本数据吗？
这是巨大的样本数据，我认为它不会有用，因为它们可以被视为具有值的两个数组。 velo_peak 乘以这个因子只是为了转换。
我在我生成的示例数据上尝试了您的代码，就像他们在示例 here 中所做的那样，它没有问题，没有 *1e-9。
您对最初的猜测有多大把握？我正在用越来越糟糕的猜测拟合数据，并且在某个时刻我收到警告，但拟合失败了。也许在 pastebin 或类似的东西上上传一个小的随机样本？

标签： python numpy plot curve-fitting gaussian

【解决方案1】：

问题似乎在最初的猜测中：

import matplotlib.pyplot as plt  # type: ignore
import numpy as np  # type: ignore
from scipy.optimize import curve_fit  # type: ignore

def double_gaussian(velo_peak, a1, mu1, sigma1, a2, mu2, sigma2):
    g1 = a1 * np.exp(-((velo_peak - mu1) ** 2) / (2 * sigma1 ** 2))
    g2 = a2 * np.exp(-((velo_peak - mu2) ** 2) / (2 * sigma2 ** 2))
    return g1 + g2

xdata = np.linspace(100, 120, 1000)
y = double_gaussian(xdata, 5, 115.2, 0.7, 4, 105.7, 0.7)
np.random.seed(1729)
y_noise = 0.2 * np.random.normal(size=xdata.size)
ydata = y + y_noise
plt.plot(xdata, ydata, "b-", label="data")

# Guess parameters:
# guess = [5, 115.2, 0.7, 4, 105.7, 0.7]  # ok (duh)
# guess = [1, 100, 1, 1, 100, 1]  # fails 10k
# guess = [5, 100, 1, 4, 100, 1]  # fails 10k
# guess = [5, 115, 1, 4, 105, 1]  # ok
# guess = [5, 115, 5, 4, 105, 5]  # warn
# guess = [5, 115, 2, 4, 105, 2]  # ok
# guess = [1, 115, 2, 1, 105, 2]  # ok
# guess = [1, 115, 3, 1, 105, 3]  # ok
# guess = [1, 110, 3, 1, 110, 3]  # ok
# guess = [1, 100, 3, 1, 120, 3]  # only one found
# guess = [1, 120, 3, 1, 100, 3]  # only one found
guess = [1, 118, 3, 1, 102, 3]  # ok
popt, pcov = curve_fit(double_gaussian, xdata, ydata, guess, maxfev=10000)
plt.plot(xdata, double_gaussian(xdata, *popt), "r-", label="fit")

plt.show()

如果无法访问数据，就很难提供帮助。您可能想尝试一些方法来找到峰值位置，并将其用作mu1/2 值，它可能会有所帮助。或者在起始参数空间上进行某种网格搜索，评估模型和数据之间的误差。

编辑：

我尝试增加优化器调用maxfev=100_000_000 的迭代次数，得到了以下结果：

# guess = [1, 100, 1, 1, 100, 1]  # bad 100kk
# guess = [1, 118, 1, 1, 102, 1]  # ok 100kk
# guess = [1, 110, 1, 1, 110, 1]  # bad 100kk
# guess = [1, 111, 1, 1, 109, 1]  # bad 100kk
guess = [1, 112, 1, 1, 108, 1]  # ok 100kk

同样，起始 mu 值似乎非常重要。

干杯！

【讨论】：