【文献阅读】

Title

Sensor placement minimizing the state estimation mean square error:Performance guarantees of greedy solutions

Autor information

引用格式：Erichson, N. Benjamin, et al. “Sparse principal component analysis via variable projection.” SIAM Journal on Applied Mathematics 80.2 (2020): 977-1002.pdf链接

Full text

Abtract

稀疏主成分分析 (SPCA) 已成为现代数据分析的强大技术，通过识别数据中的局部空间结构和消除不同时间尺度之间的歧义，改进了对低秩结构的解释。 我们通过将其表述为价值函数优化问题来展示稳健且可扩展的 SPCA 算法。
这种观点导致了一种灵活且计算效率高的算法。 此外，我们可以利用线性代数中的随机方法将该方法扩展到大规模（大数据）设置。 我们提出的创新还允许使用稳健的 SPCA 公式，尽管输入数据严重损坏，该公式仍能获得有意义的稀疏主成分。 所提出的算法使用合成数据和真实世界数据进行了演示，并显示出卓越的计算效率和诊断性能。

Introduction

PCA 分解将时变模式表示为系统状态随时间演化的主要相关空间活动的线性组合。 虽然常用，但 PCA 方法生成的全局模式通常混合或混合各种时空尺度，并且无法识别在不同尺度上起作用的潜在底层动态。 此外，经典 PCA 也倾向于过度拟合观察数小于变量数的数据 [16]。
约束或正则化矩阵分解为动态模式建模提供了更灵活的方法。具体来说，可以通过稀疏性促进正则化器引入先验信息，以获得更简约的数据近似值，这通常会提供改进的可解释性。 其中，稀疏主成分分析 (SPCA) 已成为一种流行且强大的现代数据分析技术。 SPCA 促进了模式的稀疏性，即稀疏模式只有少数active系数，而大多数系数被限制为零。 由此产生的稀疏模式通常是高度局部化的，并且比从传统 PCA 获得的全局 PCA 模式更具可解释性。 因此，PCA 的稀疏正则化允许分解策略可以专门识别数据中的局部空间结构并消除不同时间尺度之间的歧义。
虽然稀疏权重向量的想法并不新鲜，但简单的临时技术（如朴素阈值）可能会导致误导性结果。 Jolliffe 等人首先提出了一种使用 $l_1$正则化的 SPCA 正式方法[38]。 这项开创性工作导致了各种稀疏性提升算法 [70、23、24、58、57、67、39]。 稀疏 PCA 在获得可解释模式方面的成功激发了本文开发的通用方法。
具体来说，我们的方法比之前提出的 SPCA 算法提供了三个直接改进：(1) 更快且更具可扩展性的算法，(2) 对异常值的鲁棒性，以及 (3) 对非凸正则化器的直接扩展，包括 $l_0$范数。 可扩展性对于许多应用程序至关重要——例如，动力系统生成非常大规模的数据集，例如本文分析的海面温度数据。 稳健的公式允许 SPCA 部署在更广泛的环境中，否则数据污染可能会隐藏稀疏模式。 非凸正则化器目前在 SPCA 软件中不可用——我们证明了我们通过这些方法获得的模式在合成示例中更好，并且对真实数据更具解释性。
这项工作的贡献在这项工作中，我们为 SPCA 开发了一种可扩展且稳健的方法。 该方法的一个关键特征是使用变量投影来部分最小化正交约束变量。 这个想法被用在 [70] 的原始交替方法中，我们通过将问题重铸为价值函数优化来创新这个想法。 这种观点允许明显更快的算法、可扩展性和更广泛的适用性。 我们还允许对负载进行非凸正则化，这进一步提高了可解释性和稀疏性。 该方法不仅可以很好地扩展，而且可以使用线性代数的随机方法进一步加速 [28]。 此外，所提出的方法扩展到稳健的 SPCA 公式，即使输入数据严重损坏，它也可以获得有意义的主成分。 异常值被建模为对数据的扰动，就像在稳健的 PCA 模型中一样 [17,11,3]。 这些创新为现代数据分析和诊断提供了一种灵活高效的算法，能够以以前其他领先算法无法实现的规模实现广泛的关键应用。

