LK光流法与反向LK光流法

一、基本概念

光流法：光流法是一种基于灰度不变假设来跟踪图像中角点的方法
稀疏光流法：计算部分像素运动的称为稀疏光流法，以 Lucas-Kanade 光流为代表
稠密光流法：计算所有像素的称为稠密光流法，以 Hom-Schunck 光流为代表
灰度不变假设：同一个空间点的像素灰度值，在各个图像中是固定不变的。这是光流法的基本假设，同时也是一个很强的假设，在现实中很难得到保证；因此光流法得出的结果可能存在偏差
LK光流法：LK光流法基于灰度不变假设，来估计角点在两帧图像中的运动参数，然后根据运动参数寻找当前帧图像中角点在下一帧图像中的位置

二、2D中的LK光流法

2D中的LK光流法：已知$t$时刻图像中的角点$(x,y)$，求$dt$时间内该角点的运动（变换）参数$(dx,dy)$，从而得到$t+dt$时刻该角点在图像中的位置$(x+dx,y+dy)$。

1、空间点在图像中的灰度表示

假设空间中有一点
在$t$时刻该点在图像中的位置为$(x,y)$，记其此时该点在图像中的灰度为：
$$I(x,y,t)$$
设$t+dt$时刻该点在图像中的位置为$(x+dx,y+dy)$，则此时的该点在图像中的灰度为：
$$I(x+dx,y+dy,t+dt)$$
根据灰度不变假设，有式1：
$$I(x+dx,y+dy,t+dt)=I(x,y,t)——(1)$$

2、2D中的LK光流法推导

对式（1）左侧进行一阶泰勒展开：

代入式（1）可得式2：

对式2两边除以$dt$，移项可得式3:

对式3的理解,整理如下：

• $\partial I/\partial x$为$t$时刻图像在点$(x,y)$处的灰度对$x$的偏导，记作$I_x$
• $\partial I/\partial y$为$t$时刻图像在点$(x,y)$处的灰度对$y$的偏导，记作$I_y$
• $\partial I/\partial t$为$t$时刻图像在点$(x,y)$处的灰度对时间$t$的偏导，记作$I_t$
• $dx/dt$为$t$时刻图像位于$(x,y)$处的像素点在图像x方向上的运动速度，记作$u$
• $dy/dt$为$t$时刻图像位于$(x,y)$处的像素点在图像y方向上的运动速度，记作$v$

因此，式3可以记作矩阵形式式4：

在式4中，我们需要求取的变量就是像素的运动参数$u,v$，但由于这是一个二元一次方程，因此仅凭借一个像素是无法求解出$u,v$的，因此我们通常假设某一个窗口内的像素具有相同的运动，从而利用一个窗口的像素求取一组最优运动参数$u,v$，然后基于$u,v$计算下一时刻该像素在图像中的位置$(x+dx,y+dy)$

3、将2D光流法抽象成超定方程问题

考虑一个大小为 ω × ω 的窗口，它含有$w^2$数量的像素。该窗口内像素具有同样的运动，因此我们共有$w^2$个方程,如式5:

记：

则式5可以表示为式6，这是一个超定线性方程，没有精确解，但可以使用最小二乘求取该方程的解析解：

名词定义
超定线性方程：方程个数大于未知数个数的线性方程组

定理
线性方程组$Ax=b$当系数矩阵$A$超定时，最小二乘解为$x=-(A^{T}A{)}^{-1}{A}^{T}b$

因此由LK光流法求得$t$时刻图像点$(x,y)$的运动为$[u,v{]}^T=-(A^{T}A{)}^{-1}{A}^{T}b$

4、超定线性方程的最小二乘最优解定理证明

5、将2D光流法抽象为非线性优化问题

