====== 第三课 非监督机器学习方法 ======

===== 3.1 主成分分析法 =====
==== 3.1.1 主成分分析法定义 ====

主成分分析，即principal components analysis, 简称PCA，也被翻译成“主分量分析”或“主元分析”。PCA旨在利用降维的思想，把多指标转化为少数几个综合指标。在统计学中，主成分分析是一种简化数据集的技术。它是依赖于数据的线性变换。这个变换把数据变换到一个新的坐标系统中，使得任何数据投影的第一大方差在第一个坐标(称为第一主成分)上，第二大方差在第二个坐标(第二主成分)上，依次类推。主成分分析经常用减少数据集的维数，同时保持数据集的对方差贡献最大的特征。这是通过保留低阶主成分，忽略高阶主成分做到的。这样低阶成分往往能够保留住数据的最重要方面。但是，这也不是一定的，要视具体应用而定。
      

==== 3.1.2 主成分分析法基本原理 ====

主成分分析法是一种降维的统计方法，它借助于一个正交变换，将其分量相关的原随机向量转化成其分量不相关的新随机向量，这在代数上表现为将原随机向量的协方差阵变换成对角形阵，在几何上表现为将原坐标系变换成新的正交坐标系，使之指向样本点散布最开的p个正交方向，然后对多维变量系统进行降维处理，使之能以一个较高的精度转换成低维变量系统，再通过构造适当的价值函数，进一步把低维系统转化成一维系统。  

==== 3.1.3 主成分分析法的主要作用 ====

概括起来说，主成分分析主要由以下几个方面的作用。 

  -主成分分析能降低所研究的数据空间的维数。即用研究m维的Y空间代替p维的X空间(m＜p)，而低维的Y空间代替 高维的x空间所损失的信息很少。即：使只有一个主成分Y<sub>l</sub>(即m＝1)时，这个Y<sub>l</sub>仍是使用全部X变量(p个)得到的。例如要计算Y<sub>l</sub>的均值也得使用全部x的均值。在所选的前m个主成分中，如果某个X<sub>i</sub>的系数全部近似于零的话，就可以把这个X<sub>i</sub>删除，这也是一种删除多余变量的方法。
  -有时可通过因子负荷a<sub>ij</sub>的结论，弄清X变量间的某些关系。
  -多维数据的一种图形表示方法。我们知道当维数大于3时便不能画出几何图形，多元统计研究的问题大都多于3个变量。要把研究的问题用图形表示出来是不可能的。然而，经过主成分分析后，我们可以选取前两个主成分或其中某两个主成分，根据主成分的得分，画出n个样品在二维平面上的分布况，由图形可直观地看出各样品在主分量中的地位，进而还可以对样本进行分类处理，可以由图形发现远离大多数样本点的离群点。
  -由主成分分析法构造回归模型。即把各主成分作为新自变量代替原来自变量x做回归分析。
  -用主成分分析筛选回归变量。回归变量的选择有着重的实际意义，为了使模型本身易于做结构分析、控制和预报，好从原始变量所构成的子集合中选择最佳变量，构成最佳变量集合。用主成分分析筛选变量，可以用较少的计算量来选择量，获得选择最佳变量子集合的效果。 
