稀疏表达：向量、矩阵与张量

转载于：http://www.cvchina.info/2010/06/01/sparse-representation-vector-matrix-tensor-1/

稀疏表达是近年来SP, ML, PR, CV领域中的一大热点，文章可谓是普天盖地，令人目不暇给。老板某门课程的课程需要大纲，我顺道给扩展了下，就有了这个上中下三篇介绍性质的东西。遗憾的是，我在绝大多数情况下实在不算是一个勤快的人，这玩意可能充满bug，更新也可能断断续续，尽请诸位看官见谅了。顺道一提，ICCV09有一个相关的 tutorial 。

据传博文里公式数量和其人气是成反比例关系的，一个公式可以驱散50%的读者，我写完这个（上）之后点了点公式数量，觉得大约是要无人问津了。所以，在介绍稀疏表达之前，让我们先来展示下其在computer vision中的应用，吸引下眼球。

首先是图像恢复（以前有人贴过Obama还记得不），由左侧图像恢复出右侧结果

然后是类似的图像inpainting

然后是图像去模糊，左上为输入模糊图像，右下为输出清晰图像及估计的相机运动（其实是PSF），中间均为迭代过程:

再然后是物体检测（自行车），左侧输入图像，中间为位置概率图，右侧为检测结果

当然我个人还推荐Yi Ma的sparse face，这个在对抗噪声的效果上很棒，比如下图中左侧的那张噪声图像（你能辨认是哪位不？这方法可以！）

且说sparse representation这个概念，早在96-97年的时候就火了一把。最著名的大约要数Nature上的某篇文章，将稀疏性加入least square的regularization，然后得到了具有方向特性图像块（basis）。这样就很好的解释了初级视皮层（V1）的工作机理，即对于线段的方向选择特性。几乎同一时期，著名的LASSO算法也被发表在 J. Royal. Statist. Soc B。Lasso比较好的解决了least square (l2 norm) error + l1 norm regularization的问题。然而，这个时候绝大多数人没有意识到（或者没法解决）这l1 norm和稀疏性之间的联系。其实早在这之前，Osher等人提出的Total Variation （TV）已经包含了l1 norm的概念了，只不过TV原本是连续域上的积分形式。（啥？你不知道Osher…想想Level Set吧）

在进入现代的压缩感知、稀疏表示这一课题前，让我们来首先回顾下这一系列问题的核心，即线性方程组
其中矩阵，通常而言是满秩的。向量。现在已知，求解。学过线性代数的同学可能都会说：这个不难啊，因为， 故而这个方程组是欠定的，所以有无穷多组解啊，咱还可以算算基础解系啥的…

但是如果我们希望其解尽可能的稀疏：比如（即中非零元个数）尽可能的小。那么问题就会变得比较微妙了，下图给出了问题的形象示意。

换言之给定m维空间中一组过完备的基 ，如何选择最少个数的基向量，重构给定向量 ，其严格定义可以写成

时光之轮播快到2003~2004年，Donoho & Elad做了一个很漂亮的证明，如果矩阵满足某种条件，具体而言：

那么上文提及的0范数优化问题具有唯一的解。这里的是个比较诡异（请允许我使用这词）的定义：最小的线性相关的列向量集所含的向量个数（吐槽：明白了么，我做TA的时候就被这个问题问倒了）。本来想在这个概念上唠叨两句，后来发现了Elad的一个talk，清晰明了。

即便是唯一性得到了证明，求解这个问题仍然是NP难的。科研的车轮滚滚向前，转眼到了2006年，传奇性的华裔数学家Terrence Tao登场了，Tao和Donoho的弟子Candes合作证明了在RIP条件下，0范数优化问题与以下1范数优化问题具有相同的解：

其中RIP条件，即存在满足某种条件的（与N相关）常数:

RIP条件是对于矩阵列向量正交性的一种衡量（此处咱就不细说了）。其实早在1993年Mallat就提出过Mutual Coherence对于正交性进行度量，并提出了下文还要提及的matching pursuit方法。

实际上以上的1范数优化问题是一个凸优化，故而必然有唯一解，至此sparse representation的大坑初步成型。总结一下：
1. 如果矩阵满足，则0范数优化问题有唯一解。
2. 进一步如果矩阵满足RIP条件，则0范数优化问题和1范数优化问题的解一致。
3. 1范数优化问题是凸优化，故其唯一解即为0范数优化问题的唯一解。

