Scipy & Optimize: 最小化示例,如何添加约束?
Scipy & Optimize: Minimize example, how to add constraints?
我有一个数据集,我想通过最小二乘误差法找到一个混合高斯模型。
代码是这样的:
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)
,有什么问题吗?
如果求解器支持,您只能添加边界,因此仅适用于 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
对应的是没有bound。
恐怕在 scipy.minimize
.
中的 bounds
框架内,无法对您示例中的 f1+f2 <= 1
等参数组合进行约束
但是,如果您的边界被侵犯,您可以在成本函数中简单地 return np.inf
。不过,我不确定它的稳定性:
def square_error(f1,f2, other_args):
if f1+f2 <= 1:
return np.inf
# rest of the cost function
此外,我建议使用 multivariate Gaussians 的 python 实现,而不是从头开始创建它们。这将加快您的拟合速度,有助于避免错误并更具可读性。
我有一个数据集,我想通过最小二乘误差法找到一个混合高斯模型。
代码是这样的:
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)
,有什么问题吗?
如果求解器支持,您只能添加边界,因此仅适用于 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
对应的是没有bound。
恐怕在 scipy.minimize
.
bounds
框架内,无法对您示例中的 f1+f2 <= 1
等参数组合进行约束
但是,如果您的边界被侵犯,您可以在成本函数中简单地 return np.inf
。不过,我不确定它的稳定性:
def square_error(f1,f2, other_args):
if f1+f2 <= 1:
return np.inf
# rest of the cost function
此外,我建议使用 multivariate Gaussians 的 python 实现,而不是从头开始创建它们。这将加快您的拟合速度,有助于避免错误并更具可读性。