==== 3.1.4 主成分分析法的基本步骤 ====

  - 原始指标数据的标准化采集$p$维随机向量$x=(x_1, x_2, \ldots, x_p)^{\top}$。给定$n$个样本$x_i=(x_{i1}, x_{i2}, \ldots, x_{ip})^{\top}$，$i=1,2, \ldots, n$ $(n＞p)$，可构造样本阵：\[X = 
\left( \matrix{ x_{11} & x_{21} & \cdots & x_{n1}\\x_{12} & x_{22} & \cdots & x_{n2}\\... & \vdots & \ddots & \vdots \\x_{1p} & x_{2p} & \cdots & x_{np} } \right)\]
  - 对样本阵X进行标准化变换得标准化矩阵$Z = X - \overline{X}$, 其中\[\overline{x} = \frac{1}{n} \sum_{i=1}^{n}{x_i}, \overline{X} = (\overline{x}, \overline{x}, \ldots, \overline{x}) \]
  - 对标准化阵$Z$求相关系数矩阵\[R=Z^{\top} Z; \]
  - 解样本相关矩阵$R$的特征方程得$p$个从大到小排列的特征根\[ \lambda_1 \ge \lambda_2 \ge \ldots \ge \lambda_p \]将其对应的单位特征向量$u_i$确定为主成分
  - 将标准化后的指标变量转换为主成分表示,即\[ x_i = \overline{x} + \sum_{j=1}^{p}{ z_{ij} u_j }, \] 其中 \[ z_{ij} = <x_i, u_j > = x_i^{\top} u_j; \]
  - 根据需要选取前$k$个主成分构成降维的近似表示\[ x_i \rightarrow (z_{i1}, z_{i2}, \ldots, z_{ik} ) \]
 

==== 3.1.5 主成分分析法与奇异值分解 ====
奇异值分解是线性代数中一种重要的矩阵分解，在信号处理、统计学等领域有重要应用。奇异值分解在某些方面与对称矩阵或Hermite矩阵基于特征向量的对角化类似。然而这两种矩阵分解尽管有其相关性，但还是有明显的不同。对称阵特征向量分解的基础是谱分析，而奇异值分解则是谱分析理论在任意矩阵上的推广。奇异值分解在统计中的主要应用为主成分分析（PCA），它是一种数据分析方法，用来找出大量数据中所隐含的“模式”,它可以用在模式识别，数据压缩等方面。PCA算法的作用是把数据集映射到低维空间中去。 数据集的特征值（在SVD中用奇异值表征）按照重要性排列，降维的过程就是舍弃不重要的特征向量的过程，而剩下的特征向量张成空间为降维后的空间。

  * 奇异值分解
令$\widehat{X}$ 为一个居中的N×P维数据矩阵（设N>P）\[ X = \left[ \begin{array}{c|c|ccc}x_{0,0} & x_{1,0} & x_{2,0} & \cdots & x_{N-1,0}\\x_{0,1} & x_{1,1} & x_{2,1} & \cdots & x_{N-1,1}\\\vdots  & \vdots  & \vdots  & \ddots & \vdots\\x_{0,p-1} & x_{1,p-1} & x_{2,p-1} & \cdots & x_{N-1,p-1}\end{array} \right]_{N×p} = USV \] 即为$\widehat{X} $的奇异值分解，其中：
  * U是N×p 的正交阵，称作左奇异向量。
  * V是p×p 的正交阵，称作右奇异向量。
  * S是对角阵，其中，$d_1>=d_2>=\cdots >=d_p>=0$，称为奇异值。
奇异值分解总是存在的，并针对符号唯一。V的列即为主元，且$Z_j=U_jd_j$ 

  * Eckart-Young定理
令$S_q$ 为除了最前面的q个对角线元素之外的元素全为0的矩阵。则$\widehat{X_q}=US_qV^T$ 满足$\min_{rank(\widehat{X_q})=q} \big\Vert\widehat{X}-\widehat{X_q}\big\Vert$  

  * 主元分析例子-Eigenfaces(G. D. Finlayson, B. Schiele & J. Crowley.  Comprehensive colour image normalization. ECCV 98 pp. 475~490)
{{:keynote:pca1.jpg|}} {{:keynote:pca2.jpg|}} 

  * PCA与降维
通过奇异值分解进行空间变换：$X->Y$，形成从p维到q维的变换，其中p>>q，之后进行重现并计算误差。

==== 3.1.6 主成分分析法的不足 ====
1.只适用于一般的离散数据集，对于非线性数据集不适用。

