Jacobian和Hessian在`scipy.optimize.minimize`中输入

Med*_*ata 7 python optimization lambda numpy scipy

我试图了解"狗腿"方法在Python scipy.optimize.minimize函数中的工作原理.我正在调整帮助页面底部的示例.

根据笔记,狗腿洞方法需要Jacobian和Hessian参数.为此,我使用numdifftools包:

import numpy as np
from scipy.optimize import minimize
from numdifftools import Jacobian, Hessian

def fun(x,a):
    return (x[0] - 1)**2 + (x[1] - a)**2

x0 = np.array([2,0]) # initial guess
a = 2.5

res = minimize(fun, x0, args=(a), method='dogleg',
               jac=Jacobian(fun)([2,0]), hess=Hessian(fun)([2,0]))

print(res)
Run Code Online (Sandbox Code Playgroud)

编辑:

如果我按照以下帖子的建议进行更改,

res = minimize(fun, x0, args=a, method='dogleg',
               jac=Jacobian(lambda x: fun(x,a)),
               hess=Hessian(lambda x: fun(x,a)))
Run Code Online (Sandbox Code Playgroud)

我收到一个错误TypeError: <lambda>() takes 1 positional argument but 2 were given.我究竟做错了什么?

在初始猜测中计算雅可比行列式和Hessian也是正确的x0吗?

bna*_*aul 7

该错误是来自和的调用,JacobianHessian不是来自的调用minimize。替换为Jacobian(fun)Jacobian(lambda x: fun(x, a))并类似地可以解决问题Hessian(因为现在要区分的函数只有一个向量参数)。

另一件事:(a)a,如果您希望将其用作元组使用(a,)

import numpy as np
from scipy.optimize import minimize
from numdifftools import Jacobian, Hessian

def fun(x, a):
    return (x[0] - 1) **2 + (x[1] - a) **2

def fun_der(x, a):
    return Jacobian(lambda x: fun(x, a))(x).ravel()

def fun_hess(x, a):
    return Hessian(lambda x: fun(x, a))(x)

x0 = np.array([2, 0]) # initial guess
a = 2.5

res = minimize(fun, x0, args=(a,), method='dogleg', jac=fun_der, hess=fun_hess)
print(res)
Run Code Online (Sandbox Code Playgroud)


Grr*_*Grr 6

我知道这是一个玩具示例,但我想指出,使用像JacobianHessian计算衍生工具而不是推导函数本身的工具是相当昂贵的.例如,使用您的方法:

x0 = np.array([2, 0])
a = 2.5
%timeit minimize(fun, x0, args=(a,), method='dogleg', jac=fun_der, hess=fun_hess)
100 loops, best of 3: 13.6 ms per loop
Run Code Online (Sandbox Code Playgroud)

但你可以计算衍生函数:

def fun_der(x, a):
    dx = 2 * (x[0] - 1)
    dy = 2 * (x[1] - a)
    return np.array([dx, dy]

def fun_hess(x, a):
    dx = 2
    dy = 2
    return np.diag([dx, dy])

%timeit minimize(fun, x0, args=(a,), method='dogleg', jac=fun_der, hess=fun_hess)
1000 loops, best of 3: 279 µs per loop
Run Code Online (Sandbox Code Playgroud)

如你所见,它快了近50倍.它真正开始加入复杂的功能.因此,我总是尝试自己明确地推导出这些函数,无论它有多么困难.一个有趣的例子是基于内核的Inductive Matrix Completion实现.

argmin --> sum((A - gamma_u(X) Z gamma_v(Y))**2 - lambda * ||Z||**2)
where gamma_u = (1/sqrt(m_x)) * [cos(UX), sin(UX)] and
gamma_v = (1/sqrt(m_y)) * [cos(VY), sin(VY)]
X.shape = n_x, p; Y.shape = n_y, q; U.shape = m_x, p; V.shape = m_y, q; Z.shape = 2m_x, 2m_y
Run Code Online (Sandbox Code Playgroud)

与明确推导和利用这些函数相比,从该等式计算梯度和粗麻布是非常不合理的.因此,正如@bnaul所指出的,如果你的函数确实有封闭形式派生,你真的想要计算和使用它们.