matlab求kcf算法响应图_matlab求kcf算法响应图_剖析KCF

来源自我的博客

前言

核相关滤波技术是在单目标跟踪领域中占据重要地位的技术之一；而kernelized-based correlation filters(KCF)则是这一方法的核心框架。从一个新手的角度出发，我们可以逐步理解这一复杂的技术。

由于知乎的公式编辑与MathJax存在不兼容性问题，在某些情况下无法正常显示公式内容；因此，在部分页面中将这些无法正确渲染的公式转为图片格式显示以确保良好的阅读体验；如遇个别页面出现显示异常情况，请您移步至原博客查看

1.岭回归理论推导

其理论基础较为简单，在算法设计上与支持向量机中的结构风险最小化策略存在显著差异。具体而言，在保证误差风险最小的情况下尽量使得结构风险较小的优化目标下运行。相比之下，在保证误差风险最小的情况下尽量使得结构风险较小是一种更为保守的方式。值得注意的是，在该场景下的训练样本数量有限的情况下,KFC与Ridge Regression在应用上具有一定的兼容性。此外,Ridge Regression 不仅能够实现特征选择功能,还结合了对偶空间、傅里叶变换以及核函数等技术手段,从而进一步提升了算法性能。总体而言,两者之间的差异相对较小。

其中X是特征矩阵,w是权值向量,y是样本标签/响应变量,而每项都采用了L2范数的平方项,即该矩阵内所有元素的平方之和。因此,这一优化问题的关键在于确定最佳权值向量w,所采用的方法则是最直接地应用拉格朗日乘子法作为求解工具

我们假设当前的权重W和输出y都是一维向量，则矩阵的求导公式满足:

不过，由于后面要引入复频域空间，所以我们这里做一些微调:

其中，H 代表共轭转置，即在转置的同时将矩阵内所有元素变为其共轭形式，原因很简单：

2. 循环矩阵

2.1 循环矩阵的引入

为了高效确定目标中心位置，在目标跟踪中若采用滑动窗口处理，则需依次进行多次位置判断。这种方法的计算开销较大。为此作者提出了基于周期性结构的优化方法。这样我们就可以在每次迭代时不移动待选框而是通过将原始图像周期性延拓来实现移位操作。在一维情况下举例说明如下：

在每一步运算中, 矩阵中的每一行都会依次向右移位一位, 这个过程选取一个特定的例子是为了便于理解. 通过分析这一过程的具体步骤, 我们能够方便地实现任意给定矩阵X向右移动一位.

或者下移

，如：

2.2 循环矩阵的转换

循环矩阵本质上整合了各种循环移位的结果。尽管如此，在处理二维数据时

其中 F 与离散傅里叶变换中的矩阵有所差异，

就是原矩阵的傅里叶变换，

diag是将矩阵变为对角形式，后面会详细解释。

先以一维矩阵为例来证明：

Step1 定义循环矩阵 X 的多项式函数为：

这里先说明一下，单位矩阵 I 其实也是一个循环矩阵，而

其实就是将矩阵

I所有元素右移

n 个单位。

Step2 求矩阵K的特征值和特征向量：

经分析得知，在K与DFT变换所对应的基向量构成的基底具有相似性；基于此，在多项式域中进行分析可得循环移位操作对应于频域中的乘法运算这一重要性质。

， 利用矩阵与其特征值矩阵相似的特点, 可以很容易的证明该性质。

Step3 求循环矩阵 X 的特征值和特征向量：

观察到循环矩阵X的特征值等于其原始矩阵x的离散傅里叶变换结果。详细推导过程如下：

在这里我们将

替换为

DFT 变换矩阵

W ，利用矩阵对角化可知：

Step4 利用 DFT 变换矩阵 W 的性质修正 X :

通过观察可知 W 为对称矩阵，另外也可以轻松证明

，在这里呢，我们可以对W 进行适当地变换：

因此

，则

F 为酉矩阵，同时它也满足

