西瓜书学习笔记——主成分分析（公式推导+举例应用）

最新推荐文章于 2024-04-24 01:54:40 发布

Nie同学

最新推荐文章于 2024-04-24 01:54:40 发布

阅读量1.9k

点赞数 31

分类专栏：机器学习文章标签：学习笔记机器学习降维

本文链接：https://blog.csdn.net/qq_46117575/article/details/136000619

版权

机器学习专栏收录该内容

16 篇文章 8 订阅

订阅专栏

文章目录

- - 算法介绍
  - 实验分析

算法介绍

主成分分析（Principal Component Analysis，PCA）是一种常用的降维技术，用于在高维数据中发现最重要的特征或主成分。PCA的目标是通过线性变换将原始数据转换成一组新的特征，这些新特征被称为主成分，它们是原始特征的线性组合。

对于一个正交属性空间（各个属性之间是线性无关的）中的样本点，存在以下两个性质的超平面可对所有样本点进行恰当的表达：

最近重构性：样本点到这个超平面的距离足够的近。
最大可分性：样本点在这个超平面上的投影尽可能分开。

我们先从最近重构性的角度来进行推导：

假定数据样本已进行中心化，即 $\sum_{i}x_i=0$ ；
假定投影变换后得到的新坐标系为 ${w_1,2_2,...,w_d\}$ ，其中 $w_i$ 是标准正交基。
- $w_i||_2=1$
- $w_i^Tw_j=0(i \ne j)$

若丢弃新坐标系中的部分坐标，即将维度降低到 $d^\prime<d$ ，则样本点 $x_i$ 在低维坐标系下的投影是 $z_i=(z_{i1};z_{i2};...;z_{id^\prime})$ ，其中 $z_{ij}=w_j^Tx_i$ 是 $x_i$ 在低维坐标系下第 $j$ 维的坐标。若基于 $z_i$ 来重构 $x_i$ ，则会得到 $\hat x_i=\sum_{j=1}^{d^\prime}z_{ij}w_j$ 。

考虑到整个训练集，原样本点 $x_i$ 与基于投影重构的样本点 $\hat x_i$ 之间的距离为：
$\begin{aligned} \sum_{i=1}^m||\sum_{j=1}^{d^\prime}z_{ij}w_j-x_i||_2^2&=\sum_{i=1}^m\sum_{j=1}^{d^\prime}(z_{ij}w_j-x_i)^T(z_{ij}w_j-x_i)\\ &=\sum_{i=1}^m\sum_{j=1}^{d^\prime}(z_{ij}w_j^T-x_i^T)(z_{ij}w_j-x_i)\\ &=\sum_{i=1}^m\bigg((\sum_{j=1}^{d^\prime}z_{ij}w_j^T)(\sum_{j=1}^{d^\prime}z_{ij}w_j)-2\sum_{j=1}^{d^\prime}z_{ij}w_j^Tx_i+x_i^Tx_i\bigg) \end{aligned}\tag{1}$

其中 $w$ 是标准正交基，故有 $w_i||_2=1$ 和 $w_i^Tw_j=0(i \ne j)$ ，且 $x_i^Tx_i$ 是一个常数，令 $x_i^Tx_i=\text{const}$ ，所以（1）式可以简化为：
$\begin{aligned} \sum_{i=1}^m||\sum_{j=1}^{d^\prime}z_{ij}w_j-x_i||_2^2&=\sum_{i=1}^m\bigg((z_{i1}w_1^T+z_{i2}w_2^T+...+z_{id^\prime}w_{d^\prime}^T)(z_{i1}w_1+z_{i2}w_2+...+z_{id^\prime}w_{d^\prime})-2z_i^T\begin{bmatrix} w_1^T \\ w_2^T \\ ...\\ w_{d^\prime}^T \\ \end{bmatrix}x_i\bigg)+\text{const}\\ &=\sum_{i=1}^m\bigg(\sum_{j=1}^{d^\prime}z_{ij}^Tz_{ij}-2z_i^Tw^Tx_i\bigg)+\text{const}\\ &=\sum_{i=1}^m\bigg(z_i^Tz_i-2z_i^Tz_i\bigg)+\text{const}\\ &=-\sum_{i=1}^m z_i^Tz_i+\text{const}\\ &=-(z_1^Tz_1+z_2^Tz_2+...+z_m^Tz_m)+\text{const}\\ &=-\sum_{i=1}^mtr(z_i^Tz_i)+\text{const} \end{aligned}\tag{2}$

