我有一个排序随机数的直方图和高斯叠加 https://i.stack.imgur.com/dBqcU.png。直方图表示每个箱的观察值(将此基本情况应用于更大的数据集),高斯是拟合数据的尝试。显然,这个高斯并不代表直方图的最佳拟合。下面的代码是高斯的公式。
normc, mu, sigma = 30.845, 50.5, 7 # normalization constant, avg, stdev
gauss = lambda x: normc * exp( (-1) * (x - mu)**2 / ( 2 * (sigma **2) ) )
我计算了每个箱的期望值(曲线下的面积)并计算了每个箱的观察值的数量。有多种方法可以找到“最佳”配合。我关心的是通过最小化卡方来实现最佳拟合。在这个卡方公式中 https://i.stack.imgur.com/o7kGN.gif,期望值是每个 bin 的曲线下面积,观测值是每个 bin 排序数据值的出现次数。因此,我想在给定值附近波动normc、mu和sigma,以找到normc、mu和sigma的正确组合,从而产生最小的卡方,因为这些将是我可以插入到上面的代码中进行覆盖的参数我的直方图上最适合的高斯。我正在尝试使用 scipy 模块来最小化我的卡方如本例所示 https://stackoverflow.com/questions/13670333/multiple-variables-in-scipys-optimize-minimize。由于我需要波动参数,因此我将使用函数 gauss(上面定义)来绘制高斯叠加图,并将定义一个新函数来查找最小卡方。
def gaussmin(var,data):
# var[0] = normc
# var[1] = mu
# var[2] = sigma
# data is the sorted random numbers, represents unbinned observed values
for index in range(len(data)):
return var[0] * exp( (-1) * (data[index] - var[1])**2 / ( 2 * (var[2] **2) ) )
# I realize this will return a new value for each index of data, any guidelines to fix?
在此之后,我陷入困境。如何改变参数来找到产生最佳拟合的normc、mu、sigma?我最后一次尝试的解决方案如下:
var = [normc, mu, sigma]
result = opt.minimize(chi2, [normc,mu,sigma])
# chi2 is the chisquare value obtained via scipy
# chisquare input (a,b)
# where a is number of occurences per bin, b is expected value per bin
# b is dependent upon normc, mu, sigma
print(result)
# data is a list, can I keep it as a constant and only fluctuate parameters in var?
网上有很多标量函数的示例,但我找不到变量函数的示例。
PS - 到目前为止我可以发布完整的代码,但它有点长。如果您想查看它,请询问,我可以将其发布在这里或提供 googledrive 链接。