PCA

主成分$z_i\in {\mathcal{R}}^n$被写为变量的线性组合，其中$X\in {\mathcal{R}}^{n\times p}$，$n$行表示观测量，$p$列表示变量的测量。
$$z_i=X{a}_i$$
$a_i\in {\mathcal{R}}^p$是权重向量
$$Z=XA$$
$A\in {\mathcal{R}}^{p\times p}$的列通常表示为模式，基函数，主方向或者载荷。
在数学上，可以制定方差最大化问题来找到权重向量 ai。 或者，该问题可以表述为最小二乘问题，即最小化输入数据和投影数据之间的残差平方和：
$$\begin{gathered} \underset{A}{\min }f(A)=||X-XA{A}^T|{|}_F^2 \\ subjectto{A}^{T}A=I(3) \end{gathered}$$
其中PCA对权重矩阵A施加正交限制，给定$X$的SVD分解$X=U\sum V^T$
（3）的最小化由$X$的右奇异值向量$V$给出，$Z=U\sum (U,V都是酉变换，即U^{T}U=V^{T}V=I)$。可以对$\sum$截断以产生降维效果.

PCA和最小二乘问题

主成分分析（Principal Component Analysis，PCA）可以表述为一个最小二乘问题来求解。在PCA中，我们希望找到一个线性投影，将原始高维数据映射到低维空间，使得投影后的数据具有最大的方差。

假设我们有一个原始数据矩阵 $X$，其中每一列表示一个数据样本，共有 $n$个样本。我们的目标是找到一个投影矩阵 W，将原始数据 X 映射到低维空间 $Y$，使得 $Y$ 的每个维度之间的协方差为零，同时最大化投影后的数据的方差。

具体来说，我们的目标是最小化重构误差，该误差是原始数据 X 与其低维重构 $\hat{X}$之间的平方差。定义重构误差$E$ 为：
$$E=|X-\hat{X}{|}^2$$
我们可以将$\hat{X}$表示为 $Y$ 与投影矩阵 $W$ 的乘积：
$$\hat{X}=WY$$
由于投影后的数据 $Y$ 的每个维度之间的协方差为零，我们有：
$$Cov(Y)=Y{Y}^T=I$$
其中，$I$ 是单位矩阵。这意味着投影矩阵 $W$ 是正交矩阵，即 $W^{T}W=I$。

通过最小化重构误差 $E$，我们可以将 PCA 问题表述为一个最小二乘问题。我们要求解的是投影矩阵 $W$，使得 $E$ 最小。根据最小二乘的原理，我们可以通过求解以下优化问题来获得 $W$ 的解：
$$\begin{gathered} \underset{W}{\min }|X-WY{|}^2 \\ subjectto:W^{T}W=I \end{gathered}$$
这是一个约束最小二乘问题，可以使用数值优化方法（如特征值分解、奇异值分解或迭代算法）来求解。通过求解该最小二乘问题，我们可以获得最优的投影矩阵 W，从而实现主成分分析。

Variable projection

$$\underset{A,B}{\min }g(A,B)(4)$$
非线性最小二乘问题
$$\underset{A,B}{\min }\frac{1}{2}||X-\Phi (B)A|{|}^2(5)$$
变量投影指的是$X$在$\Phi (B)$范围内的最小二乘投影具有封闭形式的解，该解在迭代方法中用于明确优化$B$，与部分最小化的意思相近。将（4）改为值优化问题：
min ⁡ B { v ( B ) : = min ⁡ A g ( A , B ) }   ( 6 ) \min_B\{v(B):=\min_A g(A,B)\}\ (6) Bmin​{v(B):=Amin​g(A,B)} (6)
$v(B)$一般有一个确定的形式：
$$v(B)=\frac{1}{2}||X(I-{\mathcal{P}}_{\mathcal{R}(\Phi (B))})|{|}^2=dis{t}^2(X|\mathcal{R}(\Phi (B)))$$
令
$$A(B)=\arg \underset{A}{\min }g(A,B).$$
对于很多问题，$v(B)$的一阶和二阶导数可以求得。

