第5章-维度归约

原文参考
https://onlinecourses.science.psu.edu/stat857/node/11
http://www.doc88.com/p-315762247283.html
http://wenku.baidu.com/link?url=OrBkfJFwZ_l8soN6CRE2B1KAI-988Z0wrzG_K59W39l7ycRTNd6xeH6hFUnRugEywPz1TUsvclcjpeiS8JLJsnd97GrQdHeG4r70BtSTIsS
https://www.cs.nyu.edu/~roweis/papers/llejmlr.pdf
http://arxiv.org/pdf/cs/0212008v1.pdf
多元统计分析
等…

---

本章主要学习目标：

在第4章中，我们了解了压缩变量的算法，接下来我们熟悉一下高维数据的降维方法。

传统的降维方法包括PCA、LDA，流行学习算包括LLE、ISOMAP、LE、LTSA，如下图：

1, 主成分分析PCA

首先我们定义，输入矩阵X是n*p维的。行代表样本，列代表变量。且X是中心化的，每一列去掉均值。

1.1 奇异值分解(SVD)

由现行代数的定理可知，如果一个矩阵A是对称矩阵，那么其一定可以正交对角化。即，存在一个正交矩阵P和一个对角矩阵D使得
$A=PD{P}^{-1}=PD{P}^T$
那么，A不是对称矩阵的时候呢？
这时就引出了奇异值分解，且任何矩阵都能进行奇异值分解。

假设A是m*n维矩阵，那么$A^{T}A$一定是对称的，可正交对角化。假设v是$A^{T}A$的单位正交基构成的特征向量、λ是对应的特征值，则有：
$A^{T}Av=\lambda v$
$||Av|{|}^2=v^T{A}^{T}Av=v^T\lambda v=\lambda (v是由单位正交基构成的)$
$||Av||=\sqrt{||Av|{|}^2}=\sqrt{\lambda }=\sigma$
定义$\sigma =\sqrt{\lambda }$为矩阵A的奇异值，且奇异值是向量Av的长度。

我们令向量$u=\frac{Av}{||Av||}=\frac{Av}{\sigma }$
则我们得到，Av=σu

同理，对于n×p 输入矩阵X，可以分解为$X=UD{V}^T$

1.2 主成分

主成分分析是线性降维的主要方法：

样本X的协方差矩阵，为
$S=X^{T}X/N$
(注明：因为X已经进行了中心化的处理，所以协方差矩阵可以如上式表示)

如果对$X^{T}X$进行奇异值变换，得
$X^{T}X=VD{U}^{T}UD{V}^T=V{D}^2{V}^T$

向量v称为X的(principal component directions)主成分方向坐标基（v已经单位化），如果将X投影到该坐标基上面，就得到了对应的主成分分量。

综上，我们容易得出$z_j=X_j\ast v_j=u_j\ast d_j$
$u_j$是$x_j$的行向量沿着$v_j$向量方向上投影得到的，投影的尺度通过$d_j$来控制.

所以，

• X的主成分为$z_j=u_j\ast d_j，j=1,...,p$
• 主成分应该是按照其方差由大到小排序的，根据下面公式可知，按照奇异值的大小排序即得可得到我们想要的排序.
  $var(z)=\frac{1}{N}\ast (z-0{)}^2$
  $var(z)=\frac{1}{N}(u\ast d{)}^2=\frac{1}{N}{d}^T{u}^{T}ud=\frac{d^2}{N}$
  (因为已经中心化，所以均值为零)

1.3 主成分回归-principal components regression (PCR).

前面我们讲的PCA，是一种无监督学习方法，而PCR就是主成分与线性回归的结合，假设提取了M个主成分$z_1,...,z_m$，然后求

$y={\beta }_1{z}_1+...+...+{\beta }_m{z}_m$

有时候，PCR的效果并不一定比第6章讲的方法要好，只不过PCA是一个很好的降维方法，在PCA的基础上衍生出来的PCR也是对高维样进行回归的一个重要的参考方向。

