Ctrl+K
scipy.differentiate.

jacobian#

scipy.differentiate.jacobian(f, x, *, tolerances=None, maxiter=10, order=8, initial_step=0.5, step_factor=2.0, step_direction=0)[源代码]#

数值计算函数的雅可比矩阵。

参数:
f可调用对象

需要计算雅可比矩阵的函数。其签名必须为

f(xi: ndarray) -> ndarray

其中 xi 的每个元素都是有限实数。如果要微分的函数接受额外的参数，请包装它（例如，使用 functools.partial 或 lambda），并将包装后的可调用对象传递给 jacobian。f 不得改变数组 xi。有关矢量化以及输入和输出的维度，请参见“说明”。

xfloat array_like

计算雅可比矩阵的点。必须至少有一维。有关维度和矢量化，请参见“说明”。

tolerances浮点数字典，可选

绝对和相对公差。字典的有效键为

  • atol - 导数的绝对公差

  • rtol - 导数的相对公差

当 res.error < atol + rtol * abs(res.df) 时，迭代将停止。默认的 atol 是适当数据类型的最小正规数，默认的 rtol 是适当数据类型的精度的平方根。

maxiterint，默认值：10

要执行的算法的最大迭代次数。参见“说明”。

orderint，默认值：8

要使用的有限差分公式的（正整数）阶数。奇数将向上舍入到下一个偶数。

initial_stepfloat array_like，默认值：0.5

有限差分导数逼近的（绝对）初始步长。必须可与 x 和 step_direction 进行广播。

step_factorfloat，默认值：2.0

每次迭代中步长减小的因子；即，迭代 1 中的步长为 initial_step/step_factor。如果 step_factor < 1，则后续步骤将大于初始步骤；如果不需要小于某个阈值的步骤（例如，由于减法抵消错误），这可能会很有用。

step_direction整数 array_like

表示有限差分步长方向的数组（例如，当 x 位于函数域的边界附近时使用）。必须可与 x 和 initial_step 进行广播。其中 0（默认值），使用中心差分；其中为负数（例如，-1），步长为非正数；其中为正数（例如，1），所有步长均为非负数。

返回:
res_RichResult

类似于 scipy.optimize.OptimizeResult 实例的对象，具有以下属性。描述的编写方式就好像这些值是标量一样；但是，如果 f 返回一个数组，则输出将是相同形状的数组。

successbool 数组

算法成功终止的位置（状态 0）；否则为 False。

statusint 数组

表示算法退出状态的整数。

  • 0：算法收敛到指定的公差。

  • -1：误差估计值增加，因此终止了迭代。

  • -2：达到最大迭代次数。

  • -3：遇到非有限值。

dffloat 数组

如果算法成功终止，则为 f 在 x 处的雅可比矩阵。

errorfloat 数组

误差估计值：雅可比矩阵的当前估计值与上次迭代中的估计值之间的差的幅度。

nitint 数组

执行的算法迭代次数。

nfevint 数组

评估 f 的点的数量。

属性的每个元素都与 df 的相应元素关联。例如，nfev 的元素 i 是为计算 df 的元素 i 而评估 f 的点的数量。

另请参见

derivative，hessian

说明

假设我们希望计算函数 \(f: \mathbf{R}^m \rightarrow \mathbf{R}^n\) 的雅可比矩阵。将 \(m\) 和 \(n\) 的正整数值分别赋值给变量 m 和 n，并让 ... 表示任意整数元组。如果希望在单个点计算雅可比矩阵，则

  • 参数 x 必须是形状为 (m,) 的数组

  • 必须对参数 f 进行矢量化，以接受形状为 (m, ...) 的数组。第一个轴表示 \(f\) 的 \(m\) 个输入；其余部分用于在单个调用中评估多个点的函数。

  • 参数 f 必须返回形状为 (n, ...) 的数组。第一个轴表示 \(f\) 的 \(n\) 个输出；其余部分用于评估多个点的函数的结果。

  • 结果对象的属性 df 将是形状为 (n, m) 的数组，即雅可比矩阵。

此函数也是矢量化的，这意味着可以在单个调用中计算 k 个点的雅可比矩阵。在这种情况下，x 将是形状为 (m, k) 的数组，f 将接受形状为 (m, k, ...) 的数组并返回形状为 (n, k, ...) 的数组，并且结果的 df 属性将具有形状 (n, m, k)。

假设所需的被调用对象 f_not_vectorized 没有向量化；它只能接受形状为 (m,) 的数组。满足所需接口的一个简单解决方案是将 f_not_vectorized 包装如下

def f(x):
    return np.apply_along_axis(f_not_vectorized, axis=0, arr=x)

或者，假设所需的被调用对象 f_vec_q 是向量化的，但仅适用于形状为 (m, q) 的二维数组。为了满足所需的接口，请考虑

def f(x):
    m, batch = x.shape[0], x.shape[1:]  # x.shape is (m, ...)
    x = np.reshape(x, (m, -1))  # `-1` is short for q = prod(batch)
    res = f_vec_q(x)  # pass shape (m, q) to function
    n = res.shape[0]
    return np.reshape(res, (n,) + batch)  # return shape (n, ...)

然后将包装后的被调用对象 f 作为 jacobian 的第一个参数传递。

参考资料

[1]

雅可比矩阵和行列式, 维基百科, https://en.wikipedia.org/wiki/Jacobian_matrix_and_determinant

示例

Rosenbrock 函数将 \(\mathbf{R}^m \rightarrow \mathbf{R}\) 映射；SciPy 的实现 scipy.optimize.rosen 向量化后可以接受形状为 (m, p) 的数组，并返回形状为 p 的数组。假设我们希望在 [0.5, 0.5, 0.5] 处评估雅可比矩阵（又名梯度，因为该函数返回一个标量）。

>>> import numpy as np
>>> from scipy.differentiate import jacobian
>>> from scipy.optimize import rosen, rosen_der
>>> m = 3
>>> x = np.full(m, 0.5)
>>> res = jacobian(rosen, x)
>>> ref = rosen_der(x)  # reference value of the gradient
>>> res.df, ref
(array([-51.,  -1.,  50.]), array([-51.,  -1.,  50.]))

作为具有多个输出的函数的示例，请考虑 [1] 中的示例 4。

>>> def f(x):
...     x1, x2, x3 = x
...     return [x1, 5*x3, 4*x2**2 - 2*x3, x3*np.sin(x1)]

真正的雅可比矩阵由下式给出

>>> def df(x):
...         x1, x2, x3 = x
...         one = np.ones_like(x1)
...         return [[one, 0*one, 0*one],
...                 [0*one, 0*one, 5*one],
...                 [0*one, 8*x2, -2*one],
...                 [x3*np.cos(x1), 0*one, np.sin(x1)]]

在任意点评估雅可比矩阵。

>>> rng = np.random.default_rng()
>>> x = rng.random(size=3)
>>> res = jacobian(f, x)
>>> ref = df(x)
>>> res.df.shape == (4, 3)
True
>>> np.allclose(res.df, ref)
True

在一次调用中评估 10 个任意点的雅可比矩阵。

>>> x = rng.random(size=(3, 10))
>>> res = jacobian(f, x)
>>> ref = df(x)
>>> res.df.shape == (4, 3, 10)
True
>>> np.allclose(res.df, ref)
True
在本页上