设函数$I(x,y)$为当前帧中点$(x,y)$的灰度，函数$T(x+dx,y+dy)$为目标帧中点$(x+dx,y+dy)$的灰度。那么根据灰度不变性假设，如果$dx,dy$为当前帧角点到目标帧角点的变换关系，则$dx,dy$应使得式7取得最小值：
$$||I(x,y)-T(x+dx,y+dy)|{|}^2$$
因此，我们可用通过求取最小化式7来求取当前帧角点与目标帧角点的变换关系，将抽象为优化目标函数如式8：
$$\underset{dx,dy}{\min }F(dx,dy)={\Vert f(dx,dy)\Vert }^2$$
其中$f(dx,dy)$表示灰度差异，相当于超定方程中的$\Delta I$：
$$f(dx,dy)=I(x,y)-T(x+dx,y+dy)$$
根据高斯牛顿法，可得$f(dx,dy)$的雅克比矩阵$J$为：
$$J=\left[\frac{\partial f}{\partial dx},\frac{\partial f}{\partial dy}\right]$$
注意这里对$f$求一阶偏导后没有转置，而有的高斯牛顿法会把$J$矩阵写成列向量，最后在结果上会相差一个转置
根据高斯牛顿法，可得当前帧点到目标帧点的变换参数增量$\Delta dx,\Delta dy$为式9：

$${\left[\Delta dx,\Delta dy\right]}^T=-{\left(J^T(dx,dy)J(dx,dy)\right)}^{-1}{J}^T(dx,dy)f(dx,dy)$$

对上式两端除以两帧之间的时间差$\Delta t$，并将$J(dx，dy)$简记为$J$，同时将$(x,y)$拓展为点$(x,y)$附近的点集，即可得到式10：

$${\left[\frac{\Delta dx}{\Delta t},\frac{\Delta dy}{\Delta t}\right]}^T=-{\left(J^{T}J\right)}^{-1}{J}^T\frac{f(dx,dy)}{\Delta t}$$

• 那么式10中$\frac{\Delta dx}{\Delta t},\frac{\Delta dy}{\Delta t}$与超定方程中的$[u,v]$含义相同，均是当前帧角点$(x,y)$的运动速度
• 式10中$\frac{f(dx,dy)}{\Delta t}=\frac{\Delta I}{\Delta t}=\frac{\partial I}{\partial t}=b$

即式10可以写成：
$${\left[u,v\right]}^T=-{\left(J^{T}J\right)}^{-1}{J}^{T}b$$
因此，由高斯牛顿法求取的角点在两帧图像间的运动（变换）参数与由超定线性方程求取的结果是一致的。

6、实践中的LK光流法（多层光流）

在实际应用中，我们通常使用高斯牛顿法来实现LK光流法：

1. 我们给定一个当前帧到目标帧的初始变换参数$(dx,dy)$，通常为0或单位阵
2. 根据该初始变换参数通过式9求取初始变换参数的增量$(\Delta dx,\Delta dy)$
3. 使用该增量更新初始变换参数
4. 重复步骤1-3，直到达到最大迭代次数，或者增量$(\Delta dx,\Delta dy)$小于阈值

在使用高斯牛顿等优化方法解决光流问题的时候，我们需要假设优化的初始值靠近最优质，才能保证算法的收敛。如果相机运动较快，两张图像差异较明显，那么单层图像光流法容易达到一个局部极小值。这种情况可以通过引人图像金字塔，通过多层光流改善算法效果。

图像金字塔：图像金字塔是指对同一个图像进行缩放，得到不同分辨率下的图像。以原始图像作为金字塔底层，每往上一层，就对下层图像进行一定倍率的缩放。
单层光流：只在原始图像上进行光流法
多层光流：从图像金字塔顶层依次向图像金字塔底层做光流法，每一层的结果作为下一层光流法的初始变换，如图所示：

多层光流的好处在于当原始图像的像素运动较大时，在金字塔顶层的图像看来，运动仍然在一个很小范围内，因此光流算法是收敛的，得到了一定精度内的变换参数。然后再以该参数为初始值，再次在下一层图像上调用光流法，得到精度更高的变换参数。因此多层光流保证了光流法对大幅度变换的收敛性。

