我有一个由 f(x_0, x_1) = a(x_1 - x_0^2)^2 + (b - x_0)^2 定义的函数,其中 a 和 b 是一些参数:
def f(x):
return a*(x[1]-x[0]**2)**2+(b-x[0])**2
其中
x=np.array([x_0,x_1])
是一个 numpy 数组。那么 f 的梯度和 Hessian 矩阵都很容易找到,由下式给出:
def f_der(x):
return np.array([-4*a*x[0]*(x[1]-x[0]**2)-2*(b-x[0]), 2*a*(x[1]-x[0]**2)])
def f_hess(x):
a_11 = 12*a*x[0]**2 - 4*a*x[1] + 2
a_12 = -4*a*x[0]**2
a_21 = -4*a*x[0]**2
a_22 = 2*a
return np.array([[a_11, a_12], [a_21, a_22]])
现在我有一个 Pandas 数据框df
,它记录了 a、b 的不同值:
a b
1 2
3 7
4 12
我想创建两个名为 x_0*
和
x_1*
的新列,它们分别由 x_0 和 x_1 的值给出,这会在约束条件下最小化 f:
0 <= x_0 + ax_1 <= 1
0 <= ax_0 - x_1 <= b
0 <= x_0 <= 1
0 <= x_1 <= b
我了解 scipy 优化包,我可以一次使用一对 a、b 来完成它:
import numpy as np
import scipy
from scipy.optimize import minimize
from scipy.optimize import Bounds
from scipy.optimize import LinearConstraint
a = 1
b = 2
bounds = Bounds([0, 1], [a, b])
linear_constraint = LinearConstraint([[0, a], [a, -1]], [0, 0], [1, b])
res = minimize(f, np.array([0,0]), method='trust-constr',
jac=f_der, hess=f_hess,
constraints=linear_constraint,
options={'verbose': 1}, bounds=bounds)
res.x
但是我的原始数据框很大(〜100,000行),所以我想问是否有一种快速方法可以将最小化器一次性应用于整个数据框并快速完成。所以期望的结果是这样的:
a b x_0* x_1*
1 2 1.00000001 0.99999999
3 7 1.12649598 0.4
4 12 1.19393107 0.29411765
ab = pd.DataFrame({
'a': (1, 3, 4),
'b': (2, 7, 12),
})
x0 必须是平的,并且需要传入 a
和
b
作为参数。对于
minimize
,这看起来像
minimize(
fun=f, jac=f_der, x0=np.zeros(2*len(ab)), hess=f_hess,
args=(ab['a'].values, ab['b'].values),
)
在添加其余的界限和约束之前;成本函数看起来像
def f(
x: np.ndarray, # (n*2)
a: np.ndarray, # n
b: np.ndarray, # n
) -> float:
x0, x1 = x.reshape((2, -1))
error = a*(x1 - x0**2)**2 + (b - x0)**2
return error.dot(error)
需要最后一个最小二乘语句(或类似的语句),因为 minimize
的成本函数必须返回一个标量。这意味着您需要重写
f_der
以便雅可比行列式反映每个误差项的平方,并且返回值形状为 2n。如果 - 正如您在边界中所示 - 0
least_squares<= x0 <= 1 and 1 <= x1 <= x2, then (x1 - x0**2) will always be positive, and b - x0 will also be positive, then (a) the self-dot-product above can be replaced with a simple sum, but perhaps more importantly, (b) you can migrate to 最小化向量而不是标量。这可能比 minimize
效果更好,但不幸的是,它可能无法解决此问题,因为它不支持约束。