论文部分内容阅读
摘 要：本文用主成分分析 （PCA） 降维技术解大型超定线性方程组AX=B的最小二乘解。工业上遇到的大型方程组的求解要消耗大量的时间和内存， 系数矩阵通常是病态的， 会带来不可忽视的误差。本文设计基于PCA的算法， 并从理论上说明其可行性。实验证明该方法有效， 不仅测试误差极小， 接近原方程的误差， 而且计算时间显著减少。
关键词：PCA; SVD; 线性方程组; 最小二乘解
文章编号：2095-2163（2019）04-0091-05 中图分类号：TP399 文献标志码：A
0 引 言
很多工业问题最终会转化成解线性方程组问题。但是，通常工业级的方程组规模巨大， 直接求解会消耗大量的时间和内存。本文提出一种基于主成分分析 （PCA） [1-5]的快速降维算法求解大型线性方程组。
目前，印染行业急需解决的难题是在染色过程中， 因为事先不清楚染料光谱的精确数值， 只能先从多个布匹的染色流程中， 获取染色光谱数据， 再估计使用染料的光谱数据， 最后用于新的布匹染色。设建立一个线性模型：
N×n。问题是：使用了多种染料， 获取了大量光谱数据。方程规模非常巨大， 常规解法可能会消耗大量时间和内存。因此用机器学习的降维方法来缩小矩阵大小， 然后再解方程。
本文采用了5 000多组数据，（而这只是工业数据中极小的一部分）， 目的在于验证本文方法的有效性。
1 算法设计
解线性方程组AX=B， 其中矩阵A很大， 可能是病态的， 且一般是列满秩的 （A的行数远大于列数）， 因此可用ATA的条件数衡量方程病态程度。B有时也不只有一列。
由于各种原因， 原方程可能没有解， 只能寻求最小二乘解， 即解优化问题。
1.1 矩阵分析
为了避免多余的矩阵运算， 没有按照标准的PCA流程， 而是直接对ATA进行奇异值分解 （SVD） 或特征值分解：
由式 （6） 可知， 对误差来说， 选择哪些主成分似乎并不重要， 相对误差大致可以写成：
若对B降维，有如下公式：
其中，B1是对B的PCA重构。观察式 （7）， 显然应该选取对应特征值绝对值最大的主成分。
1.2 算法
根据上述分析设计如下算法（见表1），将一个大型方程组分解成2个小型方程组。
若使用NMF或其它分解技术， 分解A=WH， 则不能避免做最小二乘法。
理论上， 最小二乘解是方程误差最小的解， 但降维之后， 这个解就不再是原方程的最优解了。
1.3 解的处理
设最后得到的[AKX^]=V1Y1，是原方程AX=B的降维（最小二乘）解。工业中有时候默认所有量都是非负的， 但是降维解可能含有负数。为了满足非负性， 可将其修正为：
实际上， 为了方便计算， 还可对B进行降维，此时方程变为：
2 数值实验
本文算法用 Python实现， 运行于MacOS上。所用数据、源代码和运行结果都托管在Github 上， 网址为https：//github.com/Freakwill/PCA-linear-equations.
2.1 主成分选择
对矩阵A，B降维， 通过观察累计百分比曲线， 选择适当的主成分数。如图1所示。
2.2 误差分析
误差采用向量的2范数， 对矩阵来说就是Frobenius 范数。本文用相对误差公式衡量算法逼近能力。
20%的样本会被事先抽取出来， 测试误差通常比训练误差更有说服力。降维是对训练样本进行的， 而不是对测试样本进行。但对A的降维也可以应用于训练样本， 因为其是已知输入值， 而测试样本中的B作为目标值， 不参与降维。如图2所示。
由图2结果可知：
（1）原方程组没有解。图中原方程误差是理论上最小（在没有约束的情况下）， 不妨看做一个相对0点。如果降维方法的误差低于这个点， 可认为是系统误差造成的， 不是数值分析意义上的误差。
（2）降维方程组的误差衰减达到预期。“负数值置0”的处理， 对误差没有太大影响。当A的主成分数超过70时， 计算变得不稳定， 预测误差表现很夸张。
（3）从误差曲线不难看出， 解方程组时并不需要用到所有样本，样本数量对本算法没有太大影响。实际应用中， 随机挑选充分数量的样本即可。
图3是对B取前4个主成分得到的误差图， 和之前的实验相比，并没有显著影响误差。时间和主成分数基本上呈线性关系， 符合预期。为了获得较为准确的运行时间， 本次实验对每一个主成分数重复50次， 无论时间还是误差都取均值。
由图4结果可知：
（1）B主成分数对误差影响没有A大 （这是优点）， 由于B维数较低， 且主成分比较集中， 在第一个主成分处， 误差已降到0.3左右。
（2）用时与主成分数近似呈弱线性相关， B的降维也相应提高了效率。
总之， 对B降维是完全合理的。
可通过设定随机种子， 产生固定的训练-测试样本（实验重复50次）。
选用A前30个主成分，B前4个主成分。 获得实验结果見表2。
2.3 其它降维策略与拟合方法
NMF降维的效果和PCA相近， 但矩阵分解后需解较为复杂的方程。PCA的优势在于能分解出正交矩阵， 可以设计出更快捷的算法， 计算用时短。中心化PCA在主成分数较少时表现较好， 这是因为中心化PCA采用的是仿射变换， 比单纯的线性变换多了常数项， 但随着主成分数增加， 并没有表现出明显优势。目前中心化PCA只是简单重构A， 即： 其中，C1VT1是中心化后的分解， 即通常意义上的PCA， M是均值矩阵，其无法简化方程， 计算时间维持在某个常数。若能设计出简化方程的算法， 那么中心化PCA或许是不错的选择。
当主成分增加时， 负数解都是无法避免的。 数值实验表明主成分过多还有可能出现异常。几种求解方法的测试比较结果如图5所示。
显然，完全可以用其它方法求解降维后的方程组， 尤其是当人们只关心预测， 而不在乎X时。如果只为了建立A与B的联系， 完全可以采用一些非线性方法。表3中给出一些线性模型测试结果，用到了一些较复杂的计算策略， 含非线性成分， 误差无显著降低， 运行时间则较长。这些模型均由Python 第三方库scikit-learn实现[7-8]。
最后， 给出一般的计算框架， 可为解线性方程组开发新的算法：
3 结束语
PCA降维在解大型线性方程组表现的非常出色， 随着主成分增加， 误差快速递减。 这种简单的代数学技巧， 不仅计算快， 算法设计简单， 由于矩阵分解为对角矩阵和正交矩陣的乘积， 之后也无需解任何线性方程组。可以根据需要任意压缩原方程。
本文提供的方法可以应用于工业生产中。但在实验中发现过多的主成分可能使得算法不稳定。今后的工作重点：
（1）当变量有非负性或其它约束条件时， 本文还没有给出更合理的处理办法。
（2）如何实现有效的基于中心化PCA的算法。本文提到的中心化PCA没有真正起到压缩方程组的作用。
（3）应利用系数矩阵的一些特点， 如稀疏性， 设计更有针对性的算法。
参考文献
[1]HASTIE T， TIBSHIRANI R， FRIEDMAN J H. The elements of statistical learning：Data mining， inference， and prediction [M]. New York ：Springer-Verlag， 2001.
[2] HOTELLING H. Analysis of a complex of statistical variables into principal components[J]. Journal of Educational Psychology， 1993，24（6）：417-441.
[3] TURK M A， PENTLAND A P. Face recognition using eigenfaces[C]//Proceedings 1991 IEEE Computer Society Conference on ComputerVision and Pattern Recognition. Maui， HI， USA：IEEE，1991：586-591.
[4] SHALEV-SHWARTZ S， BEN-DAVID S. Understanding machine learning：from theory to algorithm[M]. New York， NY， USA ：Cambridge University Press， 2014.
[5] 焦李成，等. 稀疏学习、分类与识别[M]. 北京：科学出版社， 2017.
[6] 解锋昌， 韦博成， 林金官. 零过多数据的统计分析及其应用[M]. 北京：科学出版社， 2013.
[7] HACKELING G. Mastering machine learning with scikit-learn[M]. 2nd ed. Birmingham ：Packt Publishing， 2017.
[8] COURNAPEAU D. Scikit-learn[EB/OL]. [2019]. https：//scikit-learn.org/stable/.
关键词：PCA; SVD; 线性方程组; 最小二乘解
文章编号：2095-2163（2019）04-0091-05 中图分类号：TP399 文献标志码：A
0 引 言
很多工业问题最终会转化成解线性方程组问题。但是，通常工业级的方程组规模巨大， 直接求解会消耗大量的时间和内存。本文提出一种基于主成分分析 （PCA） [1-5]的快速降维算法求解大型线性方程组。
目前，印染行业急需解决的难题是在染色过程中， 因为事先不清楚染料光谱的精确数值， 只能先从多个布匹的染色流程中， 获取染色光谱数据， 再估计使用染料的光谱数据， 最后用于新的布匹染色。设建立一个线性模型：
N×n。问题是：使用了多种染料， 获取了大量光谱数据。方程规模非常巨大， 常规解法可能会消耗大量时间和内存。因此用机器学习的降维方法来缩小矩阵大小， 然后再解方程。
本文采用了5 000多组数据，（而这只是工业数据中极小的一部分）， 目的在于验证本文方法的有效性。
1 算法设计
解线性方程组AX=B， 其中矩阵A很大， 可能是病态的， 且一般是列满秩的 （A的行数远大于列数）， 因此可用ATA的条件数衡量方程病态程度。B有时也不只有一列。
由于各种原因， 原方程可能没有解， 只能寻求最小二乘解， 即解优化问题。
1.1 矩阵分析
为了避免多余的矩阵运算， 没有按照标准的PCA流程， 而是直接对ATA进行奇异值分解 （SVD） 或特征值分解：
由式 （6） 可知， 对误差来说， 选择哪些主成分似乎并不重要， 相对误差大致可以写成：
若对B降维，有如下公式：
其中，B1是对B的PCA重构。观察式 （7）， 显然应该选取对应特征值绝对值最大的主成分。
1.2 算法
根据上述分析设计如下算法（见表1），将一个大型方程组分解成2个小型方程组。
若使用NMF或其它分解技术， 分解A=WH， 则不能避免做最小二乘法。
理论上， 最小二乘解是方程误差最小的解， 但降维之后， 这个解就不再是原方程的最优解了。
1.3 解的处理
设最后得到的[AKX^]=V1Y1，是原方程AX=B的降维（最小二乘）解。工业中有时候默认所有量都是非负的， 但是降维解可能含有负数。为了满足非负性， 可将其修正为：
实际上， 为了方便计算， 还可对B进行降维，此时方程变为：
2 数值实验
本文算法用 Python实现， 运行于MacOS上。所用数据、源代码和运行结果都托管在Github 上， 网址为https：//github.com/Freakwill/PCA-linear-equations.
2.1 主成分选择
对矩阵A，B降维， 通过观察累计百分比曲线， 选择适当的主成分数。如图1所示。
2.2 误差分析
误差采用向量的2范数， 对矩阵来说就是Frobenius 范数。本文用相对误差公式衡量算法逼近能力。
20%的样本会被事先抽取出来， 测试误差通常比训练误差更有说服力。降维是对训练样本进行的， 而不是对测试样本进行。但对A的降维也可以应用于训练样本， 因为其是已知输入值， 而测试样本中的B作为目标值， 不参与降维。如图2所示。
由图2结果可知：
（1）原方程组没有解。图中原方程误差是理论上最小（在没有约束的情况下）， 不妨看做一个相对0点。如果降维方法的误差低于这个点， 可认为是系统误差造成的， 不是数值分析意义上的误差。
（2）降维方程组的误差衰减达到预期。“负数值置0”的处理， 对误差没有太大影响。当A的主成分数超过70时， 计算变得不稳定， 预测误差表现很夸张。
（3）从误差曲线不难看出， 解方程组时并不需要用到所有样本，样本数量对本算法没有太大影响。实际应用中， 随机挑选充分数量的样本即可。
图3是对B取前4个主成分得到的误差图， 和之前的实验相比，并没有显著影响误差。时间和主成分数基本上呈线性关系， 符合预期。为了获得较为准确的运行时间， 本次实验对每一个主成分数重复50次， 无论时间还是误差都取均值。
由图4结果可知：
（1）B主成分数对误差影响没有A大 （这是优点）， 由于B维数较低， 且主成分比较集中， 在第一个主成分处， 误差已降到0.3左右。
（2）用时与主成分数近似呈弱线性相关， B的降维也相应提高了效率。
总之， 对B降维是完全合理的。
可通过设定随机种子， 产生固定的训练-测试样本（实验重复50次）。
选用A前30个主成分，B前4个主成分。 获得实验结果見表2。
2.3 其它降维策略与拟合方法
NMF降维的效果和PCA相近， 但矩阵分解后需解较为复杂的方程。PCA的优势在于能分解出正交矩阵， 可以设计出更快捷的算法， 计算用时短。中心化PCA在主成分数较少时表现较好， 这是因为中心化PCA采用的是仿射变换， 比单纯的线性变换多了常数项， 但随着主成分数增加， 并没有表现出明显优势。目前中心化PCA只是简单重构A， 即： 其中，C1VT1是中心化后的分解， 即通常意义上的PCA， M是均值矩阵，其无法简化方程， 计算时间维持在某个常数。若能设计出简化方程的算法， 那么中心化PCA或许是不错的选择。
当主成分增加时， 负数解都是无法避免的。 数值实验表明主成分过多还有可能出现异常。几种求解方法的测试比较结果如图5所示。
显然，完全可以用其它方法求解降维后的方程组， 尤其是当人们只关心预测， 而不在乎X时。如果只为了建立A与B的联系， 完全可以采用一些非线性方法。表3中给出一些线性模型测试结果，用到了一些较复杂的计算策略， 含非线性成分， 误差无显著降低， 运行时间则较长。这些模型均由Python 第三方库scikit-learn实现[7-8]。
最后， 给出一般的计算框架， 可为解线性方程组开发新的算法：
3 结束语
PCA降维在解大型线性方程组表现的非常出色， 随着主成分增加， 误差快速递减。 这种简单的代数学技巧， 不仅计算快， 算法设计简单， 由于矩阵分解为对角矩阵和正交矩陣的乘积， 之后也无需解任何线性方程组。可以根据需要任意压缩原方程。
本文提供的方法可以应用于工业生产中。但在实验中发现过多的主成分可能使得算法不稳定。今后的工作重点：
（1）当变量有非负性或其它约束条件时， 本文还没有给出更合理的处理办法。
（2）如何实现有效的基于中心化PCA的算法。本文提到的中心化PCA没有真正起到压缩方程组的作用。
（3）应利用系数矩阵的一些特点， 如稀疏性， 设计更有针对性的算法。
参考文献
[1]HASTIE T， TIBSHIRANI R， FRIEDMAN J H. The elements of statistical learning：Data mining， inference， and prediction [M]. New York ：Springer-Verlag， 2001.
[2] HOTELLING H. Analysis of a complex of statistical variables into principal components[J]. Journal of Educational Psychology， 1993，24（6）：417-441.
[3] TURK M A， PENTLAND A P. Face recognition using eigenfaces[C]//Proceedings 1991 IEEE Computer Society Conference on ComputerVision and Pattern Recognition. Maui， HI， USA：IEEE，1991：586-591.
[4] SHALEV-SHWARTZ S， BEN-DAVID S. Understanding machine learning：from theory to algorithm[M]. New York， NY， USA ：Cambridge University Press， 2014.
[5] 焦李成，等. 稀疏学习、分类与识别[M]. 北京：科学出版社， 2017.
[6] 解锋昌， 韦博成， 林金官. 零过多数据的统计分析及其应用[M]. 北京：科学出版社， 2013.
[7] HACKELING G. Mastering machine learning with scikit-learn[M]. 2nd ed. Birmingham ：Packt Publishing， 2017.
[8] COURNAPEAU D. Scikit-learn[EB/OL]. [2019]. https：//scikit-learn.org/stable/.