又因为 $t r (A B) = t r (B A)$ 故（2）式可变为：
$\begin{aligned} \sum_{i=1}^m||\sum_{j=1}^{d^\prime}z_{ij}w_j-x_i||_2^2&=-\sum_{i=1}^mtr(z_iz_i^T)+\text{const}\\ &=-tr\bigg(\sum_{i=1}^m z_iz_i^T\bigg)+\text{const}\\ &=-tr\bigg(\sum_{i=1}^mw^Tx_ix_i^Tw\bigg)+\text{const}\\ &=-tr\bigg(w^T(\sum_{i=1}^mx_ix_i^T)w\bigg)+\text{const}\\ &\propto-tr\bigg(w^T(\sum_{i=1}^mx_ix_i^T)w\bigg) \end{aligned}\tag{3}$

根据最近重构性，式（3）应该被最小化，有：

$\begin{aligned} & \min _{\mathbf{W}}-\operatorname{tr}\left(\mathbf{W}^{\mathrm{T}} \mathbf{X} \mathbf{X}^{\mathrm{T}} \mathbf{W}\right) \\ & \text { s.t. } \mathbf{W}^{\mathrm{T}} \mathbf{W}=\mathbf{E} . \end{aligned} \tag{4}$

这就是PCA的最终优化目标。

接着我们从最大可分性出发，我们知道样本点 $x_i$ 在新空间中的超平面上的投影是 $W^Tx_i$ ，若所有的样本点的投影能尽可能的分开，则应该使投影后样本点的方差最大化。

在降维和投影中，我们通常不使用传统的方差定义，因为这里我们不涉及到随机变量和随机过程。相反，我们使用方差的术语来描述投影后样本点在新空间中的变异性，这更多地是一种几何上的直觉而非统计学中方差的定义。

方差的形式为 $x_i^T W W^T x_i$ ，可以从以下几何角度解释：

投影方向：投影矩阵 $W$ 确定了投影的方向，每一列表示新空间的基向量。
样本点变异性： $x_i$ 是原始空间中的样本点， $W^Tx_i$ 是投影到新空间后的样本点。 $x_i$ 和 $W^Tx_i$ 之间的差异，即 $x_i - W^Tx_i$ ，表示样本点在投影方向上的变异性。
方差形式：方差的形式 $x_i^T W W^T x_i$ 表示了在投影方向上的样本点变异性的量度。这个形式的计算实际上考虑了样本点在投影方向上的平方和，类似于方差的直觉。
最大化分离性：当我们尝试最大化这个“方差”时，我们实际上是在尝试选择一个投影方向，使得在新空间中的样本点在投影方向上的变异性最大化。这有助于确保在新空间中，样本点能够尽可能地分开，提高了样本点在新空间中的可分性。

请注意，这种形式的“方差”并不是统计学中方差的定义，而更多是一种在降维和投影问题中用于描述样本点分布情况的术语。这种定义的目标是为了在新空间中找到一个有利于分离的方向，以便更好地实现降维和分类等任务。

故投影后样本点的方差是

$\begin{aligned} Var(W^Tx_i)&=Var(x_i^TW)\\ &=\sum_iW^Tx_ix_i^TW\\ &=W^T(\sum_ix_ix_i^T)W\\ &=W^T\bigg( x_1x_1^T+x_2x_2^T+...+x_mx_m^T\bigg)W\\ &=W^T\bigg( \sum_itr(x_ix_i^T)\bigg)W\\ &=tr(W^TXX^TW) \end{aligned}\tag{5}$

于是优化目标可写成：
$\begin{aligned} \max _{\mathbf{W}} & \operatorname{tr}\left(\mathbf{W}^{\mathrm{T}} \mathbf{X} \mathbf{X}^{\mathrm{T}} \mathbf{W}\right) \\ \text { s.t. } & \mathbf{W}^{\mathrm{T}} \mathbf{W}=\mathbf{E}, \end{aligned} \tag{6}$