，如果我们将之前的

W 替换为

F ，那么：

所以整体来看， X 保持不变，综上可得公式(2-3)。

2.3 二维循环矩阵

在之前的推导过程中，默认假设所有操作均基于一维矩阵展开。具体而言，在处理二维循环矩阵时，则需要先将该二维结构分解为两个独立的一维循环变换序列。作为补充说明，请您稍等片刻。由于篇幅限制以及复杂性考虑，这部分较为详细的推导过程我暂时无法为您呈现。

给定一个 m×n 维矩阵 X，在构造其循环版本时，默认将其视为将原矩阵划分为块结构。这种情况下每行元素均源自上一行向右依次后移一位的结果；与此同时每列元素均源自上一列向下依次后移一位的结果。这一操作从而形成一个 mn×mn 维方阵。这一过程主要是为了使循环矩阵 X 成为一个方阵以便后续求取其实现时所需的特征值与特征向量。若原矩阵 X 已经是方阵则无需执行上述步骤以生成其循环版本即可直接采用行向量右移的方式处理；而对于列向量则采用下移方式处理以符合观察习惯这也是论文代码所采用的具体实现方式。

依次将 n×n 大小的矩阵 x 的每一行向量分别视为一个 1×n 的块矩阵，并遵循一维循环矩阵Kx所遵循的规律进行操作；通过依次向下移动这些块矩阵的位置，则可以得到一个结果。

的块矩阵。然后再将每一列向量单独看作一个

的块矩阵,按照

的方式，不断右移，最后可得一个

的块矩阵。其中二维的 DFT 变换方式为

。

3. 循环矩阵与岭回归算法的结合

构建了循环矩阵 X 后，
若判断该 (i,j) 位置处的块矩阵响应强度最大，
即目标框相较于前一个目标框向下偏移 i-1 个单位,
向右偏移 j-1 个单位。
此时对应的标签 y 的尺寸即为 n×n。
将此方法与脊回归算法相结合后可获得：

其中，

表示

x 经 DFT 变换之后的共轭形式，

表示全 1 向量，其等同于单位矩阵的特征向量，

表示矩阵元素点乘。

根据DFT时域卷积的性质：

而时域卷积多用于循环卷积的情况，在这种情况下假设原始序列被视为一个周期性序列，在这种假设下可以通过理论分析和数值验证来得出结论

可以发现：

因此可得：

利用上面的结论可以继续转换w为：

在此时此刻，在计算权重矩阵w的过程中，在傅里叶空间中被转换为一种容易计算的点乘运算；其复杂度大大降低。

4. 对偶空间的引入

我对偶空间的具体意义仅有概念不清晰的印象。记得在学习运筹学时就感觉对偶空间像是从另一个层面分析优化问题。比如，在岭回归中权重矩阵W的作用就是实现X到y的映射关系。如果将循环矩阵扩展成多行向量

个 n×n 的样本，则更直接一点就是完成从

d 维空间到一维 空间的映射/变换。关于 对 偶 空间而言，则 需 要 考虑 如何...

每个样本对问题带来的影响都会被评估，并且这个因素在优化过程中通常被用作惩罚项的系数。随后会对每个约束在不同样本上的影响进行评估。

4.1 优化角度分析

原优化问题为：

其等价为：

采用惩罚函数的方式是：

可以发现，如果将

看作

通过施加特定线性变换后，在其对应的对偶空间中表现为一种特定的形式，则我们可以将其过程简化为一个独立变量，并运用拉格朗日乘数法后得到相应的结果

将其带入原目标函数可得：

由此可得原优化问题的对偶问题，将其转换成矩阵形式为：

利用拉格朗日乘子法可得：

再将此对偶空间的最优解带入公式(4-4)可得：

4.2 矩阵变换角度

该优化方法着重源于理论根源，在实际应用时可以直接采用该理论基础，并进而对矩阵进行相应的变换。

4.3 新样本测试