1.4 偏最小二乘回归-partial least squares (PLS).

偏最小二乘回归结合了多元线性回归、主成分分析以及典型相关分析。

偏最小二乘回归与主成分回归的不同之处在于，主成分回归产生的权重矩阵反映的是预测变量X之间的协方差，而偏最小二乘回归产生的权重矩阵反映的是预测变量X与响应变量Y之间的协方差。

之所以分解矩阵变为X与Y的协方差矩阵，是因为，主成分分析中，会忽略掉那些可能对Y起决定性作用的、但是与其他变量相关性又很小的变量。且主成分只能很好的解释了X，但是并不能保证一定与Y有关。相对的，PLS中分解得到的成分(潜变量)，是抽取了与Y相关的X的主要信息。这点上，使得PLS可能比PCR得到更准确的回归。

(1) 典型相关分析

典型相关分析(Canonical Correlation Analysis)，是考察一组X变量和一组Y变量(两个及以上Y变量)之间的线性关系的方法。

1936年，Hotelling提出典型相关分析。考虑两组变量的线性组合, 并研究它们之间的相关系数p(T,U).在所有的线性组合中, 找一对相关系数最大的线性组合, 用这个组合的单相关系数来表示两组变量的相关性, 叫做两组变量的典型相关系数, 而这两个线性组合叫做一对典型变量。在两组多变量的情形下, 需要用若干对典型变量才能完全反映出它们之间的相关性。下一步, 再在两组变量的与t1,u1不相关的线性组合中, 找一对相关系数最大的线性组合, 它就是第二对典型变量, 而且p(t2,u2)就是第二个典型相关系数。这样下去, 可以得到若干对典型变量, 从而提取出两组变量间的全部信息。

令，典型变量可被表示为如下形式：
$t=a_1{x}_1+...+a_m{x}_m;$
$u=b_1{x}_1+...+b_n{x}_n;$

然后，我们可以求解t,u的相关系数: Cor(t,u)
最后，用线性规划的方法求解使得Cor(t,u)最大的参数a和b。

(2) 典型相关的思想下，如何求解PLS：

假设要对，因变量$Y=Y_1,...,Y_p$，
自变量$X=X_1,...,X_m$，的样本对象进行建模。

PLS的做法是：

首先，在自变量集中提取第一成分$T_1$
$T_1=w_{11}{X}_1+...+w_{1m}{X}_m$
同时，也在因变量集中提取第一成分$U_1$
$U_1=v_{11}{Y}_1+...+v_{1p}{Y}_p$

求解目标，为
$max(T_1)->max$
$max(U_1)->max$
$Corr(T1,U1)->max$
且，$||w_1||=||v_1||=1$

然后，建立$Y$与$T_1$的，以及$X$与$T_1$的回归方程。如果回归方程的精度达到满意程度，则终止算法；否则继续提取第二成分…(假设直到$T_r$为止)。

$X_0=t_1{\alpha }_1^{\prime }+E_1$
$Y_0=t_1{\beta }_1^{\prime }+F_1$
其中，$t_1$为n维向量，
${\alpha }_1^{\prime }=({\alpha }_{11},...,{\alpha }_{1m}),{\beta }_1^{\prime }=({\beta }_{11},...,{\beta }_{1p})$
$E_1和F_1分别是残差矩阵$
如果精度没有满足，则继续用回归得到残差矩阵$E_1,F_1$继续回归，得到$t_2,{\alpha }_1^{\prime }...$

最后，PLS通过建立Y与$T_1,...,T_r$的回归方程，进而表示为Y与原自变量X的回归方程，实现偏最小二乘回归方程。

关于PLS的R实现，可以阅读此书。

2 线性判别分析LDA

LDA的全称是Linear Discriminant Analysis（线性判别分析），是一种supervised learning。有些资料上也称为是Fisher’s Linear Discriminant，因为Ronald Fisher在1936年发明了这个算。