2.虽然如此，PCA 方法还是在许多领域（例如人脸识别）得到广泛应用，此外，也有提出使用 Kernel 方法扩展 PCA 的算法，亦即先将原始数据从特征空间映射到一个高维空间中，消除其非线性的性质，然后在这样的空间中应用 PCA ，亦即 KernelPCA 方法（参见 Scholkopf, B. and Smola, A. and Muller, K.R., //Kernel principal component analysis//, Artificial Neural Networks—ICANN'97）。

除了 KernelPCA 之外，另一种将 PCA 用在非线性数据分析的方法是使用 LocalPCA ，现在比较热的流形学习中对于数据分布在一个低维流形上的假设实际上就是假定了数据在一个足够小的局部领域内是程线性的，这也正式 PCA 用武的地方，因此，在局部范围内使用 PCA 来分析流形结构也就成了一个顺理成章的事情。比如，在 tangent space learning 中，算法的目的是得出 manifold 的切空间，而如果使用局部方法的话，在一个流行局部的切空间其实就可以直接通过对局部 sample 进行 PCA 分解而求得。

<note important>Extended by  --- //[[pluskid@gmail.com|张驰原]] 2010/04/26 14:28//</note>

3.以下是对课件的补充和扩展（介绍RPCA及其在CVPR2010人脸图像上的应用）：

不适用于非高斯噪声污染的数据集。\\
设数据$m\times n$矩阵$M=L+S$，$L$为潜在的低秩矩阵，$S$噪声矩阵。如果$S$非高斯噪声，例如稀疏且幅值不定的噪声，那么PCA将失效。\\
{{:keynote:pca_done_well.jpg|}}
{{:keynote:pca_fail.jpg|}}\\
此时，宜用改进的模型RPCA(Robust PCA)。
RPCA通过以下目标函数求解：
<jsmath>\min_{L,S}\|L\|_{*}+\lambda \|S\|_{l_{1}} \;\; s.t. \; M=L+S</jsmath>
$\lambda$是权重参数，通常设为$\lambda=\frac{1}{\sqrt{\max\{m,n\}}}$。 $\|\cdot\|_{*}$为核范数，即矩阵奇异值之和。$\|\cdot\|_{l_{1}}$为一范数，即矩阵元素绝对值之和。\\
此凸函数具有唯一最小值。使用ALM(Augmented Lagrange Multiplier)求解，最小化增强的拉格朗日函数：
<jsmath>\min_{L,S,Y}l(L,S,Y)=\min_{L,S,Y}\|L\|_{*}+\lambda \|S\|_{l_{1}}+ tr\{Y^{T}(M-L-S)\}+\frac{\mu}{2}\|M-L-S\|_{F}^{2}</jsmath>
$\mu$为另一权重，可取值$\mu=nm/(4\|M\|_{1})$。$\|\cdot\|_{F}$为Frobenius范数，即矩阵元素平方和开根号。

通过迭代计算求最优值:\\
1.初始化：$S_{0}=Y_{0}=0$，$k=0$\\
2.while $\|M-L-S\|_{F}>10^{-7}\|M\|_{F}$ do\\
3.  $L_{k+1}=D_{1/\mu}(M-S_{k}-Y_{k}/\mu)$\\
4.  $S_{k+1}=S_{\lambda/\mu}(M-L_{k+1}+Y_{k}/\mu)$\\
5.  $Y_{k+1}=Y_{k}+\mu (M-L_{k+1}-S_{k+1})$\\
6.  $k=k+1$\\
7.end while\\
8.输出$L,S$\\

对于标量$x$，函数$S_{a}(x)=sgn(x)max\{|x|-a,0\}$，即截断操作。当输入为矩阵时，对每一元素独立操作。函数$D_a(X)$也是截断操作，但作用于矩阵奇异值：设矩阵$X$的SVD分解为$X=U\Sigma V^{T}$，那么$D_a(X)=US_{a}(\Sigma)V^{T}$。

