PCA(Principal Component Analysis)，即主成分分析方法，是一种使用最广泛的数据降维算法（非监督的机器学习方法）。

其最主要的用途在于“降维”，通过析取主成分显出的最大的个别差异，发现更便于人类理解的特征。也可以用来削减回归分析和聚类分析中变量的数目。

在很多场景中需要对多变量数据进行观测，在一定程度上增加了数据采集的工作量。更重要的是：多变量之间可能存在相关性，从而增加了问题分析的复杂性。

如果对每个指标进行单独分析，其分析结果往往是孤立的，不能完全利用数据中的信息，因此盲目减少指标会损失很多有用的信息，从而产生错误的结论。

因此需要找到一种合理的方法，在减少需要分析的指标同时，尽量减少原指标包含信息的损失，以达到对所收集数据进行全面分析的目的。由于各变量之间存在一定的相关关系，因此可以考虑将关系紧密的变量变成尽可能少的新变量，使这些新变量是两两不相关的，那么就可以用较少的综合指标分别代表存在于各个变量中的各类信息。主成分分析与因子分析就属于这类降维算法。

PCA的主要思想是将n维特征映射到k维上，这k维是全新的正交特征也被称为主成分，是在原有n维特征的基础上重新构造出来的k维特征。

先假设用数据的两个特征画出散点图： 如果我们只保留特征1或者只保留特征2。那么此时就有一个问题，留个哪个特征比较好呢？ 通过上面对两个特征的映射结果可以发现保留特征1（右面）比较好，因为保留特征1，当把所有的点映射到x轴上以后，点和点之间的距离相对较大，也就是说，拥有更高的可区分度，同时还保留着部分映射之前的空间信息。 那么如果把点都映射到y轴上，发现点与点距离更近了，这不符合数据原来的空间分布。所以保留特征1相比保留特征2更加合适，但是这是最好的方案吗？

将所有的点都映射到一根拟合的斜线上，从二维降到一维，整体和原样本的分布并没有多大的差距，点和点之间的距离更大了，区分度也更加明显。 也就是说，我们要考虑的问题是：

如何找到让样本间距最大的轴？

其中，一般我们会使用方差（Variance）来定义样本之间的间距: V a r ( x ) = 1 m ∑ i = 1 m ( x i − x ˉ ) 2 Var(x)=\frac {1}{m}\sum_{i=1}^m(x_i-\bar x)^2 Var(x)=m1​i=1∑m​(xi​−xˉ)2

对于如何找到一个轴，使得样本空间的所有点映射到这个轴的方差最大。

将样本进行均值归0（demean），即所有样本减去样本的均值。样本的分布没有改变，只是将坐标轴进行了移动。 此时对于方差公式： V a r ( x ) = 1 m ∑ i = 1 m ( x i − x ˉ ) 2 , x ˉ = 0 Var(x)=\frac {1}{m}\sum_{i=1}^m(x_i-\bar x)^2,\bar x=0 Var(x)=m1​∑i=1m​(xi​−xˉ)2,xˉ=0。其中此时计算过程就少一项，这就是为什么要进行样本均值归0，可以化简的更加方便。

V a r ( X p r o j e c t ) = 1 m ∑ i = 1 m ( X p r o j e c t i − x ˉ p r o j e c t ) 2 Var(X_{project})=\frac {1}{m}\sum_{i=1}^m(X_{project}^i-\bar x_{project})^2 Var(Xproject​)=m1​i=1∑m​(Xprojecti​−xˉproject​)2 其实括号中的部分是一个向量，更加准确的描述应该是（向量的模）： V a r ( X p r o j e c t ) = 1 m ∑ i = 1 m ∣ ∣ X p r o j e c t i − x ˉ p r o j e c t ∣ ∣ 2 Var(X_{project})=\frac {1}{m}\sum_{i=1}^m||X_{project}^i-\bar x_{project}||^2 Var(Xproject​)=m1​i=1∑m​∣∣Xprojecti​−xˉproject​∣∣2

因为前面已经去均值，所以，这里只需下面的式子最大： V a r ( X p r o j e c t ) = 1 m ∑ i = 1 m ∣ ∣ X p r o j e c t i ∣ ∣ 2 Var(X_{project})=\frac {1}{m}\sum_{i=1}^m||X_{project}^i||^2 Var(Xproject​)=m1​i=1∑m​∣∣Xprojecti​∣∣2