稀疏PCA

直接将稀疏性诱导正则化器纳入优化问题（将SPCA视为正则化回归问题）：
$$\begin{gathered} \underset{A,B}{\min }f(A,B)=\frac{1}{2}||X-XB{A}^T|{|}_F^2+\psi (B) \\ subjectto{A}^{T}A=I(7) \end{gathered}$$
$B$是稀疏权重矩阵, $A$是正交矩阵，惩罚表示稀疏性诱导正则化器。使用交替部分算法最小化优化问题。
1、Update $A$. 当$B$固定，这是一个正交Procrustes问题，有封闭解$A^{\ast }=U{V}^T$,其中$X^{T}XB=U\sum V^T$.
2、Update $B$. $A$固定时，求解B，问题跨越$B$的$k$列，产生每一种情况下$b_j,j\in [1,k]$的正则化回归问题。

然后将主成分形成为观察到的变量 $Z=XB$ 的稀疏加权线性组合。 数据可以近似旋转为 $\hat{X}=Z{A}^T$。

坐标下降或最小角回归 (LARS) 用于解决 k 个子问题中的每一个 [26]。 B更新依赖于解决一个强凸问题，特别是更新是唯一的，并且通过[64]的分析所描述的算法收敛到一个固定点。 替换为非凸正则化器，例如 $\psi (B)=\alpha ||B|{|}_0+\beta ||B|{|}^2$，很难保证 B 更新的任何内容。 然而，正如我们所展示的，从可变投影的角度使用价值函数等式 6 可以产生有效的实现和直接的收敛分析。

Variable Projection角度解决SPCA

A的更新有封闭形式的解，而B的更新需要迭代方法。为了利用A的更新效率，我们投影出A并引入SPCA价值函数
$$v(B):=\underset{A}{\min }\frac{1}{2}||X-XB{A}^T|{|}_F^{2}subjectto{A}^{T}A=I$$
原始的SPCA问题的等式（7）可以看作：
$$\underset{B}{\min }v(B)+\psi (B)(8)$$
使用近端（prox)算法来最小化（8），这种方法使用单个运算符来更新B而不是使用迭代例程。而且Variable Projection的等式（8）还可以使用Huber损失函数的SPCA，而简单替换等式（7）中的$\psi (B)$会破坏a更新的有效结构。而我们可以使用Huber函数的特殊特征获得三个而不是两个变量的公式。
SPCA 框架包含一系列稀疏诱导正则化器。 稀疏性是通过在模型中引入额外信息以找到 B 中最有意义的“活动”（非零）条目来实现的，而大多数负载被限制为零。 当许多变量是冗余的，即不需要捕获底层连贯模型结构时，稀疏方法很有效。 正则化还可以防止过度拟合，并提供解决不适定问题的途径，这些问题在高维数据集分析中经常遇到。附录 A 简要讨论了一些促进稀疏性的流行正则化器。

Huber loss function