对RPCA做一些修改，可用于人脸对齐与去噪。这是RPCA研究团队一篇CVPR2010的工作：\\
Yigang Peng, Arvind Ganesh, John Wright, Wenli Xu, and Yi Ma. RASL: Robust Alignment by Sparse and Low-rank Decomposition for Linearly Correlated Images. To appear in CVPR, June 2010.\\
{{:keynote:rals.jpg|}}\\
(a)算法输入40张同一人不同光照、遮挡，姿势和表情的图片。(b)-(d)为算法输出，算法同时完成对齐(b)、提取真实人脸（c）和分离噪声（异物）(d)。算法主要思路为寻找一组几何变换$\tau$使对齐后的图像$Dº\tau$能被分解为低秩成分和稀疏噪声。

 <note important> --- //[[fancij@139.com|胡振方]] 2010/04/27 21:50//</note>

===== 3.2 非线性降维算法 =====
==== 3.2.1 局部线性映射 LLE ====
  * LLE算法的定义
LLE算法针对非线性降维问题利用线性重构的局部对称性找出高维数据空间中的非线性结构并在保持各数据点临近位置关系情况下把高维空间数据点映射为低维空间对应的数据点其计算步骤包括计算寻找数据点或邻居数据点构造数据点及计算权值矩阵并通过权值矩阵计算低维向量。
  * 与PCA和MDS(多元测度)算法的相似与区别
相同点：都是降低维度的算法。可把相距很远的数据点映射在一个平面的临近位置。
不同点：PCA和MDS都是对高维空间中的可变因素的线性建模。LLE则是非线性的降维问题，不会涉及局部最小值问题。映射后，PCA和MDS无法保留流形的基本结构，而LLE则可以保证两点之间的测量距离不会改变。
  * LLE算法的步骤
 (1)寻找每个样本点的$k-$近邻点。把相对于所求样本点距离最近的$k$个样本点规定为所求样本点的$k-$近邻点，这里$k$是一个预先给定值。由Sam T. Roweis 和 Lawrence K. Saul等人所提出的LLE算法中，采用了欧氏距离以减轻计算的复杂度。然而，若假定高维空间中的数据是非线性分布的，也可采用Diijstra距离。Dijkstra距离是一种图上的测地距离，它能够保持样本点之间的曲面特性。这一度量在ISOMAP算法中有广泛的应用。针对样本点多的情况，普通的Dijkstra算法不能满足LLE算法的要求。