映射过程如下： 红色的线是我们要找的方向 ω = ( ω 1 , ω 2 ) ω=(ω_1,ω_2) ω=(ω1​,ω2​) ；第i行的样本点 X ( i ) = ( X 1 ( i ) , X 2 ( i ) ) X^{(i)}=(X_1^{(i)},X_2^{(i)}) X(i)=(X1(i)​,X2(i)​)， X ( i ) X^{(i)} X(i)此时也是一个向量；映射到ω上做一个垂线，交点的位置就是 X p r o j e c t ( i ) = ( X p r o 1 ( i ) , X p r o 2 ( i ) ) X_{project}^{(i)}=(X_{pro1}^{(i)},X_{pro2}^{(i)}) Xproject(i)​=(Xpro1(i)​,Xpro2(i)​)对应的点；真正要求 X p r o j e c t ( i ) X_{project}^{(i)} Xproject(i)​的的模的平方，蓝色线段长度对应的平方。 把一个向量映射到另一个向量上，对应的映射长度是多少。实际上这种映射就是点乘： X ( i ) ⋅ ω = ∣ ∣ X ( i ) ∣ ∣ ⋅ ∣ ∣ ω ∣ ∣ ⋅ c o s θ X^{(i)}\cdot ω=||X^{(i)}||\cdot ||ω||\cdot cosθ X(i)⋅ω=∣∣X(i)∣∣⋅∣∣ω∣∣⋅cosθ 因为向量ω是我们要找的轴，是一个方向，因此使用方向向量就可以。因此长度为1： X ( i ) ⋅ ω = ∣ ∣ X ( i ) ∣ ∣ ⋅ c o s θ X^{(i)}\cdot ω=||X^{(i)}||\cdot cosθ X(i)⋅ω=∣∣X(i)∣∣⋅cosθ 因此，在三角形中有： X ( i ) ⋅ ω = ∣ ∣ X p r o j e c t ( i ) ∣ ∣ X^{(i)}\cdot ω=||X_{project}^{(i)}|| X(i)⋅ω=∣∣Xproject(i)​∣∣ 主成分分析法的目标是：

求ω，使得 V a r ( X p r o j e c t ) = 1 m ∑ i = 1 m ( X ( i ) ⋅ ω ) 2 Var(X_{project})=\frac {1}{m}\sum_{i=1}^m(X^{(i)}\cdot ω)^2 Var(Xproject​)=m1​∑i=1m​(X(i)⋅ω)2最大.

如果是n维数据，则有： V a r ( X p r o j e c t ) = 1 m ∑ i = 1 m ( X 1 ( i ) ω 1 + X 2 ( i ) ω 2 + … … + X n ( i ) ω n ) 2 Var(X_{project})=\frac {1}{m}\sum_{i=1}^m(X_1^{(i)}ω_1+X_2^{(i)}ω_2+……+X_n^{(i)}ω_n)^2 Var(Xproject​)=m1​∑i=1m​(X1(i)​ω1​+X2(i)​ω2​+……+Xn(i)​ωn​)2

主成分分析方法（PCA），是数据降维算法。将关系紧密的变量变成尽可能少的新变量，使这些新变量是两两不相关的，即用较少的综合指标分别代表存在于各个变量中的各类信息，达到数据降维的效果。

所用到的方法就是“映射”：将n维特征映射到k维上，这k维是全新的正交特征也被称为主成分，是在原有n维特征的基础上重新构造出来的k维特征。我们要选择的就是让映射后样本间距最大的轴。

其过程分为两步：

样本归0 找到样本点映射后方差最大的单位向量ω 最后就能转为求目标函数的最优化问题：

求ω，使得 V a r ( X p r o j e c t ) = 1 m ∑ i = 1 m ( X ( i ) ⋅ ω ) 2 Var(X_{project})=\frac {1}{m}\sum_{i=1}^m(X^{(i)}\cdot ω)^2 Var(Xproject​)=m1​∑i=1m​(X(i)⋅ω)2最大

此时，我们就可以用搜索策略，使用梯度上升法来解决。

对于最优化问题，除了求出严格的数据解以外，还可以使用搜索策略求极值。

在求极值的问题中，有梯度上升和梯度下降两个最优化方法。梯度上升用于求最大值，梯度下降用于求最小值。

