03.PCA主成分分析

参考：
讲一下numpy的矩阵特征值分解与奇异值分解
主成分分析（PCA）原理详解
源代码： https://github.com/wucng/MLAndDL

---

1.特征值分解

(1) 特征值与特征向量

如果一个向量$v$是矩阵$A$的特征向量，将一定可以表示成下面的形式：

$Av=\lambda v$

其中，$\lambda$是特征向量$v$对应的特征值，一个矩阵的一组特征向量是一组正交向量。

(2) 特征值分解矩阵

从上面我们知道了特征矩阵与特征向量的定义，即满足：$Av=\lambda v$

如果存在矩阵$A\in R^{n\mathrm{x}n}$ 有n个特征值分别为 ${\lambda }_1,{\lambda }_2,...{\lambda }_n$，对应的特征向量分别为：$v_1,v_2,...v_n$

把所有的特征向量写成矩阵形式：$V\in R^{n\mathrm{x}n}$，特征值写成对角阵 $D\in R^{n\mathrm{x}n}$，则有：$AV=VD$

如果特征矩阵$V$ 是可逆的 则有 $A=VD{V}^{-1}=VD{V}^T$ ，因为特征矩阵其实是正交矩阵 即：$V{V}^T=E$，因此 $V^{-1}=V^T$

如果矩阵$A$是正定的(可逆) 则有 $A=VD{V}^T$

import numpy as np
from numpy.linalg import eig,svd

A = np.random.random_sample([8,8])
C = np.dot(A.T, A)
# 特征值分解
vals, vecs = eig(C)

# 重构
Lambda = np.diag(vals) # 特征值对角阵
new_C = np.dot(np.dot(vecs, Lambda), vecs.T) # 与C=A.T*A相等

print(np.allclose(new_C,C)) # True

2.SVD(奇异值)分解

如果矩阵$A$可逆，则可以进行特征值分解；如果矩阵$A$不可逆，则可以是SVD分解（并且$A$可逆也能使用SVD分解）

奇异值分解是一个能适用于任意矩阵的一种分解的方法，对于任意矩阵A总是存在一个奇异值分解：

$A=U\Sigma V^T$

假设$A$是一个$m\ast n$的矩阵，那么得到的$U$是一个$m\ast m$的方阵，$U$里面的正交向量被称为左奇异向量。$\Sigma$是一个$m\ast n$的矩阵，$\Sigma$除了对角线其它元素都为0，对角线上的元素称为奇异值。是$V$的转置矩阵，是一个$n\ast n$的矩阵，它里面的正交向量被称为右奇异值向量。而且一般来讲，我们会将$\Sigma$上的值按从大到小的顺序排列。

SVD分解矩阵A的步骤：

• (1) 求$A{A}^T$的特征值和特征向量，用单位化的特征向量构成 $U$。
• 求$A^{T}A$的特征值和特征向量，用单位化的特征向量构成 $V$。
• 将$A{A}^T$或者$A^{T}A$的特征值求平方根，然后构成 $\Sigma$。

import numpy as np
from numpy.linalg import eig,svd

A = np.random.random_sample([8,6])
u, s, vh = np.linalg.svd(A)    # 这里vh为V的转置
print(u.shape,s.shape,vh.shape)
print(np.allclose(A, np.dot(u[:, :len(s)] * s, vh))) # True

3.PCA算法

• (1) 基于特征值分解协方差矩阵实现PCA算法

输入：数据集 X = { x 1 , x 2 , . . . , x n } ， X ∈ R n ∗ m X=\{x_1,x_2,...,x_n\}，X \in R^ {n*m} X={x1​,x2​,...,xn​}，X∈Rn∗m，需要降到$k$维 $X^{\prime }\in R^{n\ast k},k<m$。

1. 去平均值(即去中心化)，即每一位特征减去各自的平均值。

2. 计算协方差矩阵$X^{T}X$,注：这里除或不除样本数量$n$或$n-1$,其实对求出的特征向量没有影响。

3. 用特征值分解方法求协方差矩阵$X^{T}X$的特征值与特征向量。

4. 对特征值从大到小排序，选择其中最大的$k$个。然后将其对应的$k$个特征向量(对应的列向量)分别作为列向量组成特征向量矩阵$P$。

5. 将数据转换到$k$个特征向量构建的新空间中，即$Y=XP$。

---

• (2) 基于SVD分解协方差矩阵实现PCA算法

输入：数据集 X = { x 1 , x 2 , . . . , x n } ， X ∈ R n ∗ m X=\{x_1,x_2,...,x_n\}，X \in R^ {n*m} X={x1​,x2​,...,xn​}，X∈Rn∗m，需要降到$k$维 $X^{\prime }\in R^{n\ast k},k<m$。

1. 去平均值，即每一位特征减去各自的平均值。

2. 计算协方差矩阵$X^{T}X$。

3. 通过SVD计算协方差矩阵的特征值与特征向量。

4. 对特征值从大到小排序，选择其中最大的k个。然后将其对应的k个特征向量（右奇异矩阵）分别作为列向量组成特征向量矩阵$P$。(左奇异矩阵可以用于对行数的压缩；右奇异矩阵可以用于对列(即特征维度)的压缩。)

5. 将数据转换到k个特征向量构建的新空间中,即$Y=XP$。

"""
Author:wucng
Time:  20200109
Summary: PCA 基于特征值分解
源代码： https://github.com/wucng/MLAndDL
"""
import scipy
# from scipy.misc import imread,imshow
from imageio import imread,imsave
import numpy as np
from numpy.linalg import eig,svd
import matplotlib.pyplot as plt
from sklearn.decomposition  import PCA

class PCASelf(object):
    def __init__(self,n_components=1,mode="svd"):
        """
        :param n_components:  主成分数（压缩后的特征数）
        :param mode: "eig" 特征矩阵分解，"svd" 奇异值分解
        """
        self.n_components = n_components
        self.mode = mode

    def fit_transform(self,X:np.array):
        # 去平均值
        X = X-np.mean(X,0)
        # 协方差矩阵
        A = np.matmul(X.T,X)#/len(X)

        if self.mode == "eig":
            # 计算协方差矩阵的特征值与特征向量
            vals, vecs = eig(A)
        else:
            # u, s, vh = np.linalg.svd(A)    # 这里vh为V的转置
            _, vals, vecs = svd(A)

        # 对特征值从大到小排序，选择其中最大的k个
        index = np.argsort(vals*(-1))[:self.n_components] # 默认是从小到大排序，乘上-1 后就变成从大到小排序

        # 根据选择的K个特征值组成新的特征向量矩阵（列对应特征向量，而不是行）
        P = vecs[:,index]

        # 特征压缩后的矩阵
        return np.matmul(X,P)

if __name__=="__main__":
    # 自定义方法
    X = np.array([[-1, 1], [-2, -1], [-3, -2], [1, 1], [2, 1], [3, 2]])
    pca = PCASelf(1)
    print(pca.fit_transform(X))

    # sklearn 方法
    print(PCA(1).fit_transform(X))