CSK与KCF算法推导（三）

  本文是CSK与KCF算法推导的第三篇，主要介绍算法的主体部分推导。

• CSK与KCF算法推导（一）——DFT、相关运算和循环卷积
• CSK与KCF算法推导（二）——带核函数的岭回归问题求解
• （本文）CSK与KCF算法推导（三）——算法主体部分推导
• CSK与KCF算法推导（四）——从一维推广到二维
• CSK与KCF算法推导（五）——具体实现细节
• CSK与KCF算法推导（六）——HOG特征提取与多通道计算

循环矩阵

  循环矩阵从CSK的文章里就提出了，概念也很容易理解。一个样本$x$，可以通过循环移位生成不同的向量，组合成一个矩阵$C(x)$。下面的图里分别列出了$x$、$C(x)$和$C(x)$的转置。

  循环矩阵$C(x)$和另一个列向量$z$相乘，得到的结果就是$x$和$z$进行相关运算的结果。这个也很好理解，就是$x$在$z$上面滑动，然后进行乘加运算。同样可以写出循环卷积的矩阵表达式，循环卷积需要先将$x$倒序，然后再滑动。对比循环卷积所用的矩阵和$(C(x){)}^T$的颜色，可以发现它们俩是一样的！

  循环矩阵$C(x)$既然能和循环卷积联系起来了，那么自然可以想到把DFT矩阵也联系起来。这就能得出一个重要的结论：循环矩阵$C(x)$能够用DFT矩阵对角化！
  下面是证明过程。先根据时域卷积等于频域乘积写出等式：
$${\left(C(x)\right)}^{T}z=x\ast z={\mathcal{F}}^{-1}(\hat{x}\odot \hat{z})$$
  其中“$\odot$”是哈达玛积（Hadamard product），或者称为“element-wise product”，表示两个向量对应位置的分量直接相乘。$\hat{x}$和$\hat{z}$分别表示$x$和$z$DFT变换之后的序列。
$${\left(C(x)\right)}^{T}z=\frac{1}{s}{F}^H(\hat{x}\odot Fz)=\frac{1}{s}{F}^{H}diag(\hat{x})Fz$$
  两个列向量的哈达玛积可以转换成对角矩阵和列向量相乘。$diag(\hat{x})$表示把$\hat{x}$放到对角线上形成对角矩阵。因为上式对所有的$z$都成立，所以可以去掉$z$，得到：
$$\begin{array}{l}{\left(C(x)\right)}^T=\frac{1}{s}{F}^{H}diag(\hat{x})F\\ C(x)=\left(\frac{1}{\sqrt{s}}{F}^T\right)diag(\hat{x})\left(\frac{1}{\sqrt{s}}{F}^{\ast }\right)\end{array}$$
  又因为$F=F^T$，且$U=\frac{1}{\sqrt{s}}F$是酉矩阵。$C(x)$可以进一步表示为：
$$C(x)=Udiag(\hat{x})U^{\ast }$$
  这个结论就非常厉害了。所有循环矩阵都可以特征分解，特征值是原始序列的DFT；所有循环矩阵的特征向量都相同。这就为后面的运算带来了很大的便利。
  上面的推导是从循环卷积的角度出发的，而我在第一篇文章中还写了相关运算与DFT的关系。

  这里的相关运算相当于是$x$相对于$z$在向右滑动，所以$x$的DFT结果需要取复共轭，可以这样表示结果：
$$C(x)z=x\star z={\mathcal{F}}^{-1}({\hat{x}}^{\ast }\odot \hat{z})$$
  其中$\star$表示相关运算，经过类似的推导之后可以得到：
$$C(x)=U^{\ast }diag({\hat{x}}^{\ast })U$$

循环移位的意义

  一方面原因已经在上面推导，循环矩阵能够用DFT对角化。另一方面要清楚输入样本$x$并不都是我们关心的目标。我们关心的Target只有其中的一部分，就如下面的图所示。

  这么做就能够在训练的时候为我们关心的Target提供一些上下文（context）信息。$x$在循环移位的过程中确实起到了“密集采样”的作用。

生成标签

  在生成标签的时候我们可以给原始样本$x$一个较大的标签值，其他由循环移位产生的样本给较小的标签值，相对位移越大，标签值越小。

  因为是循环移位，所以移位越接近一个周期，这个样本越接近正样本。所以样本的标签值是一个峰值在两侧的高斯函数。