梯度下降的思想及推导问题。已知：使损失函数值变小的迭代公式： ( a k + 1 , b k + 1 ) = ( a k + 1 − η ∂ L ∂ a , b k + 1 − η ∂ L ∂ b ) (a_{k+1},b_{k+1})=(a_{k+1}-η\frac {\partial L}{\partial a},b_{k+1}-η\frac {\partial L}{\partial b}) (ak+1​,bk+1​)=(ak+1​−η∂a∂L​,bk+1​−η∂b∂L​) 我们现在要求极大值，则使用梯度上升法。梯度的方向就是函数值在该点上升最快的方向，顺着这个梯度方向轴，就可以找到极大值。即将负号变为正号： ( a k + 1 , b k + 1 ) = ( a k + 1 + η ∂ L ∂ a , b k + 1 + η ∂ L ∂ b ) (a_{k+1},b_{k+1})=(a_{k+1}+η\frac {\partial L}{\partial a},b_{k+1}+η\frac {\partial L}{\partial b}) (ak+1​,bk+1​)=(ak+1​+η∂a∂L​,bk+1​+η∂b∂L​)

现在回到对PCA算法求解的问题上来，对于 PAC 的目标函数：

求 ω 使 得 f ( X ) V a r ( X p r o j e c t ) = 1 m ∑ i = 1 m ( X 1 ( i ) ω 1 + X 2 ( i ) ω 2 + … … + X n ( i ) ω n ) 2 求ω使得f(X)Var(X_{project})=\frac {1}{m}\sum_{i=1}^m(X_1^{(i)}ω_1+X_2^{(i)}ω_2+……+X_n^{(i)}ω_n)^2 求ω使得f(X)Var(Xproject​)=m1​∑i=1m​(X1(i)​ω1​+X2(i)​ω2​+……+Xn(i)​ωn​)2 最大，可以使用梯度上升法来解决。