显然式（4）和式（6）等价，故对式（4）或（6）使用拉格朗日乘子法，引入拉格朗日乘子矩阵 $\Lambda$ ，我们的目标是：
$\mathcal{L}(\mathbf{W,\Lambda)}=tr(\mathbf{W^TXX^TW)}-tr(\mathbf{\Lambda(W^TW-E))} \tag{7}$

现在我们分别对 $W$ 和 $\Lambda$ 求偏导，并令其为零有：

$\frac{\partial \mathcal{L}}{\partial \mathbf{W}}=2 \mathbf{X} \mathbf{X}^{\mathrm{T}} \mathbf{W}-2 \Lambda^{\mathrm{T}} \mathbf{W}=0 \tag{8}$

$\frac{\partial \mathcal{L}}{\partial \boldsymbol{\Lambda}}=-\left(\mathbf{W}^{\mathrm{T}} \mathbf{W}-\mathbf{E}\right)=0 \tag{9}$

最终得到：

$\mathbf{XX^TW}=\mathbf{\Lambda^TW} \tag{10}$

$\mathbf{W^TW=E} \tag{11}$

其中 $\Lambda=diag(\lambda_1,\lambda_2,...,\lambda_d^\prime) \in \mathbb{R}^{d^\prime\times d^\prime}$ 。

将（10）式展开有：

$\mathbf{XX^T}w_i=\lambda_iw_i\tag{12}$

最后我们只需要对 $\mathbf{XX^T}$ 进行特征值分解，其中 $\lambda_i$ 是矩阵 $\mathbf{XX^T}$ 的特征值， $w_i$ 是 $\lambda_i$ 所对应的单位特征向量。

最终我们选取前 $d^\prime$ 个特征值所对应的特征向量构成的 $\mathbf{W}=(w_1,w_2,...,w_{d^\prime})$ ，这就是我们所求的最终解。

其PCA算法流程图如下图所示：

在这里插入图片描述

这就出现了一个新的问题，低维空间维数 $d^\prime$ 是由用户事先决定的，那么该如何确定 $d^\prime$ ？

我们可以从重构的角度来设置一个重构阈值，例如 $t=95\%$ ，然后选取使得下式成立的最小 $d^\prime$ 值：

$\frac{\sum_{i=1}^{d^{\prime}} \lambda_i}{\sum_{i=1}^d \lambda_i} \geqslant t$

对于式（8）和（9）求偏导过程如下：

$\begin{aligned} \frac{\partial \mathcal{L}}{\partial \mathbf{W}}&=\frac{\partial\bigg( \ tr(\mathbf{W^TXX^TW})-tr(\mathbf{\Lambda(W^TW-E))}\bigg)}{\partial \mathbf{W}}\\ &=\frac{\partial \ tr(\mathbf{W^TXX^TW})}{\partial \mathbf{W}}-\frac{\partial \ tr(\mathbf{\Lambda(W^TW-E))}}{\partial \mathbf{W}} \end{aligned} \tag{13}$

对式（13）我们可以用迹的性质， $t r (A B) = t r (B A)$ ，故我们可将 $\mathbf{W^TXX^TW}$ 和 $\mathbf{\Lambda(W^TW-E)}$ 重新排列成 $\mathbf{WW^TX^TX}$ 和 $\mathbf{\Lambda(WW^T-E)}$ ，即：

$\begin{aligned} \frac{\partial \mathcal{L}}{\partial \mathbf{W}}&=\frac{\partial \ tr(\mathbf{WW^TX^TX})}{\partial \mathbf{W}}-\frac{\partial \ tr(\mathbf{\Lambda(WW^T-E))}}{\partial \mathbf{W}}\\ &=\frac{\partial\ tr(\mathbf{W W^T X^T X)}}{\partial \mathbf{(W W^T)}} \cdot \frac{\partial\mathbf{(W W^T)}}{\partial \mathbf{W}}-2\mathbf{\Lambda W}\\ &=\mathbf{X^TX}\cdot 2\mathbf{W}-2\mathbf{\Lambda W}\\ &=2 \mathbf{X} \mathbf{X}^{\mathrm{T}} \mathbf{W}-2 \Lambda^{\mathrm{T}} \mathbf{W} \end{aligned}\tag{14}$