在训练完参数后, 当引入新的样本时, 可以直接通过 wx 方法求解响应 y, 由式（4-7）可知：

5.核函数的高维映射

5.1核函数的引入

为了缓解数据在原始空间中线性不可分的问题，在实际应用中我们通常会选择适当的变换方法。由于岭回归属于单层模型结构，并且其核心思路是从空间映射入手进行分析与优化。如常见的核函数包括线性、多项式以及高斯型（如径向基）等基础类型，在实际应用中我们通常会选择适当的变换方法。在本研究中则采用了 RBF 核函数这一特定形式作为选择依据，在理论基础上经过无限次变换后相当于嵌入了一个无穷维特征空间

如图所示，在图中显示出来的黑色与蓝色区域呈现出明显的线性不可分特征；然而通过引入一个二次函数却能够实现完美的分割；这正是核函数所发挥的重要作用。在本文的研究中采用的核函数具体为高斯型核函数：

论文中将核函数引入样本的点积，即：

5.2 核函数与循环矩阵的结合

对于任何的循环矩阵 X，还是以一维的原矩阵 x 为例，可知：

其中基于循环变换的基础作用矩阵P（即前面所述的K matrix），特别地为了区别于核函数的概念，在向量内积运算中（即两个向量各分量依次相乘），因此这两个向量同时平移一位位置。

，并不会影响结果。

由公式(5-3)可知，只要行号和列号的差值相同，其对应元素的值就相同，所以

是循环矩阵。再利用循环矩阵的特性(2-3)可得公式(5-4)：

同理可得：

虽然原论文关于w的推导错误了，但是代码是根据

来实现的，所以正确。

5.3 不同核函数的计算

从最基础的内积出发，其核函数形式就是：

对于多项式核，则有：

对于 RBF 径向基核，也就是常说的高斯核函数，有：

6. 模板图像的获取

模板图像是基于第一帧图像目标框所得到的，其具体获取过程如下：

第一步，在确保初始目标框中心固定位置的前提下，对目标框的宽度和高度进行按比例放大处理，并在论文研究中设定放大比例为2.5倍。

Step2 设定模板图像尺寸为 96，计算扩展框与模板图像尺寸的比例：

Step3 然后将 scale 同时应用于宽和高，获取图像提取区域：

第4步中，在提取HOG特征时是以cell单元的形式进行采样的。此外，在将频域中的直流分量平移到图像中心之前，则要求原始图像的高度和宽度均为cell尺寸的偶数倍。需要注意的是，在HOG特征向量化过程中，默认会忽略位于边缘附近的那些cell。因此还需要在此基础上再补充两倍于cell大小的信息。

Step5 考虑到 roi 区域可能会越界于原始图像边界，在超出部分采用与原始边界一致的像素值进行填充。

Step6 最后利用线性插值的方式将 roi 区域采样为 template 大小。

7.特征提取

7.1 f-hog特征

本段文字介绍了一种用于目标检测的特征提取方法：HOG（方向梯度直方图）。该方法通过计算图像中每个像素点的方向梯度信息，并将这些信息进行统计建模来描述目标物体的形状特性。与传统的HOG特征相比，在论文中使用的该方法具有显著的优势：它不仅能够有效捕捉目标物体在几何变换（如旋转、平移等）下的形状特性变化，并且还能较好地适应光学条件下的复杂场景变化。具体而言，在论文中使用的HOG特征相较于传统方法有何不同？主要体现在以下几个方面：首先，在计算方式上采用了更为高效的算法；其次，在统计建模的过程中引入了更多的参数变量；最后，在实验结果上表现出更高的识别率和鲁棒性。

Step1 梯度幅值计算。

对模板图像中的每个像素点，在RGB三个通道中分别计算水平方向上的dx和垂直方向上的dy梯度，并进一步计算各处的梯度幅度。取其最大值所在的通道作为基准。

Step2 梯度方向判定。