首先要求梯度。在上式中，ω是未知数 X i X_i Xi​是非监督学习提供的已知样本信息。因此对 f （ X ） f（X） f（X）中ω的每一个维度进行求导： ᐁ f = ( ∂ f ∂ ω 1 ∂ f ∂ ω 2 … ∂ f ∂ ω n ) = 2 m ( ∑ i = 1 m ( X 1 ( i ) ω 1 + X 2 ( i ) ω 2 + … + X n ( i ) ω n ) X 1 ( i ) ∑ i = 1 m ( X 1 ( i ) ω 1 + X 2 ( i ) ω 2 + … + X n ( i ) ω n ) X 2 ( i ) … ∑ i = 1 m ( X 1 ( i ) ω 1 + X 2 ( i ) ω 2 + … + X n ( i ) ω n ) X n ( i ) ) ᐁf=\begin{pmatrix} \frac {\partial f}{\partial ω_1}\\ \frac {\partial f}{\partial ω_2}\\ …\\ \frac {\partial f}{\partial ω_n}\end{pmatrix}=\frac {2}{m}\begin{pmatrix} \sum_{i=1}^m(X_1^{(i)}ω_1+X_2^{(i)}ω_2+…+X_n^{(i)}ω_n)X_1^{(i)}\\ \sum_{i=1}^m(X_1^{(i)}ω_1+X_2^{(i)}ω_2+…+X_n^{(i)}ω_n)X_2^{(i)}\\ …\\ \sum_{i=1}^m(X_1^{(i)}ω_1+X_2^{(i)}ω_2+…+X_n^{(i)}ω_n)X_n^{(i)}\end{pmatrix} ᐁf=⎝⎜⎜⎜⎛​∂ω1​∂f​∂ω2​∂f​…∂ωn​∂f​​⎠⎟⎟⎟⎞​=m2​⎝⎜⎜⎜⎛​∑i=1m​(X1(i)​ω1​+X2(i)​ω2​+…+Xn(i)​ωn​)X1(i)​∑i=1m​(X1(i)​ω1​+X2(i)​ω2​+…+Xn(i)​ωn​)X2(i)​…∑i=1m​(X1(i)​ω1​+X2(i)​ω2​+…+Xn(i)​ωn​)Xn(i)​​⎠⎟⎟⎟⎞​ 对上式进行合并，写成两个向量点乘的形式。更进一步对表达式进行向量化处理 ᐁ f = 2 m ( ∑ i = 1 m ( X ( i ) ω ) X 1 ( i ) ∑ i = 1 m ( X ( i ) ω ) X 2 ( i ) … ∑ i = 1 m ( X ( i ) ω ) X n ( i ) ) ᐁf=\frac {2}{m}\begin{pmatrix} \sum_{i=1}^m(X^{(i)}ω)X_1^{(i)}\\ \sum_{i=1}^m(X^{(i)}ω)X_2^{(i)}\\ …\\ \sum_{i=1}^m(X^{(i)}ω)X_n^{(i)}\end{pmatrix} ᐁf=m2​⎝⎜⎜⎜⎛​∑i=1m​(X(i)ω)X1(i)​∑i=1m​(X(i)ω)X2(i)​…∑i=1m​(X(i)ω)Xn(i)​​⎠⎟⎟⎟⎞​ ᐁ f = 2 m ( X ( 1 ) ω , … ， X ( m ) ω ) ⋅ ( X 1 ( 1 ) , X 2 ( 1 ) , … , X n ( 1 ) X 1 ( 2 ) , X 2 ( 2 ) , … , X n ( 2 ) … X 1 ( m ) , X 2 ( m ) , … , X n ( m ) ) ᐁf=\frac {2}{m}(X^{(1)}ω,…，X^{(m)}ω)\cdot \begin{pmatrix} X_1^{(1)},X_2^{(1)},…,X_n^{(1)}\\ X_1^{(2)},X_2^{(2)},…,X_n^{(2)}\\ …\\ X_1^{(m)},X_2^{(m)},…,X_n^{(m)}\end{pmatrix} ᐁf=m2​(X(1)ω,…，X(m)ω)⋅⎝⎜⎜⎜⎛​X1(1)​,X2(1)​,…,Xn(1)​X1(2)​,X2(2)​,…,Xn(2)​…X1(m)​,X2(m)​,…,Xn(m)​​⎠⎟⎟⎟⎞​ ᐁ f = 2 m ⋅ ( X ω ) T ⋅ X = 2 m ⋅ X T ⋅ ( X ω ) ᐁf=\frac {2}{m}\cdot (Xω)^T\cdot X=\frac {2}{m}\cdot X^T\cdot (Xω) ᐁf=m2​⋅(Xω)T⋅X=m2​⋅XT⋅(Xω) 得到梯度为： ᐁ f = 2 m ⋅ X T ⋅ ( X ω ) ᐁf=\frac {2}{m}\cdot X^T\cdot (Xω) ᐁf=m2​⋅XT⋅(Xω) 有了梯度，就可以使用梯度上升法求最大值了。

我们已经知道了PCA算法的目标函数（方差），以及如何使其最大（梯度上升）。下面使用代码来实现公式及算法流程。

下面我们以二维数据为例，将其映射到一维，即求出一系列样本点的第一主成分。

首先准备数据：

接下来对数据进行主成分分析，第一步是均值归零，定义相应的函数。

将样本进行均值归0（demean），即所有样本将去样本的均值。样本的分布没有改变，只是将坐标轴进行了移动。

其方法是：X表示一个矩阵，每一行代表一个样本，每一个样本的每一个特征都减去这个特征的均值。实现如下：

接下来就是对目标（方差）函数和梯度（导数）函数的定义。

首先定义目标函数：

然后根据梯度公式 ᐁ f = 2 m ⋅ X T ⋅ ( X ω ) ᐁf=\frac {2}{m}\cdot X^T\cdot (Xω) ᐁf=m2​⋅XT⋅(Xω)求梯度：

现在对其稍加改造，可以验证我们之前计算的梯度表达式的结果是否正确。

在这里，有一处需要注意的地方：

对于要求的轴，向量ω，实际上一个单位向量，即要让向量的模为1 因此我们使用np中的线性代数库linalg里面的norm()方法，求第二范数，也就是求算术平方根

然后就实现梯度上升代码的流程了：gradient_ascent()方法为梯度上升法的过程，在n_iters次循环中，每次获得新的ω，如果新的ω和旧的ω对应的函数f的值的差值满足精度epsilon，则表示找到了合适的 w w w：

此处要注意，对于每一次计算向量ω之前都要进行单位化计算，使其模为1。