训练（不带核函数）

  样本训练还是用的岭回归。结合第二篇文章，虽然在引入核函数之后，我们就不关心$\omega$是多少了，但在没有引入核函数的时候，能求出$\omega$自然更好。
$$\omega ={\left(X^{T}X+\lambda I\right)}^{-1}{X}^{T}y$$
  只需要将样本矩阵$X$替换成$C(x)$即可。上面的表达式里用的是转置而不是共轭转置，因为我们前面的推导是在实数的条件下得到的，因为这对于我们分析问题来说已经足够了。
  有人可能会担心，$C(x)$分解之后的三个矩阵都是复矩阵，这里是不是应该用共轭转置啊？要用共轭转置很简单，把前面岭回归的推导用共轭转置推导一遍就行，这里的转置是可以推广到共轭转置的。但其实也可以就只用转置，因为这里的$C(x)$本质上还是实数矩阵，转置是能够成立的。我们已经在前面已经有两个关于$C(x)$的分解了，可以都用来帮助化简。把$X=C(x)$代入：
$$\omega ={\left(U^{T}diag({\hat{x}}^{\ast })U^{H}Udiag(\hat{x})U^{\ast }+\lambda I\right)}^{-1}{U}^{T}diag({\hat{x}}^{\ast })U^{H}y$$
  把握几个关键点，$U^{H}U=U{U}^H=U^{\ast }U=U{U}^{\ast }=I$；对角矩阵相乘就是对角线上元素的对应相乘；对角矩阵转置就是对角线上元素求倒数；$(ABC{)}^{-1}=C^{-1}{B}^{-1}{C}^{-1}$；$Ux=\hat{x}/\sqrt{s}$；$U^{\ast }x={\hat{x}}^{\ast }/\sqrt{s}$。
$$\begin{array}{l}\omega ={\left(U^{T}diag({\hat{x}}^{\ast })U^{H}Udiag(\hat{x})U^{\ast }+\lambda I\right)}^{-1}{U}^{T}diag({\hat{x}}^{\ast })U^{H}y\\ =Udiag\left(\frac{1}{{\hat{x}}^{\ast }\odot \hat{x}+\lambda }\right)U^{\ast }Udiag({\hat{x}}^{\ast })U^{\ast }y\\ =Udiag\left(\frac{{\hat{x}}^{\ast }}{{\hat{x}}^{\ast }\odot \hat{x}+\lambda }\right)U^{\ast }y\end{array}$$
  这里$diag$里面的计算都是element-wise的。
$$U^{\ast }\omega =diag\left(\frac{{\hat{x}}^{\ast }}{{\hat{x}}^{\ast }\odot \hat{x}+\lambda }\right)U^{\ast }y$$
$${\hat{\omega }}^{\ast }=\frac{{\hat{x}}^{\ast }\odot {\hat{y}}^{\ast }}{{\hat{x}}^{\ast }\odot \hat{x}+\lambda }$$
  最后的化简结果如下。相比于原先复杂的矩阵运算，现在只需要少量DFT、IDFT和element-wise的运算即可。
$$\hat{\omega }=\frac{\hat{x}\odot \hat{y}}{{\hat{x}}^{\ast }\odot \hat{x}+\lambda }$$
  这个地方的推导也是知乎上有人提出推导有误的地方。原文中的结果里面分子的$\hat{x}$多了复共轭，因为在原文的公式（55）处落下了一个复共轭。而这部分结果作者也没有用到具体算法中，所以对算法没有影响。

训练（带核函数）

  带核函数的情况下，我们就不用求$\omega$了，只要求中间变量$\alpha$就好了。
$$\alpha ={\left(K+\lambda I\right)}^{-1}y$$
  计算Gram矩阵

  逐行观察上面的矩阵乘法，就看前两行就行。可以发现这两行的结果也是一个相对移位的关系。所以当样本矩阵是循环矩阵的时候，这一组样本的Gram矩阵也会是循环矩阵。而循环的基向量就是Gram矩阵的第一行。$K$是高维空间的Gram矩阵，这一结论自然也成立，$K$也是一个循环矩阵。