当像素点在通道m中的最大幅值出现时，则可以通过该通道的水平变化与垂直变化来确定该点的方向。研究中将0至180度划分为9个方向，在确定归属时，则会根据该像素点梯度在模板方向上的投影值来进行分类。

分开展示其方向并将其映射到这两个区间中，在[0,180 )和[0,360 )这两个区间内进行处理后，则每个像素点都会被赋予两种不同的方向

Step3 cell的分割。

设置 cell 单元尺寸为论文中所采用的4个单元长度单位（即每个 cell 的大小设为4），由此可知，在该研究中水平方向上划分为 sizeX=24（即x轴方向上有24个单元），而竖直方向上同样划分为 sizeY=24（y轴方向上也有24个单元）。在计算图像像素点梯度时采用了镜像法处理边界问题的方式进行计算操作，在这种情况下论文并未考虑边界的梯度变化情况。对于每个 cell内部所建立的方向梯度直方图而言，在本研究中应包含从9个到18个不等的具体分类统计结果（总计为27种不同的方向类别）。

Step4 cell 内像素点梯度幅值加权方式。

在论文代码中涉及cell的方向梯度分布的计算方式非常独特。对于每一个cell单元体而言，在确定其尺寸后我们创建了两个离散序列这些序列呈现出一定的规律性因此将其命名为x函数和y函数：

从图像中可以看出，在细胞内随着像素点横纵坐标位移时，相应的x和y值持续变化。此外这两个函数对称于0.5，并且满足x+y=1的关系。基于这一特点我们可以采用x和y来进行组合分解

根据公式（7-3），论文代码中将每个单元格划分为四个区域（左上、右上、左下和右下），这些区域均基于包含当前单元格在内的相邻四个单元格在同一位置的部分进行加权平均计算。其中权重值即来源于所述的公式（7-3）所对应的四个区域。举例说明时可选取左上方单元格作为参考对象。

观察到，在处理cell网格时，每个单元格的加权处理方式是基于自身左上角的位置，并结合相邻单元格的信息构建出一个虚拟单元格。随后根据附图中展示的比例关系对四个单元格的数据进行加权平均计算，并采用此比例作为权重设置

同理，在处理其余三个区域时（即其余三个区域均遵循类似的逻辑），例如，在处理右下方的区域时，则是以当前cell的右下角作为虚拟cell的右下角。具体而言，在处理左上方区域时，则是以当前cell的左上角作为虚拟cell的左上角；同样地，在处理正上方区域时，则是以当前cell的正上方作为虚拟cell的正上方；类似地，在处理正左方区域时，则是以当前cell的正左边作为虚拟cell的正左边。值得注意的是，在进行加权计算的过程中（即在构建每个新virtual cell的过程中），对于处于边界位置（即位于网格边缘）的情况（即在网格边缘的情况下），则需要特别注意只考虑实际存在的相邻单元格（即在网格边缘的情况下仅选取周边存在的部分单元格）。

Step5 方向梯度直方图计算。

对于 cell 内部的每个像素点，在将该像素点梯度幅值分别以 [0,180) 和 [0,360) 两个投影区间进行映射后累加至对应的梯度方向直方图中；在完成该 cell 特征的加权计算后，在该 cell 中最终保存了总计 27（9+18）个梯度方向信息

Step6 相对领域归一化及截断。

在每一个单元格中选取与其相邻的四个单元格（如图所示）。值得注意的是这个过程略显粗糙（好丑-_-||），因此我们不得不努力以区分这四个单元格的区别而努力

基于此, 四种不同的组合方法会分别提取其内部每个细胞中的方向梯度直方图, 并取其中前九个方向梯度进行计算, 其结果即为经过 L2 范数处理后的值 val。随后, 将该单元格内的所有27个方向上的梯度直方图进行归一化处理, 即将这些数据除以计算出的 L2 范数值 val, 这样便能够获得规范化后的 HOG 特征向量。通过这四种不同的组合方法可获得四组独立的 HOG 特征向量, 其总维度数为 (9×4) + (18×4) = 72。