Huber损失函数是回归问题中常用的损失函数，尤其是在处理异常值时。它是小误差的平方损失（均方误差）和大误差的绝对损失（平均绝对误差）的组合。 Huber 损失提供了这两个损失函数之间的折衷，允许对异常值具有鲁棒性，同时仍然惩罚大的错误。 Huber损失定义如下：
$$\begin{gathered} Huberloss=\sum [\frac{1}{2}\times (y_i-{\hat{y}}_i{)}^2]if|y_i-{\hat{y}}_i|\le \delta , \\ \delta \times (|y_i-{\hat{y}}_i|-\frac{1}{2}\times \delta )if|y_i-{\hat{y}}_i|>\delta , \end{gathered}$$
δ 是一个超参数，它确定平方损失和绝对损失之间的阈值。 当绝对差$|y_i-{\hat{y}}_i|$小于或等于阈值 δ，Huber 损失的行为类似于平方损失，并以二次方式惩罚误差。当绝对差大于阈值$\delta$时，Huber损失切换到绝对损失并对误差进行线性惩罚。 通过调整阈值$\delta$的值，我们可以控制损失函数对异常值的敏感度。较小的 $\delta$值使损失对离群值更稳健，因为它允许将更大范围的错误视为小错误。相反，较大的 $\delta$值会使损失对异常值更敏感，因为它需要将较小范围的误差视为小误差。 Huber 损失函数通常用于鲁棒回归算法，例如 Huber 回归或带有 M 估计量的鲁棒回归。它在平方损失和绝对损失之间取得平衡，为异常值提供鲁棒性，同时仍然是可微的并且适合优化。

随机化

$\hat{X}\in {\mathcal{R}}^{l\times p}$是$X\in {\mathcal{R}}^{p\times p}$的sketch。$l$被选为略大于目标秩$k$。形成样本矩阵$Y\in {\mathcal{R}}^{n\times l}$
$$Y=X\Omega$$
其中$\Omega \in {\mathcal{R}}^{p\times l}$是随机生成的测试矩阵。然后通过计算$Y=QR$的QR分解来获得标准正交基矩阵。然后将输入矩阵投影到Y的范围内形成，是低维的。
$$\hat{X}=Q^{T}X$$
只需要执行一次投影步骤以初始化（随机化的）SPCA预处理步骤。

Robust SPCA

传统上，SPCA 被表述为最小二乘问题，然而，众所周知，平方损失对异常值很敏感。 在许多现实世界的情况下，我们面临着数据因测量错误或其他影响而严重损坏的挑战。 这激发了对可以更有效地解释损坏或丢失数据的稳健方法的需求。 事实上，几位作者提出了一种稳健的 SPCA 公式，使用 $l_1$ 范数作为稳健的损失函数来处理严重损坏的数据 [47、21、36]。
对于 SPCA 的稳健公式，我们使用一个密切相关的想法，将数据矩阵分成低秩模型和稀疏模型。 这种形式的加法分解被称为鲁棒主成分分析 (RPCA)，其将高维矩阵分离为低秩和稀疏成分的卓越能力使 RPCA 成为数据分析的宝贵工具 科学 [17, 11, 3]。 具体来说，我们建议使用 Huber 损失函数 $\rho ={\rho }_H$而不是 $l_1$范数作为数据失配。 Huber 范数克服了 Frobenius 范数的一些缺点，可以用作更稳健的拟合度量 [35, 46]。

Huber loss function

