Python中的非对称高斯拟合

Question

我正在尝试为以下数据拟合非对称高斯：http://ge.tt/99iNaL53（csv 文件）。

我尝试使用 lmfit 的倾斜高斯模型和样条曲线，但我无法让高斯模型很好地拟合，并且样条曲线不是我想要的（我不希望样条曲线完全符合如下所示的数据，并且更改平滑级别无济于事）。

这是使用上述数据生成下图的代码。第二个图是我试图实现的目标，即从拟合中读取上升和衰减时间。

import numpy as np
import matplotlib.pyplot as plt
from scipy.interpolate import CubicSpline
from scipy.interpolate import UnivariateSpline
from lmfit.models import SkewedGaussianModel

data = np.loadtxt('data.csv', delimiter=',')

x = data[:,0]
y = data[:,1]

# Skewed Gaussian fit
model = SkewedGaussianModel()
params = model.make_params(amplitude=400, center=3, sigma=7, gamma=1)
result = model.fit(y, params, x=x)

# Cubic Spline
cs = CubicSpline(x, y)
x_range = np.arange(x[0], x[-1], 0.1)

# Univariate Spline
us = UnivariateSpline(x, y, k = 1)

# Univariate Spline (smoothed)
us2 = UnivariateSpline(x, y, k = 5)

plt.scatter(x, y, marker = '^', color = 'k', linewidth = 0.5, s = 10, label = 'data')
plt.plot(x_range, cs(x_range), label = 'Cubic Spline')
plt.plot(x_range, us(x_range), label = 'Univariate Spline, k = 1')
plt.plot(x_range, us2(x_range), label = 'Univariate Spline, k = 5')
plt.plot(x, result.best_fit, color = 'red', label = 'Skewed Gaussian Attempt')
plt.xlabel('x')
plt.ylabel('y')
plt.yscale('log')
plt.ylim(1,500)
plt.legend()
plt.show()

Answer 1

我为 J. Chem 写了一些东西。埃德。 [1] 涉及将非对称高斯函数拟合到数据，您可以在此处找到核心 repo [2] 但下面是关于我如何将data 设置为x = data[:,0] 和y = data[:,1] 到您正在使用的函数类型：

import numpy as np
from scipy.optimize import leastsq
from scipy.special import erf

initials = [6.5, 13, 1, 0] # initial guess

def asymGaussian(x, p):
    amp = (p[0] / (p[2] * np.sqrt(2 * np.pi)))
    spread = np.exp((-(x - p[1]) ** 2.0) / (2 * p[2] ** 2.0))
    skew = (1 + erf((p[3] * (x - p[1])) / (p[2] * np.sqrt(2))))
    return amp * spread * skew

def residuals(p,y,x):
    return y - asymGaussian(x, p)

# executes least-squares regression analysis to optimize initial parameters
cnsts = leastsq(
    residuals, 
    initials, 
    args=(
        data_set[:,1], # y value
        data_set[:,0]  # x value
        ))[0]

y = asymGaussian(data[:,0], cnsts)

最后只是绘制(y, data[:,0])。希望这会有所帮助！

[1] https://pubs.acs.org/doi/10.1021/acs.jchemed.9b00818
[2] https://github.com/1mikegrn/pyGC

Answer 2

这里有问题吗？实际上，我没有看到。

lmfit 的结果最适合偏斜的高斯模型。 您已选择在对数刻度上绘制结果。这完全改变了对合身质量或不合身的看法。

您似乎期望更合适，但不是*太好。好吧，看起来您的数据并不能完美地由单个倾斜的高斯表示。好像你没想到会这样。您可以为模型函数尝试不同的形式，比如偏斜的洛伦兹函数或其他形式。但是你的数据有那么低的x 肩膀，这绝对不像你未引用的数字。