注意到边界单元无法提供足够的方向信息，并且在实际应用中这些边缘区域的数据量往往不足以支持有效的分析或计算需求。因此需要移除这些边界单元，并重新计算网格维度大小：经过这样的处理后，则有 sizeX = 24 - 1 = 19 和 sizeY = 15 + 1 = 16

Step7 PCA降维。

该研究者从不同分辨率的大量图像中收集了多个具有36维特性的数据样本，并对这些高维数据进行了主成分分析。观察到的现象是，在前11个主成分向量所张成的空间中几乎涵盖了HOG（HOG）特征求取的所有信息。经过降维处理后的数据集在目标检测任务中的表现与使用原始36维特征时完全一致

如果用

表示第

i 组

hog 特征的第

j 个方向，则原作者代码中的降维方式分别如下：

然后将两种降维方式得到的特征进行组合，得到 27+4=31 组特征。

原论文流程示意图如下：

7.2 CN/CN2特征

该颜色特征将颜色空间划分为黑、蓝、棕、灰等共计11种类型，并将其投影至由标准正交基构成的10维子空间上，在此过程中作者列举了32768种可能的颜色向量组合。随后通过主成分分析法（PCA），采用奇异值分解提取了两个最重要的主成分作为最终的颜色特征。鉴于CN2特征在单目标跟踪领域内被广泛认可的人工特征特性，在此对其理论基础进行详细阐述：

Step1 根据如下计算公式，作者给出了模板矩阵w2c(32768×10)：

将原图中的RGB三通道数据输入其中，在图像中每个像素的位置均对应一个索引值（范围是1到32768），随后将这些索引值传递至w2c网络中。这样能够从而生成一个具有相同宽度和高度、且具有10个通道的新矩阵

x_pca，与此同时，将其形状重塑为(W×H)×10的二维矩阵；

Step2 逐帧更新外观矩阵：

Step3 开始PCA，先按列对矩阵去中心化，并计算协方差矩阵cov(10×10)：

Step4 进行奇异值分解，由于对于任何矩阵

都可以利用奇异值分解为

, 其中

,

,​

对于矩阵cov，在其中U等于V的情况下，则有AU等于UD。由此可知,U的所有列向量都是协方差矩阵对应的特征向量;而D代表协方差矩阵对应的特征值

Step5 取U的前2列特征向量，逆分解得到新的协方差矩阵：

Step6 更新协方差矩阵：

Step7 得到CN2特征：

通过查看下图中的原始图像及其灰度版本与CN2特征可视化图之间的对比关系

8. 算法实现

8.1 多通道图像特征矩阵求解

基于上述理论推导可知，在每一帧中计算可获得fhog特征。将该过程视作一个多通道系统，则可得到具有31个通道的图像特征矩阵；随后对各个通道施加二维汉宁窗以实现滤波效果。此方法旨在减少快速傅里叶转换过程中造成的频谱泄漏现象；其函数形式如下：

然后分别对各个通道求解其

，对各通道数据进行多通道处理，处理方式如下:

8.2 标签制作

利用第五章所得结论，可以求得

对于 ground truth，在论文中采用高斯分布函数来表示标签；这是因为模板函数的中心与目标框的中心重合而得出的原因。

其中,(cx,cy)代表图像特征矩阵的中心点,padding代表扩展框相对于目标框的变化幅度为2.5,

表示设定的一个值 0.125。

8.3 多尺度检测

该研究中作者采用了三种不同的尺度设定，并将这些参数组合作为统一的整体进行优化配置；随后依次应用这些参数设置并结合优化算法完成了数据处理流程

,

三种度量标准 scale 进行检测，在本节中讨论的度量标准作用于第六章所述的 ROI 矩阵

为了解决在更新操作中可能出现的目标框边缘问题（即左上角位于图像边界），作者采取了措施将这类目标框的左上角远离图像边界的边缘移动了一个单位

