0

我正在使用 lmfit 将“复杂”模型拟合到实验数据中。

我这样定义模型的参数：

def model(t, y0, h, A1, x1, s1, A2, x2, s2, A3, x3, s3, A4, x4, s4):

        y = y0 + h * (
                      A1 * (0.5 + 0.5 * erf((t - x1) / (s1 * 2**0.5))) +
                      A2 * (0.5 + 0.5 * erf((t - x2) / (s2 * 2**0.5))) +
                      A3 * (0.5 + 0.5 * erf((t - x3) / (s3 * 2**0.5))) +
                      A4 * (0.5 + 0.5 * erf((t - x4) / (s4 * 2**0.5)))
                      )
        return y


mod = Model(model)

# sys.setrecursionlimit(150000)

# Set the parameters, define their bounds and their relations
mod.set_param_hint('y0', value=self.spin_y0.value())
mod.set_param_hint('h', value=self.spin_h.value())

mod.set_param_hint('A1', value=self.spin_A1.value(), min=0, max=1),
mod.set_param_hint('x1', value=self.spin_x1.value(), min=0, max=self.x[-1])
mod.set_param_hint('s1', value=self.spin_s1.value(), min=0)

mod.set_param_hint('A2', value=self.spin_A2.value(), min=0, max=1),
mod.set_param_hint('x2', value=self.spin_x2.value(), min=0, max=self.x[-1])
mod.set_param_hint('s2', value=self.spin_s2.value(), min=0)

mod.set_param_hint('A3', value=self.spin_A3.value(), min=0, max=1),
mod.set_param_hint('x3', value=self.spin_x3.value(), min=0, max=self.x[-1])
mod.set_param_hint('s3', value=self.spin_s3.value(), min=0)

mod.set_param_hint('A4', value=self.spin_A4.value(), min=0, max=1),
mod.set_param_hint('x4', value=self.spin_x4.value(), min=0, max=self.x[-1])
mod.set_param_hint('s4', value=self.spin_s4.value(), min=0)

mod.set_param_hint('A1', expr='1-A2-A3-A4')
mod.set_param_hint('A2', expr='1-A1-A3-A4')
mod.set_param_hint('A3', expr='1-A1-A2-A4')
mod.set_param_hint('A4', expr='1-A1-A2-A3')

mod.make_params()

# Fit the data !
result = mod.fit(self.y, t=self.x)

我使用 spinBoxes 来获取将用作初始值的值（我的程序是 GUI）。

我对 A1、A2、A3 和 A4 有疑问。它们是依赖的，它们的总和必须等于 1。

如果我评论

mod.set_param_hint('A1', expr='1-A2-A3-A4')
mod.set_param_hint('A2', expr='1-A1-A3-A4')
mod.set_param_hint('A3', expr='1-A1-A2-A4')
mod.set_param_hint('A4', expr='1-A1-A2-A3

参数在 0 和 1 之间变化（我在定义参数时定义了这个范围），但它们的总和永远不会等于 1。

当这些行没有被注释时，我有以下异常：

File "/usr/lib/python3.4/site-packages/lmfit/minimizer.py", line 255, in __update_paramval
    self.__update_paramval(dep)
  File "/usr/lib/python3.4/site-packages/lmfit/minimizer.py", line 250, in __update_paramval
    if getattr(par, 'expr', None) is not None:
RuntimeError: maximum recursion depth exceeded

我真的不知道这意味着什么。我试图增加最大递归深度，但没有任何改变。

您对问题所在有任何想法吗？

4

1 回答 1

0

好的，对于未来的我和其他寻找解决方案的人，我找到了一个。我使用 scimin.fmin_slsqp 来完成这项工作。这是一个代码示例：

def model(t, params):

    y0, h, A1, x1, s1, A2, x2, s2, A3, x3, s3, A4, x4, s4 = [value for value in params]

    y = y0 + h * (
                  A1 * (0.5 + 0.5 * erf((t - x1) / (s1 * 2**0.5))) +
                  A2 * (0.5 + 0.5 * erf((t - x2) / (s2 * 2**0.5))) +
                  A3 * (0.5 + 0.5 * erf((t - x3) / (s3 * 2**0.5))) +
                  A4 * (0.5 + 0.5 * erf((t - x4) / (s4 * 2**0.5)))
                 )

    return y


def sum_residuals(params):  # the function we want to minimize

    y0, h, A1, x1, s1, A2, x2, s2, A3, x3, s3, A4, x4, s4 = [value for value in params]

    return sum((list_y - model(list_x, params))**2)


def constraints(params):

    y0, h, A1, x1, s1, A2, x2, s2, A3, x3, s3, A4, x4, s4 = [value for value in params]

    return 1 - (A1 + A2 + A3 + A4)


pfit = scimin.fmin_slsqp(sum_residuals,
                         p0,
                         f_eqcons=constraints,
                         bounds=bounds,
                         iprint=0,
                         iter=300
                        )

其中 p0 是初始参数的序列。但是， fmin_slsqp 不返回协方差矩阵，我无法计算拟合参数的置信区间。为此，我使用了引导方法。

于 2015-03-29T21:26:48.020 回答