如何用NumPy/SciPy计算数组映射函数的线性近似（雅可比矩阵）？

对于大规模非线性向量函数，手动推导雅可比矩阵既耗时又容易出错，SciPy提供了成熟的数值微分工具来自动近似计算雅可比矩阵，以下是两种常用方案：

1. 使用 scipy.optimize.approx_derivative #

这是SciPy中直接计算函数导数/雅可比矩阵的工具，支持多种差分方法（向前差分、中心差分等），可以直接处理向量输入、向量输出的函数。

结合你的示例函数，使用方式如下：

import numpy as np
from scipy.optimize import approx_derivative

def generalized_function(input_array):
    x,y,z = input_array
    output_array = [x**2+2*y, y**3-4*x+5*y, -z**2-2*x*y]
    return np.array(output_array)

# 定义要计算雅可比矩阵的输入点
input_point = np.array([1.0, 2.0, 3.0])

# 计算雅可比矩阵，使用中心差分提升精度
jacobian_matrix = approx_derivative(generalized_function, input_point, method='central')

print("雅可比矩阵：")
print(jacobian_matrix)

method='central' 指定使用中心差分，精度比默认的向前差分更高，但计算量会翻倍。你可以根据精度需求和计算成本选择合适的方法，也可通过epsilon参数手动调整差分步长。

2. 使用 scipy.optimize.jacobian #

这个函数会返回一个可调用对象，用于在任意输入点近似计算雅可比矩阵，适合需要多次计算不同点雅可比的场景：

import numpy as np
from scipy.optimize import jacobian

def generalized_function(input_array):
    x,y,z = input_array
    output_array = [x**2+2*y, y**3-4*x+5*y, -z**2-2*x*y]
    return np.array(output_array)

# 创建雅可比近似函数
jac_approx = jacobian(generalized_function)

# 在指定输入点计算雅可比
input_point = np.array([1.0, 2.0, 3.0])
jacobian_matrix = jac_approx(input_point)

print("雅可比矩阵：")
print(jacobian_matrix)

补充说明 #

• 数值微分的精度依赖于步长，默认步长通常能满足大多数场景需求，若需要更高精度可以手动调整epsilon参数。
• 如果你的函数可以用自动微分框架（如JAX）实现，自动微分的精度和效率会比数值微分更高，但这不属于NumPy/SciPy生态。

内容的提问来源于stack exchange，提问作者haricash