LDA的原理是将带有类标签的数据投影到低维空间，投影后，使得统一类标签的样本数据"更近"，不同类之间"更远"。LDA算法涉及到了降维的思想，而且其本身也是一种分类方法。

具体实现步骤如下：

1. 样本的类别数为c，下面以c=2为例；
2. w是n维空间中任一向量，x在w上的投影为$w^{T}x$
3. $m_i$为所属c类的x的样本的质心，令
   $m_1=\frac{1}{N_1}{\sum }_{n\in c_1}{w}^T{x}_n$;
   $m_2=\frac{1}{N_2}{\sum }_{n\in c_2}{w}^T{x}_n$;
4. Fisher 判别准则，为
   $J(w)=\frac{(m_2-m_1{)}^2}{s_1^2+s_2^2}$
   等价于：
   $J(w)=\frac{w^t{S}_{B}w}{w^t{S}_{w}w}$
   其中，
   $S_B=(m_2-m_1)(m_2-m_1{)}^T$
   表示类间离散度矩阵
   $S_W$
   表示类内离散度矩阵

5，求解W
W是我们要求得的投影向量。我们通过W生成新的坐标轴，使不同类别的样本数据在新Y轴的投影中能最大程度的分离。

$max: $ $J(w)=\frac{w^t{S}_{B}w}{w^t{S}_{w}w}$
$s.t. $ $W^{t}W=I$

在人脸识别中，Belhumeur等人用PCA将人脸图像投影到一个子空间，然后再实施LDA，得到最终的投影矩阵$U=U_{pca}{U}_{lda}$ 命名为"Fisherfaces"

此外，基于LDA，还有许多衍生的判别模型：
Mixture Discriminant Analysis (MDA)
Gaussian Discriminant Analysis(GDA)
Quadratic Discriminant Analysis (QDA)
Flexible Discriminant Analysis (FDA)

3 局部线性嵌入LLE

简单的讲，该方法是把高维空间数据点按维数映射到低维嵌入空间，即Xi→Yi。

整个问题最后被转化为两个二次规划问题。

LLE(ocal linear embedding)算法可以归结为三步:

1. 寻找每个样本点$X_i$的k个近邻点；
2. 由每个样本点$X_i$的近邻点计算出该样本点的局部重建权值矩阵$W_{ij}$;
3. 由该样本点的局部重建权值矩阵$W_{ij}$和其近邻点计算出对应的低维向量$Y_i$.
   

$E(W)={\sum }_i{|X_i-{\sum }_j{W}_{ij}{X}_j|}^2$ - - - - - - - -Equation(1)
其中，$X_j$表示$X_i$点的第j个近邻，且对于非近邻点$W_{i}j=0$,对于近邻点${\sum }_j{W}_{ij}=1$

在限制条件下，通过上述方程1，最小化重构错误E(W)得到的最优权值W,遵循如下对称特性，即对于特定的数据点，在其本身和其邻居数据点有旋转、缩放、平移操作时将保持其原有性质不变。

$f(Y)={\sum }_i{|Y_i-{\sum }_j{W}_{ij}{Y}_j|}^2$ - - - - - - - -Equation(2)
其中，$Y_i$是输出向量，$Y_j$表示$Y_i$的第j个近邻，
且满足$\sum y_i=0$；$\frac{1}{N}\sum y_i{y}_i^{\prime }=I$

根据方程1得到权重W后，根据上述方程2，求解得到输出向量Y。将所有的样本点映射到低维空间，即在低维空间中，根据W的几何信息进行点的重建。

根据方程2的满足条件，求解：min(f(x))

在整个求解过程中，近邻K的个数是很重要的，但是其仅仅是个经验参数，并没有很好的办法进行确定，通常算法规定K必须大于样本输出的维数。

选取最优的K

应满足，在选取K条件下，应用LLE算法，使得输入点的相对位置与输出点的相对位置尽量保持一致。通常映射的好坏，可以通过输入点与输出点之间的偶合情况来表示。（输入点之间的距离矩阵和输出点之间的距离矩阵的关系系数衡量）

