截断SVD

# 截断奇异值分解 来源：https://brashandplucky.com/2023/09/09/truncated-svd.html ## 主成分分析 *主成分分析*（PCA (https://en.wikipedia.org/wiki/Principal_component_analysis)）是我之前一篇博客的主题 (https://brashandplucky.com/2023/09/07/principal-component-analysis.html)，因此这里仅做概述。 通过PCA进行数据降维，是将数据线性变换到一个新的坐标系中完成的，该坐标系中数据的大部分变化可以用比初始数据更少的维度来描述。 简而言之，这涉及到协方差矩阵 (https://en.wikipedia.org/wiki/Covariance_matrix) 的特征分解。 主成分分析 *奇异值分解*（SVD (https://en.wikipedia.org/wiki/Singular_value_decomposition)）是一种矩阵分解技术，它将一个实矩阵`M`分解为三个矩阵`U`、`Σ`和`V`，使得`M=U*Σ*V^T`。 如果`M`是`m×n`的，那么`U`是`m×m`的，`Σ`是`m×n`的，`V`是`n×n`的。`U`和`V`都是正交矩阵，而`Σ`是矩形对角矩阵，其对角线系数非负。 这与PCA非常相似，区别在于SVD的分解是在数据矩阵上进行的，而PCA的分解是在协方差矩阵上进行的。 `Σ`的对角线系数被称为`M`的*奇异值*，通常的做法是重新排列SVD，使奇异值按降序排列。非零奇异值的数量等于`M`的秩。 SVD与PCA密切相关 (https://intoli.com/blog/pca-and-svd/)： - `V`的列是主方向/主轴（特征向量）。 - `U*Σ`的列是主成分（得分）。 - 奇异值与协方差矩阵的特征值相关。 来自维基百科： SVD矩阵 ## 截断奇异值分解 让我们用一张1024x1024的月亮灰度图像来尝试： 截断SVD（月亮） 这样的图像可以解释为1024个向量，每个向量有1024个分量。*即*，一个1024维空间中的1024个向量的集合。 如果我们在这组数据上运行PCA/SVD，三个矩阵`U`、`Σ`、`V`将是1024x1024的。特别是，`Σ`将是一个方阵对角矩阵。*即*，只有对角线上的系数可能非零。通常的做法是重新排列这三个矩阵，使得`Σ`中的对角线索引从较大值（左上角）到较小值（右下角）排序。 *截断SVD*就是将`Σ`中除了左上角`n`个以外的所有系数置零的行为。 系数截断意味着我们同时丢弃了`U`和`V`中的`1024-n`列（因为这些列乘以0后结果也为0）。 如果我们现在使用截断后的矩阵重构原始矩阵`M'=U'*Σ'*V'^T`，我们将得到`M'`，它与`M`相似。我们舍弃的系数越少，`M'`近似`M`的程度就越高。但由于PCA/SVD的信息保留特性，只保留`Σ`中最靠前的一批系数可能就足以恢复全部（或接近全部）的原始信息。 ## 重构月亮图像 让我们将上述所有方法应用于月亮图像。在单个`n`值上截断并不能很好地说明问题，因此我们将`n`在其整个范围内变动，观察近似效果如何变化。 以下图像和视频共享相同的布局： - 左半部分是重构的矩阵`M'`。 - 右半部分是重构误差`abs(M-M')`。 - 递减的黄色图表是随着越来越多的奇异值被置零，均方误差（MSE）的变化情况。 截断SVD（月亮） 此截图是仅用32个（总共1024个）分量重构的`M'`。*右键点击 + 在新标签页中打开*以查看1:1质量。 值得停下来思考一下，这32个分量实际消耗了我们什么。截断后的因子是`U'`（1024x32）、存活下来的32个奇异值以及`V'`（1024x32）。这总共是65,568个数字，而原始矩阵`M`有1,048,576个数字。*即*，实现了16倍的压缩，而图像仍然清晰可辨是月亮。 下面是一个视频，展示了使用更多分量进行重构时的相同图像对。大部分变化发生在前几帧。 截断SVD序列（月亮） 这里值得注意的是（PCA/SVD的魔力）误差图表下降的速度之快。这证明了前面的分量捕获了`M`中存在的大部分信息，而后面的分量只携带高频/低振幅的细微细节。 这让人联想到傅里叶变换、余弦/正弦变换等情况。这些变换处理的是空间与频率的对偶性，而PCA/SVD纯粹是基于方差的基变换。但类似地，所有这些方法都将信息转换为一种对偶形式，使得“信息量”以一种结构化、可管理的方式呈现出来。 FT/CT等...是*.jpeg*、*.mp3*和其他压缩系统的基础，这些系统利用了人类感知系统对亮度（相对于色度）和低频（相对于高频）更敏感的事实。 在SVD/PCA的情况下，`Σ`中的上部系数比下部系数捕获了更多的数据方差。 这些序列是使用2、4、8、16、32、64和128个系数进行的重构。 截断SVD序列（小丑） ## 较易与较难的情况 如上所述，矩阵分解可以解释为基变换到一个特殊空间，在这个空间中，数据的行之间是线性相关的。在这种情况下，较低秩的`U/V`矩阵就足以精确重构原始矩阵。实际上，`Σ`的对角线上将出现与可以丢弃而不会造成任何数据损失的行/列数量一样多的零值系数。 下面展示了一个极端情况（一个居中的正方形），其中1个系数/行/列就足够了。在这种情况下，形状内部和外部的所有行/列都是相同的。顺便说一句，方框是一个可分离的卷积滤波器（希望以后有关于低秩卷积的文章）。 截断SVD（矩形） 旋转形状会带来灾难，尽管PCA/SVD能够“自动检测”这种基变换。但这里我们处理的是离散数学，因此旋转后的形状会被“像素化”，这使得分解在数值上变得不那么纯净。 有趣的是，在视频的前几帧中，重构“坚持”要呈现为一个未旋转的方形。 截断SVD（倾斜的矩形） 下面是一个五边形形状。 截断SVD（五边形） ## 实际应用 除了数据分析中的降维，SVD截断还有一些有趣的实际应用。 ### 数据压缩 Bart Wronski (https://bartwronski.com/2020/05/21/dimensionality-reduction-for-image-and-texture-set-compression/) 有一篇非常有趣的文章，关于使用此技术进行*PBR纹理集压缩*。 *BCn纹理压缩*也基于降维。Nathan Reed (https://www.reedbeta.com/blog/understanding-bcn-texture-compression-formats/) 对此主题有很好的阐述。 ### 低秩近似 低秩近似 (https://en.wikipedia.org/wiki/Low-rank_approximation)。 Bart Wronski的另一篇有趣文章 (https://bartwronski.com/2020/02/03/separate-your-filters-svd-and-low-rank-approximation-of-image-filters/)。我希望有朝一日能自己写一篇关于*低秩卷积*的文章。 ### 实现细节 我在Maverick (https://maverickrender.com/) 的API中有一些旧的PCA/SVD C++代码，我用它来生成本文中的图像/视频。但在阅读了Atrix256 (https://blog.demofox.org/2022/07/12/calculating-svd-and-pca-in-c/) 的这篇文章后，我可能会下定决心，用Eigen (https://eigen.tuxfamily.org/) 来替换我旧的`xsvd_c`类中的实现部分。