$${\rho }_H(x;k)=\{\begin{array}{ll}k|x|-k^2/2,& |x|>k\\ x^2/2& |x|\le k\end{array}$$
$${\rho }_H(A;,k)=\sum _{i,j}{\rho }_H(A_{ij};k)$$
$${\rho }_H(x;k)=\underset{s}{\min }frac12||s-x|{|}^2+k||s|{|}_1.$$
如果不做适当改变，而是将上述鲁棒算法直接应用的话，会使等式（7）失去部分最小化A的封闭形式。等式（9）变为
$$\underset{A,B,S}{\min }{f}_H(A,B,S):=\frac{1}{2}||X-XB{A}^T-S|{|}_F^2+\psi (B)+k||s|{|}_1$$
A(B,S)由$U{V}^T$给出
$$(X-S{)}^{T}XB=U\sum V^T$$

时空SPCA

诸如 SVD 或本征正交分解 (POD) 之类的标准正交分解可能会出现过度拟合，并且生成的空间模式在空间上是密集的。 通过提升模式的稀疏性，SPCA 能够产生更易于解释的模式。
时空模态分析的目标是在空间和时间上可分离的系统分解，
$$x(t)=\sum _{j=1}^r{a}_j(t){\varphi }_j$$
${\varphi }_j$是在空间位置评估的模式，是fixed in time。有几种方法可以随时间的推移调整基础以提高性能，这些自适应方法包括增量SVD、动态正交模式、自适应h细化、最佳时间相关（OTM)模式。经典数据驱动分析寻求低秩近似，给定时间快照的数据矩阵
$$X=[x(t_1),x(t_2),...,x(t_p)]$$
适当的正交分解(POD)是高维流分析中的规范数据驱动分解，它寻求数据的最优秩 r 正交投影，该投影近似于 X 的协方差。最优低秩投影由显性 k 缩放给出 主成分 $Z=XV=U\sum$，由奇异值分解 $x=U\sum V^T$ 得到。 我们将模式定义为$\Phi =U$，导致可分离分解 $\underline{X=\Phi C^T}$, 其中 $C^T=\sum V^T$ 。 虽然 POD 模式在数值上近似于数据，但它们可能没有物理意义。
POD 模式通常不对应于持续存在的连贯结构。 另一方面，稀疏 PCA 在保持时间独立性的同时在空间模式中施加稀疏性。 在我们的框架中，空间模式由:
$$\Phi =XB=[X{b}_1,X{b}_2,...,X{b}_k]$$
表示快照（snapshot）和稀疏模式$\Phi$的稀疏线性组合，B的列是稀疏权重向量。

Conclusion

本文提出了一种用于计算稀疏主成分分析 (SPCA) 的稳健且可扩展的体系结构。
具体来说，我们将 SPCA 建模为具有正交约束的矩阵分解问题，并开发了专门的优化算法，可部分最小化变量的子集（变量投影）。 我们的 SPCA 算法具有可扩展性和鲁棒性，与当前最先进的方法相比大大提高了计算效率，同时保持了可比的性能。 更准确地说，我们已经证明： (i) 价值函数视角方法为 SPCA 提供了一个高效灵活的框架； (ii) 可以使用 Huber 损失来制定稳健的 SPCA； (iii) 可以将各种稀疏诱导正则化器纳入框架； (iv) 所提出的算法对于高维数据（即大 p）具有计算效率； (v) 线性代数的随机方法大大简化了计算需求，同时获得了低秩数据的近乎最优的近似值。
SPCA 是一种有用的数据诊断工具，具有丰富的动态特性，可在空间和时间上产生多尺度结构。 鉴于此类现象在物理、工程、生物和社会科学中普遍存在，这项工作为提高可解释性提供了一个有价值的工具，特别是在局部结构的诊断和不同时间尺度物理过程的歧义消除方面。 这项工作还为未来的发展开辟了许多途径：
方法论扩展 这种用于识别高维和多尺度数据中空间局部化空间结构的可扩展方法可以直接应用于 (1) 张量分解 [20、13、27]，它表示数据在 多维阵列结构，(2) 简约动力系统模型 [15]，识别捕获底层物理机制所需的最少非线性相互作用，以及 (3) 原位传感和控制，其中传感器和执行器通常需要 空间定位 [45]。
在工程和物理科学中的应用。 这里开发的方法将广泛适用于高维、多尺度的动力系统，并且需要可解释和简约的模型来进行预测、估计和控制。 SPCA 已经应用的具体应用包括大气化学 [66]、基因组学 [1、41] 和更广泛的生物系统 [43]。 此外，在不同领域都有巨大的进步机会，例如改进气候预测、检测和控制大脑结构以及湍流系统的闭环控制 [14]。

启发

1. Huber最小化
2. 代替POD等模式分解技术，获得解释性更强的稀疏低秩近似
3. 变量投影法。