（2）由每个样本点的近邻点计算出该样本点的局部重建权值矩阵$W$使得
\[ \min  ε(w)=\sum_{i=1}^{N}|x_{i}-\sum_{j=1}^{k}w_{j}^{i}x_{ij}|^2\]
其中$x_{ij}$ $(j=1,2,...k)$为$x_{i}$的$k-$近邻点，$w_{j}^{i}$是$x_{i}$和$x_{ij}$之间的权值，且要满足条件：\[\sum_{i=1}^{k}w_{j}^{i}=1\] 这里求取W矩阵，需要构造一个局部协方差矩阵Q<sup>i</sup> 。
\[Q_{jm}^{i}={(x_{i}-x_{ij})}^{T}(x_{i}-x_{im}) \]
将上式与∑<sub>i=1</sub><sup>N</sup>w<sub>j</sub><sup>i</sup>=1相结合，并采用拉格朗日乘子法，即可求出局部最优化重建权值矩阵：
\[w_{j}^{i}=\frac{\sum_{m=1}^{k}(Q^{i})_{jm}^{-1}}{\sum_{p=1}^{k}\sum_{q=1}^{k}(Q^{i})_{pq}^{-1}}\]
在实际运算中，$Q^i$可能是一个奇异矩阵，此时必须正则化，即令$Q^i \leftarrow Q^i + rI$。其中$r$是正则化参数，$I$是一个$k \times k$的单位矩阵。
LLE算法的最后一步是将所有的样本点映射到低维空间中。映射条件满足如下所示：
\[\min  ε(Y)=\sum_{i=1}^{N}|y_{i}-\sum_{j=1}^{k}w_{j}^{i}y_{ij}|^2\]
其中，ε(Y)为损失函数值，y<sub>i</sub>是x<sub>i</sub>的输出向量，y<sub>ij</sub>(j=1,2,...,k)是y<sub>i</sub>的k个近邻点，且要满足两个条件，即：\[\sum_{i-1}^{N}y_{i}=0,\frac{1}{N}\sum_{i-1}^{N}y_{i}y_{i}^{T}=I\]
其中$I$是$m \times n $的单位矩阵。这里的w<sub>i</sub><sup>j</sup>(i=1,2,...,N)可以存储在$N \times N$的稀疏矩阵$W$中，当$x_j$是$x_i$的近邻点时，$W_{(i,j)}=w_j^i$。否则，$W_{(i,j)}=0$。则损失函数可重写为：
\[\min  ε(Y)=\sum_{i=1}^{N}\sum_{j=1}^{N}M_{i,j}y_{i}^{T}y_{j}\]
其中$M$是一个$N \times N$的对称矩阵，其表达式为：\[M=(I-W)^{T}(I-W)\]
要使损失函数值达到最小, 则取$Y$为$M$的最小$m$个非零特征值所对应的特征向量。在处理过程中，将$M$的特征值从小到大排列，第一个特征值几乎接近于零，那么舍去第一个特征值。通常取第2~$m+1$间的特征值所对应的特征向量作为输出结果。

(3) 由该样本点的局部重建权值矩阵和其近邻点计算出该样本点的输出值。

  * LLE算法伪代码
寻找X空间中的邻点
<code>
For i=1:N
  compute the distance from Xi to every other point Xj
  find the K smallest distances
  assign the corresponding points to be neighbors of Xi
end
</code>
重构权重W
<code>
for i=1:N
  create matrix Z consisting of all neighbors of Xi subtract Xi from every column of Z
  compute the local covariance C=Z'*Z
  solve linear system C*w = 1 for w
  set Wij=0 if j is not a neighbor of I   
  set the remaining elements in the ith row of W equal to w/sum(w);
end
</code>
计算输出矩阵Y
<code>
create sparse matrix M = (I-W)'*(I-W)
find bottom d+1 eigenvectors of M (corresponding to the d+1 smallest eigenvalues)
set the q-th ROW of Y to be the q+1 smallest eigenvector 
(discard the bottom eigenvector [1,1,1,1...] with eigenvalue zero)
</code>

==== 3.2.2 等距特征映射 Isomap ====
  * Isomap算法定义\\
ISOMAP是近年来用于非线性降维的一个重要算法，它是对经典MDS的扩展，也是寻求保持数据点之间距离的低维表示。算法的基本思想是利用样本向量之间的欧氏距离Dist(i,j)计算出样本之间的测地距离D<sub>G</sub>(i,j)，然后使用经典MDS算法来获得相应的低维投影。这样能够减少样本之间的欧式距离Dist(i,j)与D<sub>G</sub>(i,j)的误差，最大程度上保持高维数据内在的非线性（等距）几何结构，通过低维数据的分析来获得相应的高维数据特性，从而达到简化分析、获取数据有效特征以及可视化数据的目标。

  * Isomap主要步骤\\
(1)选取邻域构造邻域图G。\\
计算每个样本点同其余样本点之间的欧氏距离。当X<sub>j</sub>是X<sub>i</sub>的最近的k个点中的一个时，认为它们是相邻的，即图G有边X<sub>j</sub>X<sub>i</sub>(这种邻域称为k-邻域)。设边X<sub>j</sub>X<sub>i</sub>的权为d(X<sub>j</sub>,X<sub>i</sub>)，对P=1,……,N有
\[
d(X_j,X_i)=\sqrt{|X_{i1}-X_{j1}|^2+|X_{i2}-X_{j2}|^2+...+|X_{ip}-X_{jp}|^2}
\]