$k_{opt}=argmin(1-{\rho }_{D_x,D_y}{)}^2$

其中，ρ表示两个矩阵间的标准线性相关系数。
$D_x$表示输入样本（高维）之间的距离矩阵，$D_y$表示输出样本（低维）之间的距离矩阵。

4 拉普拉斯特征映射LE

Laplacian Eigenmaps与LLE方法有一些共有的特性。都是对几何流行(manifold)进行降维方法。

拉普拉斯矩阵的一些特性已经被广泛用于聚类(如谱聚类)和划分问题。某种意义上讲，聚类和降维是一个问题的正反面，所以LE算法试图通过拉普拉斯矩阵的相关方法把问题从谱聚类中转向降维。

假设：
在$R^l$空间中，有$x_1,...,x_k$，k个点

我们的目标是
在$R^m$(m远远小于l)空间中，找到$y_1,...,y_k$
使得在低维空间中的$y_i$可以"代表"高维空间中的$x_i$

LE算法：

1). 构建邻接图( Adjacency Graph)
计算两个近邻点$x_i,x_j$之间的距离，如果认为两个点是“相近”的，则在两点之间添加一条边。如何判断"相近"，有两种方法：

a). 根据阈值ε来判断$||x_i-x_j|{|}^2<\epsilon$
凡是满足上述等式的就算作近邻
缺点：很难确定ε；优点：保留了几何上的对称性

b). 根据个数n来判断
通过计算n个最近的点作为近邻
缺点：非几何特性驱动，可能很远的第n个点也作为近邻
优点：很容易选择，不会出现断开的图

2). 计算边的权重
计算边的权重来描述两个点之间的“远近”程度，同样有两种方法来计算边的权重：

a). 利用热核函数，如果点i和点j是连接的(确定是近邻点)，则
$W_{ij}=e^{-\frac{||x_i-x_j|{|}^2}{t}}$
如果不是近邻，则，$W_{ij}=0$

b). 简单模式：t->∞ ， 即有连接的时候$W_{ij}=1$，否则为0.

3). 特征映射(eigenmaps)

拉普拉斯矩阵的定义：
假定一个图G的邻接矩阵为A，度矩阵为D，则该图的拉普拉斯矩阵L为： $L=D-A$
其中，D为度矩阵，是图中每个点的度(出度入度不重复计算)，所以D是对称矩阵，且除了对角线其余位置元素均为零。

结合前两点，此处，我们定义：
$D_{ii}={\sum }_j{W}_{ji.}$

区别：
当权重的计算用的是热核函数(t≠∞)，那么其权值可能为非整数，所以矩阵D中的元素可能是非整数的，这是与常规定义的度矩阵的差异。

所以，LE算法中的拉普拉斯矩阵定义为：$L=D-W$

计算其特征值与特征向量：$Lf=\lambda Df$
可以得到k个特征向值，并对其进行排序，得到
${\lambda }_0,...,{\lambda }_{k-1}$，且
$0={\lambda }_0\le {\lambda }_1,...\le ...,{\lambda }_{k-2}\le {\lambda }_{k-1}$

不考虑0值，我们取较大的特征值λ
$x_i->(f_1(i),...,f_m(i))$
作为m维欧几里得空间

5 局部切空间排列算法LTSA

LTSA算法的思想是找出每个数据点的邻近点，用邻域中低维切空间的坐标近似表示局部的非线性几何特征；再通过变换矩阵将各数据点邻域切空间的局部坐标映射到统一的全局坐标上；最后，经过一些列数学推导，将求解整体嵌入坐标问题转换为求解矩阵的特征值问题，从而实现高维数据的维度约简。

具体实现步骤如下：

1. 选取邻域。
   对样本点$x_i$，选取包含其自身在内的k个近邻点作为邻域；$x_{ij}(j=1,..,k)$是其k个邻近点；${\bar{x}}_i$ 为上述k个邻近值的均值；
