在许多研究领域当中，通常都需要对多个属性进行描述刻画，收集大量的数据并找到他们的分布规律和内在联系。每个属性可能彼此独立，也可能彼此关联。有些属性可能发挥的作用更大，而有些属性又不会起作用。因此如果只是传统地单独对每个属性进行分析，可能会得到完全错误的结论。
因此我们需要寻找到一种合理的方法，去减少所需要分析的属性——即去掉大量冗余的属性或作用不大的属性，但同时尽量减少原对象信息的丢失，这也是我们常说的“数据降维”。
而 PCA算法，即主成分分析，就可以做到我们想要的数据降维。将高维度中发挥作用的数据保留下来，去掉噪音，从而达到快速计算的目的。
接下来我们就会通过 PCA 在人脸识别上的应用来解释这个算法。
主要思路
其实具体的思路非常简单。此处我会以 olivetti faces 人脸数据集来作为例子，在这个数据集当中有400张人脸，共40人。
数据地址：https://cs.nyu.edu/~roweis/data/olivettifaces.gif
对于单张的图像来说，“属性”实际上就是组成每个图片的像素点。对于 olivetti faces，每张人脸的大小为47 * 57，即一共有 2679 个属性。
在这 2679 个属性当中，实际上也是有很多冗余的噪音存在的。我们需要将这 2679 维特征映射到k维上，这 k 维是全新的正交特征也被称为主成分，是在原有 2679 维特征的基础上重新构造出来的 k 维特征。PCA 的工作就是从原始的空间中顺序地找一组相互正交的坐标轴，新的坐标轴的选择与数据本身是密切相关的。这相当于只保留包含绝大部分方差的维度特征，而忽略包含方差几乎为 0 的特征维度，实现对数据特征的降维处理。
PCA 算法原理
那么问题来了：我们如何寻找到这些主成分向量呢？
实际上，只要通过计算数据的协方差矩阵，再寻找协方差矩阵的特征值和特征向量，选择特征值大的 k 个特征向量构建新矩阵，再将原始图像数据转换到由新矩阵构建的新空间当中，从而实现数据降维。
计算协方差矩阵有很多种方法，一般都是使用特征值分解矩阵和奇异值分解矩阵（SVD）。SVD 相较于前者，性能更好，也是普遍使用的算法（例如 sklearn 的 pca 函数就是使用 svd 算法实现的）。
协方差矩阵的特征向量分解
我们可以通过下面的公式来求得原始矩阵 X 的协方差矩阵 C：
C=\frac{X^{T}X}{n-1}
接着对协方差矩阵 C 进行特征值分解：
C=W \lambda W^{-1}
其中 \lambda 为由特征值组成的对角阵，W 为特征向量组成的矩阵。
实际上，在 numpy 当中，可以直接使用 np.linalg.eig(C) 来求解特征值与特征向量。
最后是将原始数据 X 导入到新的向量空间当中：
X_k=XW_k
SVD 奇异值分解
当数据量很大的时候，也会带来很多问题：
- 求解协方差矩阵以及求解其特征向量，都是一个计算量非常大的工作。
- 因为原始矩阵与其转置相乘，如果矩阵当中存在非常小的数，也会在平方当中丢失。
- 只能在方阵上使用
因此我们需要一个更加稳健的算法，也就是我们的SVD算法。
奇异值分解是一个能适用于任意矩阵的一种分解的方法，对于任意矩阵A总是存在一个奇异值分解：
A=U \Sigma V^{T}
其中 ∑ 为奇异值。
由此，我们根据协方差公式，将 X 进行奇异值展开：
C=\frac{X^{T}X}{n-1}= \frac{V\Sigma U^{T}U\Sigma V^{T}}{n-1}=\frac{V\Sigma^2 V^T}{n-1}
而由于 V 是一个酉矩阵（即V的转置等于V的逆矩阵），因此
C=V\frac{\Sigma ^2}{n-1}V^{-1}
我们在这里与这个式子对比一下：C=W \lambda W^{-1}
就会发现，实际上 \lambda 是与原始矩阵 X 的奇异值是有关联的：
\lambda=\frac{\Sigma ^2}{n-1}
而 C 的特征向量，就是 X 的奇异值分解之后的矩阵 V。因此之后我们只需要 X 右乘 V 即可。
在 SVD 当中，只需要对原始矩阵做奇异值分解即可，因此大大减少了计算量。
在人脸识别当中使用 PCA
对于我们使用的人脸数据，我们需要将所有灰度数值排成一列组成一个数组，这样有助于我们之后的分析。因此，我们就会得到一个 2679 * 400 的矩阵（2679 为图像的长乘高），每一列代表一个人脸信息。
数据预处理
在应用 PCA 之前，我们需要对数据进行归一化。
我们需要求解出一个“平均脸”，实际上就是对每个属性求平均值。然后求解每个数据与均值的距离。
# data 即为人脸数据
expec = np.mean(data, axis=0)
expec = np.reshape(expec, (1, len(expec)))
data = data - expec
我们可以来观察一下求得的平均脸
PCA
大家可以直接调用 sklearn 的 PCA() 函数来获取 PCA。
不过我在这里使用了 numpy 的 SVD 函数，但实际上二者是一样的。
# svd 计算并求取空间向量
u, s, vh = np.linalg.svd(temp_data, full_matrices=False)
eig_vecs = temp_data @ vh.T # 直接右乘 v，由于 vh 为 v 的转置因此需要转置回来
# 选出前 k 个数据
features = eig_vecs[:, :k]
return features, expec, vh
这里大家可以关注一下 k 值的选择。
numpy 得 SVD 函数默认已经对奇异值进行了从大到小的排列，而每一个奇异值则代表了对应特征向量表达数据特征的能力。因此我们只需求解前 k 个奇异值的平方和在总奇异值平方和当中的占比，即可解决 k 的取值问题。例如我需要让我之后的向量空间能表达 90% 的原始数据，则只需让前 k 个奇异值平方和在总奇异值平方和的占比达到 90% 即可。在这组数据里我选择的 k=90，相对于以前的 2679 个属性，可谓是大大减少了。
人脸匹配
获得了 PCA 之后，我们也获得了降维后的向量空间。
之后对测试样本要进行匹配，则只需将我们的测试样本投影到降维后的向量空间，再求取与各个训练样本的欧氏距离，欧氏距离最短的正是与之相匹配的。
这里与 KNN 算法类似，我们也可以寻找欧氏距离最近的几组数据，从中以所占个数最多的训练样本作为匹配结果。不过我这里没有使用。
这里随机从数据中抽取一个图像后显示结果。
# 抽取测试样本，组建训练样本
i = np.random.randint(0, data.shape[1])
training_data = np.concatenate((data[:, :i], data[:, (i+1):]), axis=1)
test_data = data[:, i:i+1].T
features, expec, vh = pca(training_data) # 将上面 pca 的代码打包后组成 pca 函数
eig_face = np.array((test_data - expec) @ vh.T)[:, :k] # 将测试数据放入降维的向量空间当中
# 求解测试数据与训练集各个样本的欧氏距离
dist = np.sqrt(np.sum(np.square(eig_face - features), axis=1))
# 返回找到的数据
ins = list(dist).index(np.min(dist))
res_face = training_data[:, ins]
# 画图
fig, ax = plt.subplots(1, 2, figsize=(4, 4), sharex=True, sharey=True)
ax[0].set_title("Test Sample")
ax[1].set_title("Result Sample")
displayFace(test_data, ax[0])
displayFace(res_face, ax[1])
运行结果：
交叉验证
由于数据量并没有特别大，我选择使用交叉验证来对模型进行评估
res_list = []
for i in range(400):
training_data = np.concatenate((data[:, :i], data[:, (i+1):]), axis=1)
test_data = data[:, i:i+1].T
training_label = label[:i] + label[i+1:]
test_label = label[i]
features, expec, vh = pca(training_data)
eig_face = np.array((test_data - expec) @ vh.T)[:, :k]
dist = np.sqrt(np.sum(np.square(eig_face - features), axis=1))
ins = list(dist).index(np.min(dist))
res_list.append(1 if training_label[ins] == test_label else 0)
最后所得到的准确度 ACC = 0.98，使用 kNN 可能可以得到更好的数据
PCA 算法所存在的缺陷
PCA 算法极度依赖于所训练出的数据集，也正因如此，一旦需要向数据集当中添加数据，PCA 算法又需要重新对所有数据进行训练来生成新的向量空间。
在这一点上，PCA 与 KNN 是一致的，只不过 PCA 更加快速。
所以目前商用的人脸识别大多采用更加复杂的卷积神经网络来解决上述问题。
参考资料
数字图像处理（第三版） 胡学龙
https://zhuanlan.zhihu.com/p/37777074
Comments | NOTHING