【问题标题】:Scipy LeastSq errorbarsScipy LeastSq 误差线
【发布时间】:2012-10-31 16:31:12
【问题描述】:

我正在使用来自 SciPy 的 LeastSq 将实验光谱拟合到理论期望值。当然存在与实验值相关的误差。我怎样才能把这些喂给 LeastSq 或者我需要一个不同的例程?我在文档中找不到任何内容。

【问题讨论】:

  • 它们是数据格式错误(例如“null”)还是outliers?在后一种情况下,没有标准的方法来处理它。
  • 能否详细说明错误内容和错误分布情况?
  • 误差基本上是通道计数的平方。

标签: python scipy least-squares


【解决方案1】:

scipy.optimize.leastsq 函数没有内置的方法来合并权重。但是,scipy.optimize.curve_fit 函数确实有一个 sigma 参数,可用于指示每个 y 数据点的方差。

curve_fit 使用1.0/sigma 作为权重,其中sigma 可以是长度为N 的数组,(与ydata 长度相同)。

因此,您必须以某种方式根据误差线的大小推测每个 ydata 点的方差,并使用它来确定sigma

例如,如果您声明误差条的一半长度代表 1 个标准差,那么方差(curve_fit 称为 sigma)将是标准差的平方。

sigma = (length_of_error_bar/2)**2

参考:

【讨论】:

  • 我在残差函数中加入了一个权重函数,它工作正常。我可能会看一下 curve_fit 例程,但由于我实际上使用的是 lmfit(它处理参数很棒),所以我的第一种方法是加权。
  • 哦,是的,我想你可以这样做——在这种情况下,你可以使用leastsq
【解决方案2】:

我自己也在做这件事,所以我会分享我所做的,也许我们可以从社区获得一些 cmets。我有一组在确定的时间间隔内采集的数据点,我从中计算了标准偏差。我想用 sin 函数来拟合这些点。 Leastsq 通过基于一组参数 p 最小化残差或数据点与拟合函数之间的差异来做到这一点。我们可以通过将残差除以方差或标准差的平方来加权残差。

如下:

from scipy.optimize import leastsq
import numpy as np
from matplotlib import pyplot as plt
def sin_func(t, p):
    """ Returns the sin function for the parameters:
        p[0] := amplitude
        p[1] := period/wavelength
        p[2] := phase offset
        p[3] := amplitude offset
    """
    y = p[0]*np.sin(2*np.pi/p[1]*t+p[2])+p[3]
    return y

def sin_residuals(p, y, t, std):
    err = (y - p[0]*np.sin(2*np.pi/p[1]*t+p[2])-p[3])/std**2
    return err

def sin_fit(t, ydata, std, p0):
""" Fits a set of data, ydata, on a domain, t, with individual standard 
deviations, std, to a sin curve given the initial parameters, p0, of the form:
        p[0] := amplitude
        p[1] := period/wavelength
        p[2] := phase offset
        p[3] := amplitude offset
    """
    # optimization # 
    pbest = leastsq(sin_residuals, p0, args=(ydata, t, std), full_output=1)
    p_fit = pbest[0]

    # fit to data #
    fit = p_fit[0]*np.sin(2*np.pi/p_fit[1]*t+p_fit[2])+p_fit[3]

    return p_fit

【讨论】:

    猜你喜欢
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 2014-11-21
    • 2011-11-27
    • 2012-04-10
    • 1970-01-01
    • 2017-07-22
    • 2015-10-04
    相关资源
    最近更新 更多