##主成分分析(PCA)是一个很普通的数据将维方法,但是也是很常用的一个方法。虽然不难理解,但是时间一长容易和LDA搞混,于是决定自己重新梳理一遍。
####数据集的0维代表
假
设
你
有
n
个
d
维
样
本
x
1
,
x
2
,
⋯
 
,
x
n
,
如
果
我
们
希
望
用
一
个
d
维
的
向
量
代
表
这
n
个
样
本
,
或
者
说
找
一
个
代
表
向
量
和
其
他
样
本
的
总
体
差
距
越
小
越
好
(
体
现
代
表
性
)
,
定
义
平
方
误
差
准
则
函
数
:
假设你有n个d维样本\mathbf {x_1,x_2,\cdots, x_n},如果我们希望用一个d维的向量代表\\ 这n个样本,或者说找一个代表向量和其他样本的总体差距越小越好(体现代表性)\\,定义平方误差准则函数:
假设你有n个d维样本x1,x2,⋯,xn,如果我们希望用一个d维的向量代表这n个样本,或者说找一个代表向量和其他样本的总体差距越小越好(体现代表性),定义平方误差准则函数:
J
0
(
x
0
)
=
∑
k
=
1
n
∥
x
0
−
x
k
∥
2
J_0(\mathbf x_0)=\sum_{k=1}^{n}\|\mathbf{x_0-x_k}\|^2
J0(x0)=k=1∑n∥x0−xk∥2
我们要找到能够使上面
J
0
(
x
0
)
J_0(\mathbf x_0)
J0(x0)最小化的那个
d
d
d维向量
x
0
\mathbf x_0
x0.容易猜到这个问题的解答是
x
0
=
m
\mathbf{x_0=m}
x0=m,其中
m
\mathbf m
m是样本均值,
m
=
1
n
∑
k
=
1
n
x
k
\mathbf m=\frac 1 n\sum_{k=1}^n\mathbf x_k
m=n1k=1∑nxk
我们直接对
x
0
\mathbf x_0
x0求偏导有:
∂
J
0
(
x
0
)
∂
x
0
=
∑
k
=
1
n
2
(
x
0
−
x
k
)
=
0
\frac {\partial J_0(\mathbf x_0)} {\partial\mathbf x_0}=\sum_{k=1}^n2(\mathbf{x_0-x_k})=0
∂x0∂J0(x0)=k=1∑n2(x0−xk)=0
2
n
x
0
=
∑
k
=
1
n
2
x
k
2n\mathbf x_0=\sum_{k=1}^n2\mathbf x_k
2nx0=k=1∑n2xk
样本均值是数据集的零维表达,它非常简单,以至于丢到太多信息,无法反映出样本之间的不同,完全可以假设所有样本都是
m
\mathbf m
m.所以我们想能不能用一维数据来代表样本集,显然这一维需要时能体现差别的,如我们有100名小学生数据,假设有年龄、性别,成绩,身高这4个记录,显然如果我们选择身高作为代表数据,显然不能去判断这个学生是不是优秀学生,性别也是如此,只有选取成绩这一项,才能体现他们的差距。那么接下里就引入PCA方法了。
####PCA分析
一个
n
维
向
量
,
如
何
变
为
一
维
呢
?
先
给
个
直
观
的
几
何
图
形
n维向量,如何变为一维呢?先给个直观的几何图形
n维向量,如何变为一维呢?先给个直观的几何图形
如图显示一个二维数据样本,
假
设
x
和
y
不
相
关
假设x和y不相关
假设x和y不相关,如果我们从
x
和
y
x和y
x和y中选取一维代替原先的数据,那么应该选择哪一维呢?很显然应该选择
x
x
x维,因为我们的直观感受是,
x
x
x维上的数据区分比较大,更能代表数据的特征。在将这种直观感受转化为数学语言前,我们还要考虑一个小问题,如何将
d
d
d维数据转化为一维?
从上面的图可以看出,我们通过投影方法可以将数据映射到坐标轴上。实际上是执行下面这个过程:以投影到x轴为例
z
=
(
1
,
0
)
T
(
x
,
y
)
z=(1,0)^T(x,y)
z=(1,0)T(x,y)
用向量的形式写为:
z
=
w
T
x
\mathbf{z=w^Tx}
z=wTx
于是乎我们就回忆起中学学过的向量了
点
P
P
P在向量
A
B
⃗
上
的
投
影
点
P
′
满
足
\vec {AB}上的投影点P'满足
AB上的投影点P′满足
∣
A
P
′
∣
=
∣
A
P
∣
c
o
s
θ
=
A
P
⃗
⋅
A
B
⃗
∣
A
B
⃗
∣
|AP'|=|AP|cos\theta=\frac{\vec{AP}\cdot\vec{AB}}{|\vec {AB}|}
∣AP′∣=∣AP∣cosθ=∣AB∣AP⋅AB
取
∣
A
B
⃗
∣
=
1
|\vec{AB}|=1
∣AB∣=1,若A为原点,则P点的投影坐标为:
z
=
A
P
⃗
⋅
A
P
⃗
z=\vec{AP}\cdot\vec{AP}
z=AP⋅AP
显然AB不用是x轴,A点不是原点也没有关系,因为它不影响数据集的分布,总是可以通过去中心化转化为一维坐标。
####数学推导一
好了我们差不多可以进行数学分析了。假设
e
\mathbf e
e为数据集要投影的直线上的单位向量,那么这条直线方程可以表示为:
x
=
m
+
a
e
\mathbf{x=m+}a\mathbf e
x=m+ae
其中
a
a
a为一个实数标量,表示直线上的某个点离
m
\mathbf m
m的距离。
假
设
原
数
据
点
x
k
在
直
线
上
的
投
影
为
m
+
a
k
e
,
那
么
通
过
最
小
化
平
方
误
差
准
则
函
数
,
我
们
能
够
得
到
一
系
列
最
优
a
k
的
集
合
,
其
过
程
如
下
:
假设原数据点\mathbf x_k在直线上的投影为\mathbf m+a_k\mathbf e,那么通过最小化平方误差准则函\\数,我们能够得到一系列最优a_k的集合,其过程如下:
假设原数据点xk在直线上的投影为m+ake,那么通过最小化平方误差准则函数,我们能够得到一系列最优ak的集合,其过程如下:
J
1
(
a
1
,
⋯
 
,
a
n
,
e
)
=
∑
k
=
1
n
∥
(
m
+
a
k
e
)
−
x
k
∥
2
J_1(a_1,\cdots,a_n,\mathbf e)=\sum_{k=1}^n\|(\mathbf m+a_k\mathbf e)-\mathbf x_k\|^2
J1(a1,⋯,an,e)=k=1∑n∥(m+ake)−xk∥2
=
∑
k
=
1
n
∥
a
k
e
−
(
x
k
−
m
)
∥
=\sum_{k=1}^n\|a_k\mathbf {e}-(\mathbf x_k-\mathbf m)\|
=k=1∑n∥ake−(xk−m)∥
=
∑
k
=
1
n
a
k
2
∥
e
∥
2
−
2
∑
k
=
1
n
a
k
e
T
(
x
k
−
m
)
+
∑
k
=
1
n
∥
x
k
−
m
∥
2
=\sum_{k=1}^n a_k^2\|\mathbf e\|^2-2\sum_{k=1}^n a_k\mathbf e^T(\mathbf{x}_k-\mathbf m)+\sum_{k=1}^n\|\mathbf x_k-\mathbf m\|^2
=k=1∑nak2∥e∥2−2k=1∑nakeT(xk−m)+k=1∑n∥xk−m∥2
由于
∥
e
∥
=
1
\|\mathbf e\|=1
∥e∥=1,通过对
a
k
a_k
ak求偏导,并令结果为0,我们可以得到:
a
k
=
e
T
(
x
k
−
m
)
a_k=\mathbf e^T(\mathbf x_k-\mathbf m)
ak=eT(xk−m)
从几何意义上说,只要把向量
x
k
\mathbf x_k
xk向通过样本均值的直线
e
\mathbf e
ez作投影就能够得到最小方差结果。而减去样本均值就是对应于前面所说的去中心化:
####寻找投影方向
现在剩下的问题是寻找投影方向$ \mathbf e$从前面的图形例子中我们可以感受到:我们要投影的方向实际上就是数据集方差最大的方向。根据香农的信息论得知,数据越混乱,其包含的信息量就越大,熵越大。这里我们可以用方差来表示数据的波动性,也就是说:方差就是信息方差越大,信息量就越大。
对于高维数据,我们定义离散度矩阵:
S
=
∑
k
=
1
n
(
x
k
−
m
)
(
x
k
−
m
)
T
\mathbf S=\sum_{k=1}^n(\mathbf x_k-\mathbf m)(\mathbf x_k-\mathbf m)^T
S=k=1∑n(xk−m)(xk−m)T
这个离散度矩阵是半正定的,为了证明我们只需要证明该矩阵的特征值不小于0即可。具体过程见补充部分一。
投影点的均值:
μ
=
1
n
∑
k
=
1
n
(
m
+
a
k
e
)
=
m
+
e
n
∑
k
=
1
n
a
k
\mathbf \mu=\frac 1 n\sum_{k=1}^n(\mathbf m+a_k\mathbf e)=\mathbf m+\frac{\mathbf e} n\sum_{k=1}^n a_k
μ=n1k=1∑n(m+ake)=m+nek=1∑nak
=
m
+
e
n
∑
k
=
1
n
e
T
(
x
k
−
m
)
=
m
=\mathbf m+\frac{\mathbf e} n\sum_{k=1}^n \mathbf e^T(\mathbf x_k-\mathbf m)=\mathbf m
=m+nek=1∑neT(xk−m)=m
那么方差为:
V
=
1
n
∑
k
=
1
n
∥
m
+
a
k
e
−
m
∥
2
\mathbf V=\frac 1 n\sum_{k=1}^n\|\mathbf m+a_k\mathbf e-\mathbf m\|^2
V=n1k=1∑n∥m+ake−m∥2
=
1
n
∑
k
=
1
n
∥
e
(
e
T
(
x
k
−
m
)
)
∥
2
=
∑
k
=
1
n
e
t
(
x
k
−
m
)
(
x
k
−
m
)
T
e
=\frac 1 n\sum_{k=1}^n\|\mathbf e(\mathbf e^T(\mathbf x_k-{\mathbf m}))\|^2=\sum_{k=1}^n\mathbf e^t(\mathbf x_k-\mathbf m)(\mathbf x_k-\mathbf m)^T\mathbf e
=n1k=1∑n∥e(eT(xk−m))∥2=k=1∑net(xk−m)(xk−m)Te
=
e
t
S
e
=\mathbf e^t\mathbf{Se}
=etSe
好了我们终于得到一个简单的式子,最大化
e
t
S
e
\mathbf e^t\mathbf{Se}
etSe,别忘了我们还有一个约束条件
∥
e
∥
=
1
\|\mathbf e\|=1
∥e∥=1,对于有约束条件的极值问题我们一般用拉格朗日乘子法,于是得到
L
(
e
,
λ
)
=
e
t
S
e
−
λ
(
e
t
e
−
1
)
L(\mathbf e,\lambda)=\mathbf e^t\mathbf{Se}-\lambda(\mathbf e^t\mathbf e-1)
L(e,λ)=etSe−λ(ete−1)
对这个式子求偏导得:
∂
L
∂
e
=
2
S
e
−
2
λ
e
=
0
\frac{\partial L}{\partial \mathbf e}=2\mathbf{Se}-2\lambda\mathbf e=0
∂e∂L=2Se−2λe=0
S
e
=
λ
e
\mathbf{Se}=\lambda\mathbf e
Se=λe
发现什么了没有?没错,
λ
和
e
\lambda和\mathbf e
λ和e分别对应离散度矩阵的特征值和特征向量,我们选取特征值最大的特征向量作为投影方向。显然结论不局限于二维平面,对于高维数据也是成立的,对于高维数据我们往往是将其投影到超平面。实际操作步骤见补充二。
####小结
用PCA将维,只需要求出离散度矩阵的特征值和特征向量即可,对于普通数据,直接调用matlab函数便可立即求解,而对于上百万维数据,则需要用迭代的方法求解特征值,我们不必求出全部特征值,只需要选取前面若干个最大特征值即可。特征值在矩阵分析中具有非常重要的地位,很多数据挖掘算法都会用到特征值这一概念,所以请务必深刻理解特征值。
####补充一
我们在求
e
t
S
e
\mathbf e^t\mathbf{Se}
etSe没有考虑最大的存在性,由数学原定知道,当
S
\mathbf S
S是半正定矩阵时,
e
t
S
e
\mathbf e^t\mathbf{Se}
etSe存在最大值。定义为
S
=
∑
k
=
1
n
(
x
k
−
m
)
(
x
k
−
m
)
T
\mathbf S=\sum_{k=1}^n(\mathbf x_k-\mathbf m)(\mathbf x_k-\mathbf m)^T
S=k=1∑n(xk−m)(xk−m)T
我们可以做数据预处理:将样本的每一维减去对应维的样本均值,那么离散度矩阵(实际为样协方差矩阵的n-1倍)可以写为:
S
=
∑
k
=
1
n
x
k
x
k
T
=
x
1
x
1
T
+
⋯
+
x
n
x
n
T
\mathbf S=\sum_{k=1}^n\mathbf x_k\mathbf x_k^T=\mathbf x_1\mathbf x_1^T+\cdots+\mathbf x_n\mathbf x_n^T
S=k=1∑nxkxkT=x1x1T+⋯+xnxnT
=
(
x
1
,
x
2
,
⋯
 
,
x
n
)
(
x
1
T
,
x
2
T
,
⋯
 
,
x
n
T
)
T
=(\mathbf{x_1,x_2,\cdots,x_n)}(\mathbf x_1^T,\mathbf x_2^T,\cdots,\mathbf x_n^T)^T
=(x1,x2,⋯,xn)(x1T,x2T,⋯,xnT)T
令
X
=
(
x
1
,
x
2
,
⋯
 
,
x
n
)
\mathbf X=(\mathbf{x_1,x_2,\cdots,x_n)}
X=(x1,x2,⋯,xn)
则
S
=
X
X
T
\mathbf {S=XX}^T
S=XXT
设
λ
\lambda
λ为其特征值
ξ
\xi
ξ为对应的特征向量,则有:
S
ξ
=
λ
ξ
=
X
X
T
ξ
\mathbf S\xi=\lambda\xi=\mathbf{XX}^T\xi
Sξ=λξ=XXTξ
ξ
T
X
X
T
ξ
=
λ
∥
ξ
∥
2
=
∥
X
T
ξ
∥
2
≥
0
\xi^T\mathbf{XX}^T\xi=\lambda\|\xi\|^2=\|\mathbf X^T\xi\|^2\ge 0
ξTXXTξ=λ∥ξ∥2=∥XTξ∥2≥0
所以
λ
≥
0
\lambda\ge0
λ≥0
证毕。
####补充二
总结一下PCA的算法步骤:
设有m条n维数据。
1)将原始数据按列组成n行m列矩阵X
2)将X的每一行(代表一个属性字段)进行零均值化,即减去这一行的均值
3)求出协方差矩阵 C = 1 m X X T C=\frac{1}{m}XX^\mathsf{T} C=m1XXT
4)求出协方差矩阵的特征值及对应的特征向量
5)将特征向量按对应特征值大小从上到下按行排列成矩阵,取前k行组成矩阵P
6)Y=PX即为降维到k维后的数据
####推荐博客
###关于LDA
本来是想再自己写一下关于LDA的东西的,但是发现已经有很多博客,而且都写得很好,所以这里搬一些连接过来就好了。LDA与PCA最主要区别是:LDA是针对多个类别做投影,而PCA是对原始数据做投影,对多个类别做投影不能只考虑最大化投影方差,还要考虑类之间的间隔,记住这一点就不会搞混了。
下面是参考文章,注意文章后面给的参考资料也是很好的