进一步可以考虑含噪声情况，即

可以得到相似的结果，有兴趣的同学可以查阅相关文献。理论坑只有大牛能挖，但一般人也能挖挖这个优化算法啊，于是SP、ML、CV邻域里都有做这个优化算法的，这个出招可就真是五花八门了。据我所知，大致可以分为三大流派：
1. 直接优化

一般的方法是greedy algorithm，代表有Matching Pursuit, Orthogonal Matching Pursuit
2. 优化

还记得上面提到的LASSO么，这就是它的模型。
3. 如果已知拉格朗日乘子，优化无约束凸优化问题

解这个的方法现在基本上soft thresholding的方法一统天下，常见的有coordinate descent, Bregman Iteration (又是Osher)等
4. 如果未知拉格朗日乘子，优化

这类方法又叫Homotopy，可以认为是3的扩展。核心出发点是objective function是的分段线性函数。

除此之外，还有利用p范数逐次逼近0范数的方法等等，此处不再赘述。顺道说一句，稀疏表示在不同的领域连名称都不同，搞信号的管这个叫basis pursuit，搞统计的叫l1 regularization….然后，让我们把话题拉回到Nature的那篇文章：如果我们不知道矩阵，只知道一堆向量。我们应当如何构造，使得在这一字典（矩阵）下的表示最稀疏？类比以上过程，这个问题被称为Dictionary Learning，可以写成以下优化问题：

这个东西可就相对麻烦了，最关键的是这个优化不是凸的（优化变量相乘）。所以一般的想法是block descent：首先固定，优化（相当于多个独立的1范数优化问题）；其次将计算出的固定，优化，这就是一个（可能带约束）的least square问题。如此反复，直到算法收敛到某个（局部）极小值。实际上解这个问题的方法目前有三种：efficient sparse coding algorithm NIPS 06; K-SVD tsp 06; Online dictionary learning for sparse coding, ICML 09 & JMLR 10。前两种都是batch的方法，后一种是online的，据个人测试最后一种的方法比前两者要快很多很多….下面这个是我利用ICML09的方法从1200张彩色图像中训练出一组过完备基，具有比较好的方向特性。

最后，还记得本文开头的那些demo么？INRIA做了一个sparse representation的matlab工具包SPAMS，虽然不开源，但其效率（大部分时候）是现有公开工具包之冠（底层用了intel的MKL），利用这个工具包，几行简单的matlab代码就可以几乎实现以上提及的所有demo了….大家有兴趣的话，欢迎尝试^_^

下期预告：借着collaborative filter的东风，Candes在08年又挖出了matrix completion的新坑。于是，当向量的1范数推广到矩阵的迹范数（trace norm）之后…..

---

（中）

在开始正文之前，咱首先得说明一下，这篇东西偏向于理论，各位看官可以自行跳过某些部分。这方面的工作奠基人同样也是compressive sensing的大牛之一E.J Candes（Donoho的得意门生），以及Candes的学生Ben Recht，前者刚从caltech被挖到stanford，后者目前刚到wisconsin做AP。Candes大牛，stanford统计系出生，师从Donoho。Candes原来的主要工作集中在小波分析上（实际上C牛非常多产），比如著名的curvelets以及ridgelets，04年左右开始和Tao合作从事compressive sensing的理论工作，这里有他的简要介绍。

继续唠叨，上回说到借着collaborative filtering的东风，矩阵的稀疏表示受到了广泛的关注。说到矩阵的稀疏性，大部分看官可能有所误解。这个矩阵稀疏表示严格而言可以分为两种：
1. 矩阵元素的稀疏性，即矩阵非0元个数相对较少。参照向量的范数，同样可以定义矩阵的0范数，并将其松弛到矩阵的1范数的优化问题。
2. 矩阵奇异值的稀疏性，即矩阵奇异值中非0元的个数（即矩阵的秩）相对较少。仿照向量情况下0范数与1范数的关系，同样可以将其松弛的到迹范数（trace norm）的优化问题。