每一种尺度都可以计算出 f(z) 的相应频率分布情况；尽管该频率分布情况对应于循环卷积中每个子块的频率特性特征值计算结果；其中位于第 i 行第 j 列的子块对应的频率特性特征值等于目标框向右移动 i−1 个单位、向下移动 j−1 个单位后的频率特性的平均值；进而可以通过分析该 response matrix 确定其峰值 peak_value 及其对应的峰值位置

当最大响应位置位于图像边界之外时，在其左右两侧进行响应强度对比：若右侧强度高于左侧或下侧强度高于上侧，则将最大响应位置相应地向强度较高的那一侧移动一定距离

然后计算此位置与图像中心的距离 res。

在各个不同尺度过量层级上，在各个不同规模层次上（即不同的尺寸级别），都存在一个对应的惩罚权重参数 scale_weight。通过将此惩罚权重参数与各自维度下最大的特征响应值相乘的方式计算出每个维度的真实最大特征响应值，并选取对应于这些真实最大特征响应值所处的最大尺寸级别那个维度作为最佳尺寸级别（即 best_scale）。基于这一原则对目标框参数 T(x,y,w,h) 进行更新

8.4 模板更新

首先从当前帧中提取原始尺度下的fhog特征矩阵z，在作者提供的代码中，并未对这一区域应用汉宁窗进行滤波处理。这可能是出于在更新过程中认为该区域的重要性相对较低考虑，在这种情况下与其对应的模板信息可能也会随之发生变化

的更新公式如下：

9.总结

整体而言，本文利用循环矩阵方法显著降低了循环移位运算的计算复杂度。接着借助傅里叶变换将主要的矩阵乘法转换为点乘运算。随后进一步利用对偶空间和核函数的技术成功提升了岭回归分类器的表现能力。总体而言该算法设计思路清晰结构简洁这与其算法结构简单明了直接相关因而达到了较好的跟踪效果

不过其中也暴露了很多问题:

虽然模板和

该部分对当前帧的重要性较低；然而它对整个矩阵进行了调整。因此，在出现一段遮挡的情况时, 该模板将无法恢复; 这与 siam-fc 的模板更新不同, 深度网络通常只是更新目标框的位置信息以及大小调整;

在论文中采用了三种度量标准，在这些度量标准的应用范围相对有限的情况下表现出了一定的作用，在后来Martin Danelljan提出的改进算法DSST中对该方法的性能表现有了显著提升；然而，在DSST算法中采用了33种度量标准以实现更高的分类精度目标的同时也导致该方法的整体运行速度有所下降

当目标发生形变时（即形状发生变化），其结果会出现较大程度的负面效果；这是因为 KCF 的本质即为基于模板的匹配机制，在这种情况下（目标变形），自然难以实现精确匹配。

在论文中引入的汉宁窗虽能有效降低 FFT 的频谱泄漏效果，在其特定的空间分布特征下会导致边缘像素cell 的数值趋近于零的情况出现。这种现象从而导致大量数据信息的缺失。此外，在实际应用中可能会出现主响应峰值并非精确对齐理论峰值的情况，并且有可能导致边缘区域的目标信号被抑制或遗漏

参考资料

High-speed tracking technology utilizing kernel-based correlation filters was introduced by João F. Henriques in the year 2015, as published in IEEE TRANSACTIONS ON PATTERN ANALYSIS AND MACHINE INTELLIGENCE.

Forsyth D. Discriminatively Trained Part-Based Models for Object Recognition[J]. 2014.

https://www.jianshu.com/p/69a3e39c51f9

https://www.cnblogs.com/torsor/p/8848641.html

原作者C++源码：https://github.com/joaofaro/KCFcpp

QiangWang复现C++源码：https://github.com/foolwood/KCF

这篇博士学位论文的链接为：http://www.robots.ox.ac.uk/~joao/publications/henriques_phd.pdf