Scipy & Optimize：最小化示例，如何添加约束？

Question

我有一个数据集，我想通过最小二乘误差法找到一个混合高斯模型。

代码是这样的：

from sklearn.neighbors import KernelDensity
kde = KernelDensity().fit(sample)

def gaussian_2d(x,y,meanx,meany,sigx,sigy,rho):
    # rho <= 1
    part1 = 1/(2*np.pi*sigx*sigy*sqrt(1-0.5**2))
    part2 = -1/2*(1-rho**2)
    part3 = (((x-meanx)/sigx)**2-2*rho*(x-meanx)*(y-meany)/(sigx*sigy)+((y-meany)/sigy)**2)
    return part1*exp(part2*part3)

def square_error(f1,f2, u1,v1,sigu1,sigv1,rho1, u2,v2,sigu2,sigv2,rho2, u3,v3,sigu3,sigv3,rho3):
    # 1. Generate Mixed Gaussian Model
    def gaussian1(x,y):
        return gaussian_2d(x,y,u1,v1,sigu1,sigv1,rho1)
    def gaussian2(x,y):
        return gaussian_2d(x,y,u2,v2,sigu2,sigv2,rho2)
    def gaussian3(x,y):
        return gaussian_2d(x,y,u3,v3,sigu3,sigv3,rho3)
    mixed_model = f1*gaussian1(x,y)+f2*gaussian2(x,y)+(1-f1-f2)*gaussian3(x,y)
    # 2. Calculate the sum of square error
    sum_error = 0
    for point in sample:
        error = (exp(mixed_model(point)) - exp(kde.score(point)))**2
        sum_error += error
    return sum_error

# How can I add constraints:
# f1+f2 <= 1
# rho1,2,3 <= 1
result = sp.optimize.minimize(square_error)

但我不知道如何在minimize方法中添加约束。http://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.minimize.html#scipy.optimize.minimize 中的示例很难理解。

更新：这就是我最终的结果，

result = sp.optimize.minimize(
    square_error,
    x0 = [0.2,0.5,
          1,1,1,1,0.3,
          1,1,1,1,0.3,
          1,1,1,1,0.3,],
    bounds = [(0., 1.),(0., 1.),
              (None, None),(None, None),(0., None),(0., None),(0., 1.),
              (None, None),(None, None),(0., None),(0., None),(0., 1.),
              (None, None),(None, None),(0., None),(0., None),(0., 1.),])

但它给了我TypeError: square_error() takes exactly 17 arguments (1 given)，有什么问题？

Answer 1

如果解算器支持，则只能添加添加边界，因此仅适用于method='L-BFGS-B',TNC和SLSQP。

边界通过元组序列传递(min, max)，其长度对应于参数的数量。拟合 3 个参数的示例如下：

result = sp.optimize.minimize(
                        square_error,
                        method='L-BFGS-B',
                        bounds=[(0., 5.), (None, 1.e4), (None, None)])

这里，None对应于无界。恐怕在inf1+f2 <= 1的框架内不可能对您的示例中的参数组合进行约束。boundsscipy.minimize

但是，np.inf如果违反了界限，您可以简单地返回成本函数。不过，我不确定它的稳定性：

def square_error(f1,f2, other_args):
    if f1+f2 <= 1:
        return np.inf
    # rest of the cost function

此外，我建议使用多元高斯函数的 python 实现，而不是从头开始创建它们。这将加快您的安装速度，有助于避免错误并且更具可读性。