咱下面会分别聊聊这两个问题。首先，咱的出发点是machine learning中的collaborative filtering，这个概念并不是啥新东西了，最早大约可以追朔到1992的某篇同名文章。这玩意是做啥的呢，通俗的说，每次你在淘宝上闲逛的时候，下面都会有一行推荐商品。这些个网络服务商（淘宝，Amazon, Ebay）就在想了，如果这个推荐系统做的足够好，那么消费者（比如你我）的购物欲望就会得到刺激，这个销量也就上去了。实际上，这和超市里玲琅满目的货架是一个道理。

这里就得提提Netflix Prize这件事了，话说netflix是家在线dvd租赁公司，这公司就抱了同样的想法。不过这家公司想了个主意：该公司提供数据，出资100万美刀，奖励研发这个推荐系统算法的小组，并要求这些算法发表在学术会议或期刊之上。这可以算是现实版的百万富翁了（学术和money两不误），于是collaborative filtering着实火了一把（比如SIGKDD上的不少文章）。最终历时两年，由AT&T实验室成员组成的BellKor’s Pragmatic Chaos赢得了这100万刀。顺到一提，国内也有不少家伙参与了这个Prize，比如排名第二的Ensemble组里就能看到中科院某所学生的身影。

这个推荐系统咋做呢？我们先从简单的模型开始。以netflix为例，netflix有个影评系统，在你租完DVD以后会让你打分（1-5分）。当然不是所有人都会认真去打，实际上只有少数家伙会给打分（这世界上懒人何其之多）。同样，对每个用户而言，他也只可能给部分看过的DVD打分。假设现在有个用户和部电影，如果把所有评分列成一张大表，可以得到矩阵。其中，每一行对应一个用户的评分，每一列对应一部电影的用户评价。可以想象，这个矩阵中只有少部分元素是已知的（图1）。

从现有的用户数据，来预测未知的用户数据，这就是collaborative filtering了。那么这个东西怎么实现呢？解释起来难，做起来容易，这个模型放在在topic model里叫做Probabilistic latent semantic analysis （PLSA），放在代数里叫做矩阵分解（Matrix Fatorization）或者矩阵填充（Matrix Completion），这里就只能形象的解释下。虽然用户千奇百怪、电影成千上万，但总可以归结为若干类型：比如有腐女向、宅男向电影之分，再比如有悲剧也有喜剧。如果把这些latent factor画成一个空间，那么不同的用户群体应当位于这个latent factor空间的不同位置，体现了不同用户的喜好。如果可以把用户喜好连同潜在的latent factor一同计算出来，预测也自然水到渠成了。从某种角度来看，奇异值分解过程也就是上述的剥离latent factor和用户喜好的过程，这其中的philosophy可以参见这篇文章。

咱首先要谈的是矩阵奇异值的稀疏性，为此先来回忆下奇异值分解。
1. 奇异值非负，即
2. 奇异值非0元的个数即为矩阵的秩（rank）
如果把奇异值写成对角矩阵的形式（比如SVD分解的标准形式），其对角元为。进一步，矩阵的迹范数（trace norm）定义为矩阵奇异值之和，即有

现在我们可以把collaborative filtering的基本问题回顾一下，给定一张推荐数据表，已知其下标子集中的元素（也就是有评分的部分），如何恢复这个矩阵？这就是matrix completion的问题了…

乍眼一看，这基本就是mission impossible了，即使只有一个元素未知，这个矩阵也不可能唯一。但是如果我们加一些限制条件，这个问题就变得有趣起来了。Candes考虑的是这么一个问题：

其中表示在子集上的投影（即只取子集上的对应元素）。实际上，同样的问题可以有不同的表达形式，如果把这个优化问题稍作调整，可以得到相对容易解释的模型：

其中Frobenius范数也就是矩阵的2范数。从这个式子来看，我们希望找到这么一个矩阵，使得其在已有的数据上和用户评分尽可能的一致（2范数意义下），同时具有比较低的秩（受到上限的约束）。这里对于秩的约束，很多时候是为了降低模型自身的复杂度（比如collaborative filtering，multiple instance learning）。当然，这里也可以看成是一个fidelity term + regulariztion term的经典形式。

