如何在Python中使用minuit.Minuit拟合高斯曲线？

Question

我正在尝试用minuit.Minuit函数来给一些简单的数据拟合一个高斯曲线，但我的参数没有任何变化。如果有人能帮我，我会非常感激。

import numpy as np
import minuit

xCurve = np.array([0,1,2,3,4,5,6,7,8,9])
yCurve = np.array([0,1,2,3,4,5,4,3,2,1])


def Gaus(a,b,c):
    return a*np.exp(-((xCurve-b)**2)/(2*c**2))

m = minuit.Minuit(Gaus,a=4.5,b=5,c=0.4)
m.printMode=1
m.migrad()
m.printMode=0
m.values()

a = m.values['a']
b = m.values['b']
c = m.values['c']
d = m.values['d']
print a
print b
print c
print d

它报了一个错误：minuit.MinuitError: 协方差不是正定的。

Answer 1

Minuit 是一个用来最小化的工具，但你给它的是一个拟合函数，而不是目标函数。（这个函数实际上不是正定的，所以这个错误提示是合理的。）

要得到你真正想要的结果，可以这样做：

def gauss(x, a,b,c):
    return a*np.exp(-((x-b)**2/(2*c**2)))

def minimizeMe(a,b,c):
    return sum((gauss(x, a,b,c) - y)**2 for x, y in zip(xCurve, yCurve))

m = minuit.Minuit(minimizeMe, a=4.5, b=5, c=0.4)
m.printMode = 1
m.migrad()

这样做虽然没有充分利用你的 Numpy 数组，但如果你把最小化和拟合的函数结合起来，应该可以通过 ufunc 来实现。

PyMinuit 旨在提供更底层的拟合技术访问。如果你只对普通的最小二乘法感兴趣，可能会觉得直接使用 Minuit 的接口有点麻烦。另一方面，如果你打算用套索回归来限制一些参数，提供非二次或甚至非对称的损失函数，或者如果你计划进行一些无法简单表示为函数拟合的优化，那么这个底层接口就会很有帮助。