Python lmfit在加权合适后缩小了卡方太小

Question

我正在使用一些测试数据并使用以下代码在Python 2.7中运行适配器lmfit。我需要权重为1/y（使用Leven-Marq。例程）。我已经定义的权重，并在这里使用它们：

from __future__ import division 
from numpy import array, var 
from lmfit import Model 
from lmfit.models import GaussianModel, LinearModel 

import matplotlib.pyplot as plt 
import seaborn as sns 

xd = array([1267, 1268, 1269, 1270, 1271, 1272, 1273, 1274, 1275, 1276, 
    1277, 1278, 1279, 1280, 1281, 1282, 1283, 1284, 1285, 1286, 1287, 1288, 
    1289, 1290, 1291, 1292, 1293, 1294, 1295, 1296, 1297, 1298, 1299, 1300, 
    1301, 1302, 1303, 1304, 1305, 1306, 1307, 1308, 1309, 1310, 1311, 1312, 
    1313, 1314, 1315, 1316, 1317, 1318, 1319, 1320, 1321, 1322, 1323, 1324, 
    1325, 1326, 1327, 1328, 1329, 1330, 1331, 1332, 1333, 1334]) 
yd = array([238, 262, 255, 271, 270, 281, 261, 278, 280, 254, 289, 285, 304, 314, 
    329, 342, 379, 450, 449, 564, 613, 705, 769, 899, 987, 1043, 1183, 1295, 1298, 
    1521, 1502, 1605, 1639, 1572, 1659, 1558, 1476, 1397, 1267, 1193, 1016, 951, 
    835, 741, 678, 558, 502, 480, 442, 399, 331, 334, 308, 283, 296, 265, 264, 
    273, 258, 270, 262, 263, 239, 263, 251, 246, 246, 234]) 

mod = GaussianModel() + LinearModel() 
pars = mod.make_params(amplitude=25300, center=1299, sigma=7, slope=0, intercept=450) 
result = mod.fit(yd, pars, method='leastsq', x=xd, weights=1./yd) 
rsq = 1 - result.residual.var()/var(yd) 
print(result.fit_report()) 
print rsq 

plt.plot(xd, yd,   'bo', label='raw') 
plt.plot(xd, result.init_fit, 'k--', label='Initial_Guess') 
plt.plot(xd, result.best_fit, 'r-', label='Best') 
plt.legend() 
plt.show() 

的输出是：

[[Model]] 
    (Model(gaussian) + Model(linear)) 
[[Fit Statistics]] 
    # function evals = 27 
    # data points  = 68 
    # variables  = 5 
    chi-square   = 0.099 
    reduced chi-square = 0.002 
    Akaike info crit = -434.115 
    Bayesian info crit = -423.017 
[[Variables]] 
    sigma:  7.57360038 +/- 0.063715 (0.84%) (init= 7) 
    center:  1299.41410 +/- 0.071046 (0.01%) (init= 1299) 
    amplitude: 25369.3304 +/- 263.0961 (1.04%) (init= 25300) 
    slope:  -0.15015228 +/- 0.071540 (47.65%) (init= 0) 
    intercept: 452.838215 +/- 93.28860 (20.60%) (init= 450) 
    fwhm:  17.8344656 +/- 0.150037 (0.84%) == '2.3548200*sigma' 
    height:  1336.33919 +/- 17.28192 (1.29%) == '0.3989423*amplitude/max(1.e-15, sigma)' 
. 
. 
. 
. 
0.999999993313 

最后一行（刚好高于此处，或立即plt.plot(xd, yd, 'bo', label='raw')之前）是R^2，将所得配合附在这里。 。

R^2和输出的视觉检查表明这是一个合理的拟合。我期待1.00的订单减少卡方（source）。但是，降低的卡方值的返回值比1.00小几个数量级。

由于默认值是no weightslmfit我需要一个加权拟合，我已经定义了权重，但我认为我需要以不同的方式指定它们。我的怀疑是这种重量的规格可能会导致减小的卡方非常小。

是否有不同的方式来指定权重或其他参数，使得曲线拟合后的减少的卡方接近或等于1.00的相同量级？

Answer 1

lmfit中的权重是残差在最小二乘意义上被最小化的乘法因子。也就是说，它取代

residual = model - data 

与

residual = (model - data) * weights 

一种常见的做法，一个是我想你可能会打算，是说的权重应为1.0/variance_in_data，因为这是通常意味着要在1左右减少卡路里，以达到最佳效果，因为您要链接的优秀文章将对此进行讨论。

正如那里所讨论的那样，问题是确定数据中的方差。对于很多情况，例如当信号由统计统计支配时，数据的差异可以估计为sqrt(data)。这忽略了许多噪音源，但通常是一个很好的起点。碰巧，我相信使用

result = model.fit(..., weights=np.sqrt(1.0/yd)) 

会导致您的情况下约0.8的卡方减少。我认为这可能是你想要的。

此外，为了阐明相关的一点：您链接的写法讨论了缩小卡方远离1时拟合参数的不确定性.Lmfit在默认情况下进行缩放（scale_covar选项可关闭此功能），因此改变权重的比例不会改变参数sigma,center等中的不确定性的比例。不确定性（和最佳拟合值）的值将改变一些，因为权重的改变改变每个数据点的重点，但最佳拟合值不会有太大变化，并且即使您对数据方差的估计（以及如此reduced chi-square）少数几个数字，估计的不确定性也应保持相同的数量级数量级。

也就是说，将脚本更改为使用weights=1.0/np.sqrt(yd)会使卡方的卡方值接近1，但它不会很大程度上改变拟合变量的不确定性。