$$K=C(k^{xx})k^{xx}=\left[\begin{array}{ccccc}\kappa (x,P^{0}x)& \kappa (x,P^{1}x)& \kappa (x,P^{2}x)& \cdots & \kappa (x,P^{n}x)\end{array}\right]$$
  这里的$P$矩阵的作用是让列向量$\phi (x)$循环移位。比如$P^0$是不移位，$P^1$是循环移位1次，$P^2$是循环移位2次，以此类推。
$$\alpha ={\left(C(k^{xx})+\lambda I\right)}^{-1}y={\left(Udiag({\hat{k}}^{xx})U^{\ast }+\lambda I\right)}^{-1}y$$
$$U^{\ast }\alpha =diag\left(\frac{1}{{\hat{k}}^{xx}+\lambda }\right)U^{\ast }y$$
$${\hat{\alpha }}^{\ast }=\frac{{\hat{y}}^{\ast }}{{\hat{k}}^{xx}+\lambda }$$
  用循环矩阵的形式表示$K$，然后特征分解，利用对角矩阵和酉矩阵的特性即可完成上面的化简。这时候计算量相对比较大的地方就是求$k^{xx}$了，它是向量$x$和$x$的循环移位的向量在高维空间的内积构成的向量。
  然后对上面的等式两边取复共轭，就可以得到：
$$\hat{\alpha }=\frac{\hat{y}}{{{\hat{k}}^{xx}}^{\ast }+\lambda }$$
  但这和论文里的结论不一样，论文里的结论，分母的${\hat{k}}^{xx}$没有取复共轭。但实际上这里取复共轭或者不取复共轭都行，因为${\hat{k}}^{xx}$是实数。
  什么样的实数序列做DFT之后还能是实数？观察旋转因子$W_N^{kn}$，旋转因子在复平面上的分布就是有对称性的，当实数序列满足一定条件之后，虚部就能够抵消。

  观察可以发现，只要序列满足$x[i]=x[N-i]i=1,2,...,N-1$这个条件，这些旋转因子的虚部就都能够被抵消。对于二维的图像来说，就是我们不关心它的第0行和第0列，第0行和第0列可以是任意实数，然后剩余的数据是中心对称的就可以。
  序列$x$的自相关$k^{xx}$自然就满足这一条件。所以更一般的结果还可以是：
$$\hat{\alpha }=\frac{\hat{y}}{real({\hat{k}}^{xx})+\lambda }$$
  我们还可以让结果更简单，就是如果$y$给定的标签同样满足上面的条件，那么$\hat{y}$就也是一个实数。这样整个$\hat{\alpha }$也是实数了！！
$$\hat{\alpha }=\frac{real(\hat{y})}{real({\hat{k}}^{xx})+\lambda }$$
  通过上面的计算，我们还可以想到一个结论，$K=C(k^{xx})$不光是一个循环矩阵，也还是一个对称矩阵，$K=K^T$。结合这一结论我们可以重新进行一次推导：
$$\alpha ={\left(K^T+\lambda I\right)}^{-1}y$$
$$\alpha ={\left(C(k^{xx}{)}^T+\lambda I\right)}^{-1}y={\left(U^{\ast }diag({\hat{k}}^{xx})U+\lambda I\right)}^{-1}y$$
$$U\alpha =diag\left(\frac{1}{{\hat{k}}^{xx}+\lambda }\right)Uy$$
$$\hat{\alpha }=\frac{\hat{y}}{{\hat{k}}^{xx}+\lambda }$$
  这样就能直接得到论文中给出的结论。

带核函数的检测

  检测就是我们的训练完成之后送进一个新的样本$z$来预测输出。不带核函数的情况下因为直接计算出了$\omega$，所以只需要求$\omega$和$z$的内积就能得到结果。带核函数的情况在上一篇文章中介绍了，不需要显式地求出$\omega$，只需要将输入样本$z$和所有训练样本做内积，得到一个维度是样本个数的列向量，然后把这个列向量和$\alpha$做内积即可得到一个预测结果。
  逐个样本$z$进行检测太慢了，我们把$z$也扩展成循环矩阵$C(z)$，然后就可以直接做一批样本的检测，这一批样本相互之间是空间上的位移关系，最后得到$f(z)$是一个列向量，表示一批样本的检测结果。

  在没有引入核函数的情况下，上面的$C(z)$和$(C(x){)}^T$相乘得到的也是一个循环矩阵。引入核函数之后，每个样本都经过$\phi (x)$映射到高维空间。我们把高维空间的$C^{\prime }(z)$和$(C^{\prime }(x){)}^T$相乘的结果记作$K^z$
