不同拟合类型的 scipy.odr 输出不一致？答案

【问题标题】：Inconsistency in scipy.odr output for different fit types?不同拟合类型的 scipy.odr 输出不一致？
【发布时间】：2020-10-24 14:53:54
【问题描述】：

对于odr.set_job(fit_type=2)（普通最小二乘法）odr.set_job(fit_type=0)（显式 ODR），对于相同数据的两种拟合类型，不确定性 sd_beta 和协方差矩阵 cov_beta 存在显着差异 - 这在y 部分，仅。

有人知道这是从哪里来的吗？下面是一个重现效果的代码示例。另一方面：对于sd_beta 和cov_beta 之间的隐藏比例因子，我在这里找到了一些帮助：How to compute standard error from ODR results?

import numpy as np
import scipy.odr as spodr

def linear_func(pars, x):
    """Define a (linear) function to fit the data with."""
    a, b = pars
    return b * x + a

x = np.asfarray( [1.01585124, 2.03965699, 2.89148684] )
y = np.asfarray( [1.16201599, 2.15275179, 2.87248892] )
dy = np.asfarray( [0.15520856, 0.33160654, 0.45461028] )

real_data = spodr.RealData(x, y, sy=dy)
linear_model = spodr.Model(linear_func)
odr = spodr.ODR(real_data, linear_model, beta0=[0., 1.])
# fit_type=2 means least square fit
odr.set_job(fit_type=2)
out = odr.run()
print('fit_type=2 -------------------------------------------')
out.pprint()
print('\n')

# fit_type=0 means explicit ODR
odr.set_job(fit_type=0)
out = odr.run()
print('fit_type=0 -------------------------------------------')
out.pprint()

给出以下结果：

fit_type=2 -------------------------------------------
Beta: [0.22131257 0.92950773]
Beta Std Error: [0.0444006 0.03004  ]
Beta Covariance: [[ 0.10678885 -0.06586652]
 [-0.06586652  0.04888189]]
Residual Variance: 0.01846085213372428
Inverse Condition #: 0.02685695006039764
Reason(s) for Halting:
  Sum of squares convergence


fit_type=0 -------------------------------------------
Beta: [0.24510628 0.91688192]
Beta Std Error: [0.07210707 0.03482789]
Beta Covariance: [[ 2.29687687 -1.03022278]
 [-1.03022278  0.53584146]]
Residual Variance: 0.0022636956774060857
Inverse Condition #: 0.017439127336256896
Reason(s) for Halting:
  Sum of squares convergence

【问题讨论】：

标签： python numpy scipy regression curve-fitting

【解决方案1】：

事实证明，这不是 buggy 不一致，而更像是一个陷阱。当然，对于普通最小二乘类型 (OLS)，x 误差实际上为零，因为它们只是没有出现在数学模型中。

现在，认为在 ODR 案例中省略 sx 意味着 sx → 0 是不正确的。实际上对于sx=None，sx默认为unity，见here。与 OLS 案例相比，这当然会产生不同的结果。

在这里，如下所示显式设置sx=data.x*ResEps 会得到相同的结果。

至少，应该在文档中更明确地指出这个上下文，我会尝试为此准备一个问题。有趣的是，在代码中，sx 和 covx 的同时传递会引发警告（这是有道理的），但如果两者都没有给出，则会发生上述情况。

# obtain the machine precision for floating point operations
MachineEps = np.finfo('float').eps
ResEps = np.sqrt(MachineEps)

real_data = spodr.RealData(x, y, sy=dy, sx=x*ResEps)
linear_model = spodr.Model(linear_func)
odr = spodr.ODR(real_data, linear_model, beta0=[0., 1.])

# fit_type=0 means explicit ODR
odr.set_job(fit_type=0)
out = odr.run()
print('fit_type=0 -------------------------------------------')
out.pprint()

产量

fit_type=0 -------------------------------------------
Beta: [0.22131258 0.92950773]
Beta Std Error: [0.04440071 0.03004007]
Beta Covariance: [[ 0.10678941 -0.06586688]
 [-0.06586688  0.04888212]]
Residual Variance: 0.01846085213371618
Inverse Condition #: 0.02685687854396314
Reason(s) for Halting:
  Sum of squares convergence

【讨论】：