4

我有以下由变量参数化的优化代码n。

from symfit import parameters, Eq, Ge, Fit
import numpy as np
n = 3
xdata = np.sort(np.random.choice(range(1, 4*n), n)) # Make fake data
print(xdata)
p1, p2, p3 = parameters('p1, p2, p3')
model = p1*p2*p3
constraints = [
    Eq(xdata[0]*p1+(xdata[1]-xdata[0])*p2+(xdata[2]-xdata[1])*p3, 1),
    Ge(p1, p2),
    Ge(p2, p3),
    Ge(p3, 0)
    ]

fit = Fit(- model, constraints=constraints)
fit_result = fit.execute()
print(fit_result)

我想将它用于更大的值，n但我不知道如何更改线条

 p1, p2, p3 = parameters('p1, p2, p3')
 model = p1*p2*p3

而且constraints要应付一个任意大的n。

该代码正在使用symfit库。该链接显示了如何使用的示例以及parameters文档的链接。

怎么能做到这一点？

4

3 回答 3

4

Numpy与库的交互非常好。symfit您试图概括的所有操作在使用它时都相当简单。

设置

n = 3
_data = np.sort(np.random.choice(np.arange(1, 4 * n), n))

  1. 字符串格式

    tuple您可以使用简单的迭代器动态创建 a of 参数 and str.join，然后您可以将其传递给parameters构造函数以获取 a tupleof 您的参数。

params = parameters(', '.join(f'p{i}' for i in range(1, n+1)))
                ^^
# p1, p2, p3 = parameters('p1, p2, p3')

  1. np.prod

    这个操作非常简单。 np.prod计算：

    给定轴上数组元素的乘积

    当应用于一个tuple参数symfit时，会产生您想要的p1*p2*...pn

model = np.prod(params)
        ^^
# model = p1*p2*p3

  1. np.concatenate+np.diff

    可能是最复杂的概括线，但仍然不太难理解。您想将数据数组中连续元素的差异乘以您的参数，并对结果求和。由于第一个元素与前一个元素没有区别，您可以使用np.concatenate重新添加它。

u = np.concatenate((_data[:1], np.diff(_data)))
c1 = Eq(np.sum(params * u), 1)
              ^^
# Eq(xdata[0]*p1+(xdata[1]-xdata[0])*p2+(xdata[2]-xdata[1])*p3, 1)

  1. np.column_stack

    您希望将参数滚动查看为约束：p1-p2, p2-p3, ... pn, 0。这只是堆叠一个一次性元组，用原始tuple参数填充零，然后使用列表推导解包到您的Ge构造函数中。

ges = [Ge(*group) for group in np.column_stack((params, params[1:] + (0,)))]

  1. 合身！

    我没有在这里改变任何东西！

constraints = [c1, *ges]
fit = Fit(- model, constraints=constraints)
fit_result = fit.execute()
于 2019-04-08T20:54:17.463 回答
3

参数字符串需要以动态方式从n

paramstr = ', '.join(['p{}'.format(i) for i in range(1, n)])
# when n=1, paramstr = 'p1, p2, p3'

用作函数paramstr的参数parameters。

paramvals = parameters(paramstr)

model可以通过减少paramvals其产品来重构。


from functools import reduce

model = reduce(lambda x, y: x * y, paramvals, 1)

现在到了有趣的地方！constraints可以重构为：

eqs = xdata[0] * paramvals[0] + sum(
    (xdata[i] - xdata[i-1]) * paramvals[i]
    for i in range(1, n)
)

ges = [
    Ge(paramvals[i-1], paramvals[i])
    for i in range(1, n)
]

ges.append(
    Ge(paramvals[-1], 0)
)

constraints = [
    Eq(eqs, 1),
    *ges
]
于 2019-04-08T20:45:19.843 回答
2

我对 Symfit 一无所知，但如果您只是想将上述代码推广到任意 N，那么：

  • 您可以任意生成一个看起来像"p1, p2, p3"任意 N 的字符串，并将其解析为参数列表：
params_string = ", ".join("p{}".format(i + 1) for i in range(n))
params = parameters(params_string)

解析字符串以获取参数列表的概念对我来说听起来很糟糕，我确信有更好的方法来以编程方式声明一堆参数，但这将尽可能接近您的原始代码正在做的事情。

编辑：查看 Symfit 文档，它似乎parameters(s)只是一个捷径，你实际上可以这样做：

params = [Parameter("p{}".format(i + 1)) for i in range(n)]

这不需要您构建自己的所有参数名称的连接字符串，以便 Symfit 可以将它们拆分回单独的参数名称。这还允许您为参数定义其他属性，例如它们的初始值或它们的最小/最大界限。

  • 您可以概括Eq约束：
coeffs = [xdata[0]] + [(xdata[i+1] - xdata[i]) for i in range(n-1)]
eq_constraint = Eq(sum(param * coeff for param, coeff in zip(params, coeffs), 1)

或者，作为另一个答案，使用 numpy 操作：

coeffs = np.concat(xdata[:1], np.diff(xdata))
eq_constraint = Eq(np.sum(params * coeffs), 1)
  • 您可以概括Ge约束：
ge_constraints = [Ge(params[i + 1], params[i]) for i in range(n - 1)] + [Ge(params[-1], 0]

constraints = [eq_constraint] + ge_constraints

这可以再次使用 numpy 操作来完成，但我将把它留给@user3483203 的回答。

  • 您可以使用以下方法将所有参数相乘reduce：
model = reduce(lambda x, y: x * y, params, 1)

或使用numpy.prod：

model = np.prod(params)

这应该足以将上述内容推广到任意 N。

于 2019-04-08T20:53:41.660 回答