(2)计算任意两个样本向量之间的最短路径。\\
当图G有边X<sub>j</sub>X<sub>i</sub>时，设最短路径如d<sub>G</sub>(X<sub>i</sub>,X<sub>j</sub>)=d(X<sub>j</sub>,X<sub>i</sub>)，否则设d<sub>G</sub>(X<sub>i</sub>,X<sub>j</sub>)=$\infty$。对k=1,2,...,N，\\
d<sub>G</sub>(X<sub>i</sub>,X<sub>j</sub>)=min {d<sub>G</sub>(X<sub>i</sub>,X<sub>j</sub>),d<sub>G</sub>(X<sub>i</sub>,X<sub>k</sub>)+d<sub>G</sub>(X<sub>k</sub>,X<sub>j</sub>)}

(3)计算D维映射\\
将经典MDS应用到距离矩阵D<sub>G</sub>。
令
\[ D_Y(i,j)=\|y_i-y_j\|，D_G(i,j)=d_G(i,j), \]
最小化代价函数
\[E=\|\tau(D_G)-\tau(D_Y)\|_{L^2}\]
求解出D<sub>G</sub>的最大d个特征值以及对应的特征向量即可。


==== 3.2.3 BoostMap ====

  * 设计目标-大幅度降低图像数据库系统的获取时间
    
    嵌入（Embedding）被描述为一个机器学习任务，AdaBoost方法用来把多个简单的1维嵌入结合成为一个d维的嵌入；按照与一个查询对象的相似程度获得所有DB对象的排序。
  * 主要想法
    数据集S，距离函数D：S×S→ $R^1$ 
    
    {{:keynote:boostmap1.jpg|}} 
     
    要点：D可能计算复杂；D可能不满足三角不等式。
  * 问题定义
    嵌入被看做分类器，如果a与b或c接近，给出对a,b,c的估计。其中X为对象集，$D\le x$ 为距离函数。\[ \begin{equation} P_x(q,x _1,x_2)=\left\{\begin{array}{cl}1 & \textrm{if }D_x(q,x_1)<D_x(q,x_2)\\0 & \textrm{if }D_x(q,x_1)=D_x(q,x_2)\\-1 & \textrm{if }D_x(q,x_1)>D_x(q,x_2)\end{array}\right. \end{equation} \] \[ (即近似函数P(r:x,y) = sgn(D(y,r)-D(x,r)))\] 找到一个嵌入$F:X→R^d$ 和一个测度$D_{R^d}$ \[ \tilde{F}(q,x_1,x_2)=D_{R^d}(F(q),F(x_2))-D_{R^d}(F(q),F(x_1)),\] 用来评价任意三元组。
  * 例子：一些关于"北美城市”的近似函数 
     {{:keynote:boostmap2.jpg|}} 
     
  * BoostMap的输出
输出是一个分类器：$ H = /sum_{j=1}^dα_j\tilde{F}_j$ 最终输出是一个嵌入$ F:X→R^d$ 和一个距离测度$ D:R^d×R^d→R^d$ \[ \begin{equation} F(x)=(F_1(x),\ldots,F_d(x))\end{equation}\]  \[ D_{R^d}\left((u_1,\ldots,u_d),(v_1,\ldots,v_d)\right)=\sum_{j=1}^d(α_j\big| u_j-v_j\big|)\] 
    

==== 3.2.4 非线性降维算法总结 ====

  * Isomap
    * 不仅仅考虑欧氏几何结构
    * 多项式时间解

  * LLE
    * 不需要估计成对距离
    * 优化不牵涉局部最小化

  * BoostMap
    * 主要用于数据库中相似度检索
    * 离线训练，适合在线使用