三、光流法的应用拓展

光流法不仅是一种算法，更是一种思想，他能够用来解决一类问题：
设有一组数据$A$，经过一个变换$p()$，得到了一组数据$B$，函数$W()$可以提取该类型数据的一个特征，且数据$p(A)$和数据$B$的该特征是相同的，即满足：
$$W(B)=W(p(A))$$
那么可以构造优化目标函数，通过最小化特征差异求取最优变换$p$:
$$\underset{p}{\min }F(p)={\Vert f(p)\Vert }^2,f(p)=W(B)-W(p(A))$$
每次迭代时，最优变换$p$的增量$\Delta p$为，$J$为$f(p)$的雅克比矩阵：
$$\Delta p=-{\left(J^{T}J\right)}^{-1}{J}^{T}f(p)$$

四、逆向光流法(inverse compositional)

在光流法中，由于每次迭代过程中特征差$f(x)$都会发生改变，因此每次迭代过程中都要重新计算特征差$f(x)$的雅克比矩阵$J$,而逆向光流法则可以解决这个问题

1、逆向光流法思想

在每次迭代前，已知数据$A$当前变换为$p$，我们不再求取数据$A$当前变换$p$的增量，而是求取一个针对数据$B$的变换量$q$，使得变换后的$q(B)$与当前的$p(A)$特征差最小，即优化目标函数为：
$$f(q)=W(q(B))-W(p(A))$$
其中，对于每一次迭代，$q$的初始值为0，即没有变换；$p$则是一个与$q$无关的变换，在$f(q)$中为一个常量

2、逆向光流法推导

根据高斯牛顿算法，求取$f(q)$的雅克比矩阵，式中$q$为初始值，如式4-1：

$$J(q)=\frac{\partial f(q)}{\partial q}=\frac{\partial (W(q(B))-W(p(A)))}{\partial q}=\frac{\partial W(q(B))}{\partial q}-\frac{\partial W(p(A))}{\partial q}$$
对该式子有如下简化：

• 由于$q$的初始值为0，因此式中$q(B)=B$，因此$\frac{\partial W(q(B))}{\partial q}=\frac{\partial W(B)}{\partial q}$
• 由于$W(),p(),A$均为与$q$无关的常量，因此$\frac{\partial W(p(A))}{\partial q}=0$

综上，式子4-1可以化简为式4-2：
$$J(q)=\frac{\partial f(q)}{\partial q}=\frac{\partial W(B)}{\partial q}$$
即此时对特征差求取关于$q$的雅克比矩阵就等于对目标数据$B$求关于变换参数$q$的雅克比矩阵
根据高斯牛顿法，可知当前当前迭代的$q$变换增量为式4-3：
$$\Delta q=-{\left(J^{T}J\right)}^{-1}{J}^{T}f(q)$$
其中，式4-3中$q$为初始值0，因此式4-3中的$f(q)$为：
$$f(q)=W(q(B))-W(p(A))=W(B)-W(p(A))$$

3、逆向光流法迭代更新

在逆向光流法推导中，我们获得了针对$B$的变换增量$\Delta q$，我们需要把这个变换增量的逆变换$(\Delta q{)}^{-1}$更新到针对$A$的变换$p$中；而因为针对$B$的变换$q$的初值为0，且在本次迭代中没有更新，因此下一次迭代开始时候，针对$B$变换$q$的值认为0，因此式子4-1仍可以化简位式4-2；且数据$B$也没有变换，因此下次迭代时候的雅克比矩阵$J$与本次迭代时的雅克比矩阵$J$相同，即在迭代过程中只需要计算一次雅克比矩阵$J$。

这里包含了一个逆向光流法迭代更新的前提条件：数据$A$与数据$B$间的变换时可逆的，这样才能 使用目标数据$B$的变换增量的逆变换更新源数据$A$