实际上矩阵的rank是一个不那么友好的函数，rank自身是非凸、不连续的，最后的结果就是对于rank的优化问题是NP难的。类比0范数与1范数的关系，矩阵的秩（rank）相当于这个对角阵的0范数；矩阵的迹范数（trace norm）相当于这个对角矩阵的1范数。为此，如果这个对角矩阵足够稀疏，即矩阵的秩，那么可参照向量的稀疏表示，利用矩阵的迹范数(trace norm)代替矩阵的秩（rank）。

同样，由于迹范数(trace norm)是凸的，上式是一个凸优化问题，故而必有唯一的最优解。如果这种近似是可以接受的，那么这个问题自然也就解决了。

这种近似靠谱么？这就是Candes和Recht回答的关键问题。Candes从random orthogonal model出发，证明了在此假设下从某个秩为的真实矩阵中均匀抽取个元素，且满足(这里不妨设，反之只需要转置即可)

则凸优化问题的唯一最优解至少以概率逼近原始矩阵，即有

其中均为某常数。更进一步，如果矩阵的秩足够小，对于元素数量的要求会进一步降低。

咱来聊聊这个结果，这说明在random orthogonal model假设成立的条件下，如果相对于比较小，那么只需要知道这个矩阵中约个元素，就可以很高的概率恢复出这个矩阵。举例而言，如果我们有一个秩为10的矩阵，那我们大致只需要从中随机抽取约270万个元素就可以（以很高概率）恢复出原始矩阵了(当然270万貌似也是一个很大的数，但原始矩阵约含有1700万个元素…）。实际上，这是一个相对保守的界，Recht在此基础上还进行了一系列的理论工作。自从出现了这个之后，under mild condition，大家都把rank直接放成trace norm了…从实用的角度来说，Candes告诉我们用凸优化去近似一个NP问题，可能得到很好的解。从实验结果来看（代码见此），这种近似有时候效果一流，但有时候也根本不work（违背了假设条件），故而具体问题还得具体对待。

虽然早在04年NIPS上，就有人提出了类似的优化方法（MMMF），用trace norm代替rank，并且ML领域中也确实有不少类似的工作。但是，Candes的工作解决了根本的理论问题，并为一系列的rank minimization的问题指了一条出路。这里有一个比较有意思的地方是，MMMF是从构造最大间隔线性分类器的角度出发来考虑matrix factorization的问题，并且用的是low norm，但和matrix completion的模型本质上是差不多的，两者关系大家可以自行推导下。

咱接着要讨论的是矩阵元素的稀疏性，这个工作也和Candes有着很大的关系。咱先把上面的公式照着copy一遍：

如果咱已知矩阵的全部元素，这个东西类似很常见的PCA了：

这样问题就变成了去噪+降维。进一步把F范数（2范数）改写为0范数：

为啥是0范数呢，这是基于这么一种假设：误差相对于总体样本而言总是稀疏的。于是，我们可以引入辅助变量表示误差，并把上式稍作改写：

这里的用于平衡矩阵的秩和误差的稀疏性。同样，rank和0范数什么的都是相当讨厌的东西，于是咱松弛一下，就有

这就是Robust Principle Component Analysis (RPCA) 或者Principle Component Pursuit 的核心模型了。这幅图很好的说明了RPCA和PCA的区别（转自Yi Ma主页）。

PCA	RPCA

(libo：上图注释有误，两图中蓝色为RPCA，绿色为PCA，两图分别展示了在有噪声，和有错误数据corruption的情况下RPCA和PCA的结果)

说起RPCA，这里岔开两句，这个东西原来是Yi Ma的学生John Wright发在NIPS09上的一篇文章。结果接收之后，被Candes指出了一个bug（审稿人没看出来），于是Candes对这个问题进行了考虑，从而就有了一篇叫做《Robust Principal Component Analysis?》的文章（preprint）。Candes证明了在同matrix completion基本相同的假设下，这种近似以很高的概率恢复精确结果（详细结果可见RPCA的论文）。特别的，此时可以简单选择参数。Matrix Completion（MC）和 RPCA在Yi Ma的主页上有一个简单的介绍，上面列举了相关文献与代码的链接。

MC和RPCA在computer vision上有啥用呢？John Wright在NIPS的文章里做了两个实验：背景建模，人脸阴影去除。大家有兴趣可以查查cvpr 10的paper, 有用MC来做video denoising的，有用RPCA来做人脸对齐的…还有那篇best paper也是紧密相关。咱本来还想聊聊这些模型的优化算法，鉴于篇幅所限，就只能留到（下）篇去了。

看了Ma Yi的工作后，才开始对Sparse represenation开始感兴趣。从最开始的惊奇，到中间的觉得无用，再到现在觉得还是有点用，其对这整套想法的把握真可谓一曲三折。不过总得来说，都是在项目和实验的驱动下进行的。

和ML的人不同，我们关注的点在CV中的application。解决图像分类中重要的一个方法是bag of features，也就是形成很多words图像patch后，去fit新的图，然后送去投票编码，最后分类。

其中一个很关键的步骤是怎么得到这些words的图像patch？一种方法是K-means做聚类。当然这个方法有很多naive的假设，不过总得来说效果凑合。所以最近几年大家都在死命搞features extraction，直到发现最后无论如何都打不败的是SIFT,HOG之流。

当然，科学总在进步中。大家发现feature这块都做不动了，自然开始动其他的念头，要不转行去做其他应用，要么得用ML的新工具了。最近在ICCV上发现在K-means上的重大改进，其核心就是放松了原来kmeans聚类后分配新样本去哪个类的分配方式，加上考虑到什么locally constraints，最后取得了相当好的效果。仅仅和最基本的dense sift+linear svm就打败了multiple features+nonlinear svm。哈，仔细研究了论文，发现解那个式子，可以用lasso搞搞。没准回头可以实现一下。

sparse很好，但是并不是什么新问题和新方法。总结一下这个熟悉又陌生的新思潮吧。这个总结有很多版本，看过Prfans上的，不过现在prfans好像down掉了。

最最开始要算Terrance Tao在IT上那篇best papers，因为他是华人，长得特别像电影明星，然后20+就当UCLA的教授。。。所以很吸引人。不过关于这个方面在CV还有ML开始火还在于以下几篇文章作者的工作吧。

首先是这篇文章，
John Wright, Allen Yang, Arvind Ganesh, Shankar Shastry, and Yi Ma, Robust face recognition via sparse representation. (To appear in IEEE Trans. on Pattern Analysis and Machine Intelligence)

大致做法就是用training data表出新的data，即新的x是training data的线性组合，如果training data数大大小于feature数，那么这就相当于CS里的x_0=Ax，A的每列是个sample，x是每个sample的权重。于是给出新的 x_0，可以通过优化x的l_1找到最sparse的解，即用最少数量的sample表出新的点。然后还可以在Ax前乘个measurement matrix，用来做特征选择，其实就是进一步降维。我在yaleB上的实验也证明了这个方法还是比较好的。但是中间的小trick是：
如果我用来做遮挡问题，那么训练库中的样本出现问题后，测试集上的结果就将变得相当差。比如用张猴子脸训练，检测人脸的时候用许多脸来拼凑，最后得到的是猴子脸而不是人脸。。。
比较囧的是，好多CVPR和ICCV上2009的文章都来follow这种motivation的工作。

然后这篇文章也挺有意思的
Robert Calderbank, Sina Jafarpour, and Robert Schapire, Compressed learning: Universal sparse dimensionality reduction and learning in the measurement domain(Preprint, 2009)

其中， 证明了measurement space里train个svm和原空间里train的svm是同一个问题。关键或许是 投影到多个space里去，然后在每个space里 train，然后做boosting，类似于multi-view learning。感觉最近muliti-kernel learning是个大趋势呀。
当然一些研究还得结合自己问题的数据上的特殊性解决才行。

给一些关于Compression Sensing(Sparse Representation方面的链接：
Compressive Sensing资源主页：
Compressive Sensing Resources （最权威最全面的Compressive Sensing资源主页，几乎什么都能找的到）；
Compressive Sensing （和上面的差不多）；
Compressive Sensing Listing； 马毅的课程主页
Compressive Sensing Videos；Compressed Sensing Codes （还有Compressive Sensing Resources 的Software一栏中）；
Nuit Blanche；Compressive Sensing: The Big Picture；Terence Tao‘s What’s new;

理论方面的代表人物：
David Donoho; Emmanuel Candes;

转载自视觉计算研究论坛【http://www.sigvc.org/bbs/thread-1951-1-1.html】