2. 局部线性拟合。
   计算点$x_i$处d维切空间的正交基$Q_i$和每一个$x_{ij}$在切空间上的正交投影${\theta }_j^i=Q_i^T(x_{ij}-{\bar{x}}_i)$
3. 局部坐标整合
   根据n个局部投影${\theta }^i=[{\theta }_1^i,{\theta }_2^i,...,{\theta }_k^i]$,计算全局坐标 T= { τ i } i = 1 n \{\tau_{i}\}^n_{i=1} {τi​}i=1n​，
   ${\tau }_i=[{\tau }_{i1},{\tau }_{i2},...,{\tau }_{ik}]$

6 MDS

MDS是通过对样本数据的相似性矩阵(相关矩阵或距离矩阵)进行适当的谱分解，得到能够近似表示高维数据的相对构造的、在低维空间中的投影，进而对原数据进行降维的一种方法。其思想和主成分分析（PCA）非常相似。在欧氏距离之下，这种算法求得的主坐标（principal coordinates）相当于PCA中前t个主成分的得分。

假设输入矩阵为
$X=x_1,x_2,...,x_n\in R^p；$

X的距离阵为$D=(d_{ij})$
$d_{ij}^2=(x_i-x_j{)}^T(x_i-x_j).$

算法步骤：

1), 根据距离阵D计算矩阵B

定理： 令
$A=(a_{ij}),a_{ij}=-\frac{1}{2}{d}_{ij}^2$,
$H=I_n-\frac{1}{n}{1}_n{1}_n^{\prime }$ 是中心化的变化矩阵
$B=HAH$
则，$D$ 是欧式型的距离阵的充要条件：
$B$ 是半正定矩阵，即$B\ge 0$

必要性证明：

设D是欧式型的，则存在$x_i$，使得两点的距离:
$d_{ij}^2=(x_i-x_j{)}^T(x_i-x_j)=-2a_{ij}$

由矩阵B的定义，得
$B=HAH=A-\frac{1}{n}AJ-\frac{1}{n}JA+\frac{1}{n^2}JAJ$
,其中$J=1_n{1}_n^{\prime }$

可得，
$b_{ij}=(a_{ij}-{\overline{a}}_{i.}-{\overline{a}}_{.j}+{\overline{a}}_{..})$
其中，
${\overline{a}}_{i.}=\frac{1}{n}{\sum }_{j=1}^n{a}_{ij}$
${\overline{a}}_{.j}=\frac{1}{n}{\sum }_{i=1}^n{a}_{ij}$
${\overline{a}}_{..}=\frac{1}{n^2}{\sum }_{i=1}^n{\sum }_{j=1}^n{a}_{ij}$

再带入$-2a_{ij}=(x_i-x_j{)}^T(x_i-x_j)$

得，
$b_{ij}=(x_i-\overline{x}{)}^T(x_i-\overline{x})$

可推得
$B=(HX)(HX{)}^T\ge 0$

必要性证毕.

充分性证明：

当B≥0时，记Λ = diag(${\lambda }_1,...,{\lambda }_k$)
为B的特征根对角矩阵，
Γ 为对应的n*k 维特征向量矩阵

则，$B＝\Gamma \Lambda {\Gamma }_T$
定义 $B＝\Gamma {\Lambda }^{1/2}$,则 $B＝Z{Z}^T$

2. 对矩阵B进行谱分解

7 ISOMAP

Isomap以MDS为计算工具，创新之处在于计算高维流形上数据点间距离时，不是用传统的欧式距离，而是采用微分几何中的测地线距离（或称为曲线距离），并且找到了一种用实际输入数据估计其测地线距离的算法（即图论中的最小路径逼近测地线距离)。

8 KPCA

入svm一样，KPCA运用了核方法对技巧，将原空间中的数据映射到高维空间.
$K(x_i,x_j)=\phi (x_i{)}^T\phi (x_j)$
利用核技巧，我们可以非线性空间$\phi (x_i)$使用pca方法