从局部导数看场 · 第 1 章
雅可比与梯度:导数是局部线性映射
把一元导数推广到多元,导数就从一个数变成一个矩阵,即雅可比;标量场的雅可比转置后是梯度。本章推导梯度为何指向最快上升方向、为何垂直于等值面,并给出相机投影和位姿优化中的雅可比。
系列目录第 1 / 3 章
问题是什么
多元函数的”导数”到底是一个什么样的对象?为什么标量场的导数是一个向量,而向量场的导数是一个矩阵?
这个问题看起来只是记号上的区别,却决定了后面所有算子的位置。散度、旋度、海森、拉普拉斯都是从同一个对象,即雅可比矩阵 (Jacobian Matrix),出发得到的。本章先把这个对象讲清楚。
直觉怎么理解
一元函数的导数回答的是:输入改变一点点,输出会改变多少。在一点附近放大函数图像,曲线看起来像一条直线,导数就是这条直线的斜率。
多元情形完全相同,只是”斜率”换成了矩阵。把一个映射 放大到某一点附近,弯曲的映射看起来就像一个线性变换:一个小方格被映射成一个小平行四边形。描述这个线性变换的矩阵,就是该点的雅可比。
雅可比的每一列回答”沿某一个输入方向移动时,所有输出怎样变化”;每一行回答”某一个输出分量对所有输入方向的敏感程度”。
标量场 只有一个输出,雅可比只有一行。把这一行立起来写成列向量,就是梯度 (Gradient)。
数学 / 几何表达
从一元导数到雅可比
一元函数在 附近的一阶展开为
对 ,如果存在矩阵 使
就称 在 处可微, 为它的雅可比。取 并令 ,可以读出每个元素:
雅可比的形状是”输出维数 × 输入维数”,不一定是方阵。只有 时才讨论它的行列式,即局部体积缩放比例。
链式法则就是矩阵乘法
设 ,。对 的小变化 ,有 ,再有 。于是
复合映射的局部线性近似,就是各段局部线性近似的乘积。机器人中的长链式求导(例如”位姿 → 三维点 → 像素 → 残差”)都是这个公式的反复应用。
梯度:标量场雅可比的转置
对标量场 ,雅可比是一个 的行向量。定义梯度为它的转置:
于是一阶展开可以写成 。沿单位方向 的方向导数 (Directional Derivative) 为
梯度为何指向最快上升方向
由 Cauchy–Schwarz 不等式(见内积、协方差与相关系数),对单位向量 有
等号成立当且仅当 与 同向。所以在 的点,沿 方向上升最快,最大变化率为 ;沿相反方向下降最快。这就是梯度下降的依据。
例如 ,有 。在点 ,最快上升方向是 的单位向量,而不是位置向量 的方向。各方向的”陡峭程度”不同,梯度会偏向更陡的坐标。
梯度垂直于等值面
设 是等值面 上过 的任意光滑曲线,。由于 恒成立,对 求导,并用链式法则:
可以是等值面在 处的任意切向量,所以 与整个切平面正交,即梯度是等值面的法向量。
导数与梯度之间隔着一个内积
严格地说,导数 是一个线性函数:它接收一个方向 ,输出变化率 。梯度则是借助内积,把这个线性函数表示成一个向量:对所有 ,。
在笛卡尔坐标和标准内积下,两者的分量完全相同,所以这一步很容易被忽略。换成一般的内积 ( 对称正定)时,梯度变为 ,与偏导数组成的向量不再相同。曲线坐标和流形上的梯度公式,都源于这一点。
机器人中的例子
相机投影的雅可比
对不含畸变的针孔相机,相机坐标系中的点 ()投影到像素
逐项求偏导:,,; 同理。于是
这是一个 矩阵,它预测三维点的微小移动会引起多大的像素位移。有两点可以直接读出来:
- 所有元素都含 ,同样的三维位移,近处的点在图像上移动得更多;
- 第三列描述沿光轴方向的移动,它使点沿着从主点出发的射线方向在图像上移动。第 2 章推导光流散度时会再用到这一点。
隐式表面的法线
如果表面由隐式函数 定义,由”梯度垂直于等值面”可知,在 的光滑点,单位法线为
对精确的欧氏符号距离函数 (Signed Distance Function, SDF),在光滑区域内有 (程函方程,Eikonal Equation),所以梯度本身就是单位法线。
需要注意的是,RGB-D 融合得到的截断符号距离函数 (TSDF) 通常只是沿视线方向的投影距离近似,还带有噪声、插值误差和截断饱和。它不能在全域当作精确的 SDF 使用,用梯度估计法线时,应只在表面附近、观测充分的体素中进行。
位姿优化中的梯度
设残差向量 依赖于位姿的局部参数 ,目标函数为 。由链式法则,
是 矩阵,梯度是 向量。Gauss–Newton 法用 近似海森,每一步求解 ;第 3 章会说明这个近似从何而来。
这里的求导变量是位姿的局部参数,而不是三维世界坐标。六维位姿参数的雅可比还取决于扰动是左乘还是右乘,不同约定下的公式不能混用。ICP、RGB-D 里程计和光束法平差 (Bundle Adjustment) 中的 Jacobian 和 gradient,都属于这一套语言。
可视化 / 代码
下面的代码依次验证四件事:用 SymPy 求投影雅可比;与有限差分结果对照;最小二乘目标的梯度等于 ;SDF 的梯度范数为 1。
import numpy as np
import sympy as sp
# 1. 针孔投影的雅可比(符号推导)
X, Y, Z, fx, fy, cx, cy = sp.symbols("X Y Z f_x f_y c_x c_y", positive=True)
proj = sp.Matrix([fx * X / Z + cx, fy * Y / Z + cy])
J_sym = proj.jacobian([X, Y, Z])
sp.pprint(J_sym)
# 2. 与中心差分对照
J_fn = sp.lambdify((X, Y, Z, fx, fy), J_sym)
p_fn = sp.lambdify((X, Y, Z, fx, fy, cx, cy), proj)
P0 = np.array([0.4, -0.2, 3.0])
intr = (500.0, 500.0, 320.0, 240.0)
eps = 1e-6
J_num = np.zeros((2, 3))
for j in range(3):
d = np.zeros(3); d[j] = eps
J_num[:, j] = ((p_fn(*(P0 + d), *intr) - p_fn(*(P0 - d), *intr)) / (2 * eps)).ravel()
print(np.abs(J_num - np.array(J_fn(*P0, 500.0, 500.0))).max()) # 约 1e-8
# 3. 最小二乘:grad E = J^T r,R^6 -> R^100
rng = np.random.default_rng(0)
A = rng.standard_normal((100, 6))
b = rng.standard_normal(100)
r = lambda xi: np.sin(A @ xi) - b
J_r = lambda xi: np.cos(A @ xi)[:, None] * A
E = lambda xi: 0.5 * r(xi) @ r(xi)
xi0 = rng.standard_normal(6)
g_analytic = J_r(xi0).T @ r(xi0)
g_numeric = np.array([(E(xi0 + eps * e) - E(xi0 - eps * e)) / (2 * eps) for e in np.eye(6)])
print(J_r(xi0).shape, g_analytic.shape, np.abs(g_analytic - g_numeric).max())
# (100, 6) (6,) 约 1e-8
# 4. 单位球的 SDF:f(p) = |p| - 1,梯度 p/|p| 的范数恒为 1
p = rng.standard_normal((5, 3))
grad = p / np.linalg.norm(p, axis=1, keepdims=True)
print(np.linalg.norm(grad, axis=1)) # [1. 1. 1. 1. 1.]
参考资料
- Hartley, R., & Zisserman, A. (2004). Multiple View Geometry in Computer Vision (2nd ed.). Cambridge University Press.
- Barfoot, T. D. (2017). State Estimation for Robotics. Cambridge University Press.
- Strang, G. (2016). Introduction to Linear Algebra (5th ed.). Wellesley-Cambridge Press.
- OpenStax. Calculus Volume 3, Chapter 4 (Differentiation of Functions of Several Variables). https://openstax.org/details/books/calculus-volume-3
- Khan Academy. Derivatives of multivariable functions(含雅可比矩阵与梯度的视频)。https://www.khanacademy.org/math/multivariable-calculus/multivariable-derivatives