对于式（7）约束条件为什么可以写成迹的形式？

在这个问题中，我们要证明在优化问题中 $tr(\mathbf{\Lambda(W^TW-E)})$ 和 $\mathbf{\Lambda(W^TW-E)}$ 是等价的，即它们对于最优解是相同的。下面我们来给出详细的证明！

为此，我们可以比较两者的梯度。考虑对拉格朗日函数 $\mathcal{L}(\mathbf{W,\Lambda})$ 是关于 $W$ 的偏导数。首先计算 $\mathcal{L}$ 对 $\mathbf{W}$ 的梯度得到式（10），解出最优解 $\mathbf{W^\star}$ 。现在我们比较两种形式对最优解的梯度。

对于 $tr(\mathbf{\Lambda(W^TW-E)})$ ：
$\frac{\partial tr(\mathbf{\Lambda(W^T W-E))}}{\partial \mathbf{W}}=2 \Lambda \mathbf{W}$
对于 $\mathbf{\Lambda(W^TW-E)}$ ：
$\frac{\partial \mathbf{\Lambda(W^T W-E)}}{\partial \mathbf{W}}=2 \Lambda \mathbf{W}$

从上面的梯度计算中，我们可以看到，两种形式的梯度在最优解处是相同的，因此在优化问题中 $tr(\mathbf{\Lambda(W^TW-E)})$ 和 $\mathbf{\Lambda(W^TW-E)}$ 是等价。这意味着它们会产生相同的最优解。

我们将矩阵运算转化为迹（trace）形式在某些情况下具有优势，尤其是在涉及到矩阵的导数计算和优化问题时。以下是一些使用迹的好处：

简化导数计算：迹运算的导数计算相对简单，尤其是对于迹的线性性质，这可以简化复杂矩阵函数的导数计算。这种简化对于深度学习等领域中的梯度下降等优化问题很有用。
矩阵的不可交换性：在矩阵乘法中，矩阵的乘法是不可交换的，即 $A B$ 一般不等于 $B A$ 。但是，通过转化为迹形式，有时可以使计算更加简单。
整合约束条件：在优化问题中，通过引入拉格朗日乘子法，迹运算可以用于整合约束条件，使得问题的形式更加一致，更容易处理。
向量化和矩阵形式：将矩阵运算表示为迹形式有助于使用向量化和矩阵形式，从而更好地利用硬件加速（如GPU）来加速计算。

实验分析

数据集如下图所示：
在这里插入图片描述

读入数据集：

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt

# 读取数据
data = pd.read_csv('data/correlated_dataset.csv')

# 提取特征列
features = data.drop('Target', axis=1)

计算协方差矩阵以及特征值分解：

# 计算协方差矩阵
cov_matrix = features.cov()

# 特征值分解
eigenvalues, eigenvectors = np.linalg.eig(cov_matrix)

选择最小 $d^\prime$ 值：

# 选择重构阈值
reconstruction_threshold = 0.95

# 计算特征值的累积和
eigenvalue_cumsum = np.cumsum(eigenvalues) / np.sum(eigenvalues)

# 找到满足阈值条件的最小d'
selected_d_prime = np.argmax(eigenvalue_cumsum >= reconstruction_threshold) + 1

# 输出选择的最小d'
print("选择的最小降维维度 d' =", selected_d_prime)

选择的最小降维维度 d' = 7

重新计算特征值和特征向量：

# 根据选择的最小d'重新计算选取的特征值和特征向量
selected_indices_reconstructed = np.argsort(eigenvalues)[-selected_d_prime:][::-1]
selected_eigenvalues_reconstructed = eigenvalues[selected_indices_reconstructed]
selected_eigenvectors_reconstructed = eigenvectors[:, selected_indices_reconstructed]

# 投影数据到新的空间
reduced_data_reconstructed = features.dot(selected_eigenvectors_reconstructed)

可视化数据降维后的数据：

plt.scatter(reduced_data_reconstructed.iloc[:, 0], reduced_data_reconstructed.iloc[:, 1], c=data['Target'], cmap='viridis')
plt.xlabel('Principal Component 1')
plt.ylabel('Principal Component 2')
plt.title('PCA of Correlated Dataset (Reconstructed)')
plt.show()