$$K^z=C(k^{zx})k^{zx}=\left[\begin{array}{ccccc}\kappa (z,P^{0}x)& \kappa (z,P^{1}x)& \kappa (z,P^{2}x)& \cdots & \kappa (z,P^{n}x)\end{array}\right]$$
  把结果用循环矩阵化简：
$$f(z)=C(z){\left(C(x)\right)}^T\alpha =K^z\alpha =C(k^{zx})\alpha =Udiag({\hat{k}}^{zx})U^{\ast }\alpha$$
$$U^{\ast }f(z)=diag({\hat{k}}^{zx})U^{\ast }\alpha$$
$${\hat{f}}^{\ast }(z)={\hat{k}}^{zx}\odot {\hat{a}}^{\ast }$$
  这里的结果和论文中的结果有点不一样，这是因为$K^z$的表示方式不太一样，原文中用的是$f(z)=(K^z{)}^T\alpha$，所以$K^z=C(x){\left(C(z)\right)}^T$

$$K^z=C(k^{xz})k^{xz}=\left[\begin{array}{ccccc}\kappa (x,P^{0}z)& \kappa (x,P^{1}z)& \kappa (x,P^{2}z)& \cdots & \kappa (x,P^{n}z)\end{array}\right]$$
  把结果用循环矩阵化简：
$$f(z)=(K^z{)}^T\alpha =(C(k^{xz}){)}^T\alpha =U^{\ast }diag({\hat{k}}^{xz})U\alpha$$
$$Uf(z)=diag({\hat{k}}^{xz})U\alpha$$
$$\hat{f}(z)={\hat{k}}^{xz}\odot \hat{a}$$
  这两个结论的区别在于求$x$和$z$的互相关时，两个序列的相对滑动方向不一样。结合我的第一篇文章中关于相关运算用DFT加速的内容：
$${\hat{k}}^{zx}=\hat{z}\odot {\hat{x}}^{\ast }$$
$${\hat{k}}^{xz}=\hat{x}\odot {\hat{z}}^{\ast }$$
  我更倾向于使用后者的结论，那样在没有引入核函数的情况下计算$f(z)$可以表示为：
$$f(z)={\mathcal{F}}^{-1}(\hat{x}\odot {\hat{z}}^{\ast }\odot \hat{\alpha })$$

核函数计算

  有些核函数的计算也可以利用循环矩阵来加速。

多项式核

$$\kappa (u,v)=(u^{T}v+a{)}^b$$
  在训练的时候计算$k^{xx}$，$u$和$v$一个始终是原始样本$x$，另一个是循环移位的样本。这就是自相关运算，所以可以用DFT来加速。
$$k_i^{xx}=\kappa (x,P^{i}x)=(x^T{P}^{i}x+a{)}^b={\left({\mathcal{F}}^{-1}(\hat{x}\odot {\hat{x}}^{\ast })+a\right)}^b$$
  在检测的时候计算$k^{xz}$。这时候是互相关运算，互相关运算要注意它们相对滑动的方向，然后再确定哪个需要取复共轭。
$$k_i^{xz}=\kappa (z,P^{i}x)=(z^T{P}^{i}x+a{)}^b={\left({\mathcal{F}}^{-1}(\hat{x}\odot {\hat{z}}^{\ast })+a\right)}^b$$

高斯核

$$\kappa (u,v)=exp\left(-\frac{1}{{\sigma }^2}||u-v|{|}^2\right)$$
  在训练的时候计算$k^{xx}$。
$$k_i^{xx}=\kappa (x,P^{i}x)=exp\left(-\frac{2}{{\sigma }^2}\left(||x|{|}^2-{\mathcal{F}}^{-1}(\hat{x}\odot {\hat{x}}^{\ast })\right)\right)$$
  在检测的时候计算$k^{xz}$。
$$k_i^{xz}=\kappa (x,P^{i}x)=exp\left(-\frac{1}{{\sigma }^2}\left(||x|{|}^2+||z|{|}^2-2{\mathcal{F}}^{-1}(\hat{x}\odot {\hat{z}}^{\ast })\right)\right)$$