Paper Reading《Optimal Whitening and Decorrelation》

Paper Reading《Optimal Whitening and Decorrelation》

@author：Heisenberg

paper url

code url

论文发表自2016《 The American Statistician》，主要涉及了PCA、PCA-cor、ZCA、ZCA-cor以及Cholesky分解五种白化方法的理论推导证明以及在鸢尾花数据集上的测试表现。最终结论是ZCA-cor白化得到与原始变量最大相似的球形变量；PCA-cor白化能得到与原始变量最大压缩的球形变量。

1. Introduction and Notation

对于数据矩阵$x=(x_1,\cdots ,x_d{)}^T$，数学期望为$E(x)=\mu =({\mu }_1,{\mu }_2,\cdots ,{\mu }_d)$，协方差矩阵$var(x)=\Sigma$。通过线性变换$z=(z_1,z_2,\cdots ,z_d)=Wx,s.t.var(z)=I$的过程称为白化或者球化。其中$W\in {\mathbb{R}}^{d\times d}$ , $var(z)=I=W^T{x}^{T}xW=W^T\Sigma W$，也即
$$W^{T}W={\Sigma }^{-1}$$
下面介绍几个推导过程中常用的矩阵分解：

• 根据$\rho =\frac{Cov(X,Y)}{\sqrt{Var(X)Var(Y)}}$，我们令$P$为皮尔逊自相关系数矩阵则有$\Sigma =V^{\frac{1}{2}}P{V}^{\frac{1}{2}}$。 其中$V$是形如$V=\text{diag}({\sigma }_1^2,\cdots ,{\sigma }_d^2)$的对角阵。
• 由特征值分解$\Sigma =U\Lambda U^T$，$\Lambda$是$\Sigma$的特征值按照大小顺序组成的对角阵，形如$\Lambda =\text{diag}({\lambda }_1,\cdots ,{\lambda }_r,\cdots ,0),r\le d$。$U$是对应的特征值单位化后的特征向量组成的正交阵，$s.t.U^{T}U=I,U^T=U^{-1}$；同理，${\Sigma }^{-\frac{1}{2}}=U{\Lambda }^{-\frac{1}{2}}{U}^T$
• 同样地，对自相关矩阵$P$做特征值分解有$P=G\Theta G^T$和$P^{-\frac{1}{2}}=G^T{\Theta }^{-\frac{1}{2}}{G}^T$，其中$G^{T}G=I$且$\Theta$是$P$的特征值按照大小顺序组成的对角阵。

2. Rotational freedom

2.1 非标准化数据

根据约束**(1)**其实我们无法唯一确定线性变换$W$。我们可以通过极分解将$W$分解为$W=Q_1{\Sigma }^{-\frac{1}{2}}$，其中$Q_1$是正交阵，满足于$Q_1^T{Q}_1=I_d$。

为什么可以这么分解呢？

$W^{T}W=({\Sigma }^{-\frac{1}{2}}{)}^T{Q}_1^T{Q}_1{\Sigma }^{-\frac{1}{2}}=({\Sigma }^{-\frac{1}{2}}{)}^T{I}_d{\Sigma }^{-\frac{1}{2}}={\Sigma }^{-1}$.

相当于对原数据进行了${\Sigma }^{-\frac{1}{2}}$的放缩，和一个$Q_1$的旋转变换。（正交阵的几何意义即是保持两点的欧式距离不变的线性变换，可以实现旋转以及镜像变换等）。

对${\Sigma }^{-\frac{1}{2}}$进行特征值分解可得：$W=Q_1(U{\Lambda }^{-\frac{1}{2}}{U}^T)=(Q_{1}U){\Lambda }^{-\frac{1}{2}}{U}^T$，可以视为对矩阵$W$进行的SVD分解。$Q_{1}U$即左奇异矩阵，可以理解原矩阵进行了$U^T$的旋转之后，进行了${\Lambda }^{-\frac{1}{2}}$的放缩后又进行了由$Q_1$决定的旋转。

2.2 标准化数据

很多情况下我们要处理标准化数据$i.e.V^{-\frac{1}{2}}x$.

基于此对W进行矩阵分解可得$W=Q_2{P}^{-\frac{1}{2}}{V}^{-\frac{1}{2}}$，其中$Q_2$是正交阵，满足于$Q_2^T{Q}_2=I_d$。

为什么可以这么分解呢？

$W^{T}W=(V^{-\frac{1}{2}}{)}^T(P^{-\frac{1}{2}}{)}^T{Q}_2^T{Q}_2{P}^{-\frac{1}{2}}{V}^{-\frac{1}{2}}=(V^{-\frac{1}{2}}{)}^T(P^{-\frac{1}{2}}{)}^T{I}_d{P}^{-\frac{1}{2}}{V}^{-\frac{1}{2}}=(V^{-\frac{1}{2}}{)}^T{P}^{-1}{V}^{-\frac{1}{2}}={\Sigma }^{-1}$

对$P^{-\frac{1}{2}}$进行特征值分解可得：$W=Q_{2}G{\Theta }^{-\frac{1}{2}}{G}^T{V}^{-\frac{1}{2}}$，可以视为对矩阵$W$进行的SVD分解。$Q_{2}G$即左奇异矩阵，可以理解原矩阵先进行方差矩阵$V^{-\frac{1}{2}}$的放缩，而后经历了$G^T$的旋转，再经过${\Theta }^{-\frac{1}{2}}$的放缩，又进行了由$Q_2$决定的旋转。

2.3 两种旋转自由度的关系

所以矩阵$W$可以有两种自由旋转的方式，分别由$Q_1,Q_2$决定。

由$Q_1{P}^{-\frac{1}{2}}{V}^{-\frac{1}{2}}=Q_2{\Sigma }^{-\frac{1}{2}}$，可得$Q_1=Q_{2}A,whereA=P^{\frac{1}{2}}{V}^{\frac{1}{2}}{\Sigma }^{-\frac{1}{2}}$

3 互协方差与互相关

Cross-covariance and cross-correlation

互协方差矩阵定义为：
$$\begin{array}{rl}\mathbf{\Phi }=\left({\varphi }_{ij}\right)& =\mathrm{cov}(\mathbf{z},\mathbf{x})=\mathrm{cov}(\mathbf{W}\mathbf{x},\mathbf{x})=\mathbf{W}\mathbf{\Sigma }=Q_1{\mathbf{\Sigma }}^{1/2}\end{array}$$
互自相关矩阵定义为：（互自相关矩阵对角线元素不一定为1）
$$\begin{array}{rl}\mathbf{\Psi }=\left({\psi }_{ij}\right)& =\mathrm{cor}(z,\mathbf{x})=\mathbf{\Phi }{\mathbf{V}}^{-1/2}=Q_{2}A{\Sigma }^{1/2}{\mathbf{V}}^{-1/2}=Q_2{\mathbf{P}}^{1/2}\end{array}$$

互相关函数是信号分析里的概念，表示的是两个时间序列之间的相关程度，即描述信号x(t),y(t)在任意两个不同时刻t1，t2的取值之间的相关程度。具体计算方法可以参考互相关（cross-correlation）及其在Python中的实现

可以看到互协方差和互自相关分别由旋转变量$Q_1,Q_2$决定，这为我们通过适当地选择或限制$\Phi$或$\Psi$来选择和区分白化变换提供了衡量标准。

对$\Psi$的讨论：

由$x=W^{-1}z$对于列维度的向量$x_j$都分别由不相关的$z_1,\cdots ,z_d$线性组合而成，所以$SMC(x_j,z)=1$.

SMC，即squared multiple correlations：一个变量与其他变量的多重相关性平方，作为对一个变量公共性的评估。可以理解为在OLS中的$\rho (y_i,\hat{\mathbb{Y}})$

用矩阵表示则为$diag({\Psi }^T\Psi )=diag(P^{1/2}{Q}_2^T{Q}_2{P}^{1/2})=diag(P)=(1,\cdots ,1{)}_d^T$。即$z,x$的组内关系是确定的，$x_j$这一列完全由$z_j$这一列线性变换而来，故$\rho (z_j,x_j)=1$.

而组间关系$\rho (x_i,z_i)$则要用$diag(\Psi {\Psi }^T)$来表示，不同的白化方法对应着不同的值。

4 不同的白化方法

由上文所述方法我们确定了范式后可知，一种白化方法我们确定了旋转矩阵$Q_1,Q_2$,也就可以推出白化矩阵$W$.反之亦然。

4.1 ZCA

ZCA 也即zero-phase components analysis

$W^{ZCA}={\Sigma }^{-1/2}$，$Q_1=I$.

4.2 PCA

${\mathbf{W}}^{\mathrm{P}\mathrm{C}\mathrm{A}}={\mathbf{\Lambda }}^{-1/2}{\mathbf{U}}^T$

4.3 Cholesky分解

cholesky分解即将一个方阵分解为一个上三角形和一个下三角形矩阵的乘积，两者互为转置。

${\mathbf{W}}^{\mathrm{C}\mathrm{h}\mathrm{o}\mathrm{l}}={\mathbf{L}}^T$

4.4 ZCA-cor分解

由Zummer（2009）等人首先用于变量重要性以及变量筛选统计上。

${\mathbf{W}}^{\mathrm{Z}\mathrm{C}\mathrm{A}-\mathrm{c}\mathrm{o}\mathrm{r}}={\mathbf{P}}^{-1/2}{\mathbf{V}}^{-1/2}$，相当于按照标准差进行了标准化，然后又做了自相关系数$P^{-1/2}$的放缩。

4.5 PCA-cor分解

${\mathbf{W}}^{\mathrm{P}\mathrm{C}\mathrm{A}-\mathrm{c}\mathrm{o}\mathrm{r}}={\Theta }^{-1/2}{G}^T{\mathbf{V}}^{-1/2}$

五种白化方法总结：
$$\begin{array}{lrrrll}& \begin{array}{r}\text{Sphering matrix}\\ W\end{array}& \begin{array}{r}\text{Cross-covariance}\\ \mathbf{\Phi }\end{array}& \begin{array}{r}\text{Cross-correlation}\\ \mathbf{\Psi }\end{array}& \begin{array}{l}\text{Rotation mat}\\ {\mathbf{Q}}_1\end{array}& \begin{array}{l}\text{Rotation mat}\\ {\mathbf{Q}}_2\end{array}\\ \text{ ZCA }& {\mathbf{\Sigma }}^{-1/2}& {\mathbf{\Sigma }}^{1/2}& {\mathbf{\Sigma }}^{1/2}{\mathbf{V}}^{-1/2}& \mathbf{I}& \mathbf{I}\\ \text{ PCA }& {\mathbf{\Lambda }}^{-1/2}{\mathbf{U}}^T& {\mathbf{\Lambda }}^{1/2}{\mathbf{U}}^T& {\mathbf{\Lambda }}^{1/2}{\mathbf{U}}^T{\mathbf{V}}^{-1/2}& {\mathbf{U}}^T& {\mathbf{U}}^T{\mathbf{A}}^T\\ \text{ Cholesky }& {\mathbf{L}}^T& {\mathbf{L}}^T\mathbf{\Sigma }& {\mathbf{L}}^T\mathbf{\Sigma }{V}^{-1/2}& {\mathbf{L}}^T{\mathbf{\Sigma }}^{1/2}& {\mathbf{L}}^T{\mathbf{V}}^{1/2}{\mathbf{P}}^{1/2}\\ \text{ ZCA-cor }& {\mathbf{P}}^{-1/2}{\mathbf{V}}^{-1/2}& {\mathbf{P}}^{1/2}{\mathbf{V}}^{1/2}& {\mathbf{P}}^{1/2}& \mathbf{A}& \mathbf{I}\\ \text{ PCA-cor }& {\Theta }^{-1/2}{G}^T{\mathbf{V}}^{-1/2}& {\Theta }^{1/2}{\mathbf{G}}^T{\mathbf{V}}^{1/2}& {\Theta }^{1/2}{\mathbf{G}}^T& {\mathbf{G}}^T\mathbf{A}& {\mathbf{G}}^T\end{array}$$

5 优化白化过程

优化目标是通过最小的额外调整来消除相关性，目的是使转换后的变量z保持尽可能类似于原始向量x。PCA侧重于前者，有利于数据压缩降维；ZCA侧重于后者，有利于映射回向量源空间。

5.1 ZCA

ZCA中我们采用最小化原始变量$x$和白化变量$z$之间总平方距离的白化变换。在下列变换中$\mathrm{E}(z_c)=\mathrm{E}(x_c)=0$。
$$\begin{array}{rl}\mathrm{E}\left({\left({\mathbf{z}}_c-{\mathbf{x}}_c\right)}^T\left({\mathbf{z}}_c-{\mathbf{x}}_c\right)\right)& =\mathrm{E}(z_c^T{z}_c^T)-2\mathrm{E}({\mathbf{z}}_c{\mathbf{x}}_c^T)+\mathrm{E}(x_c^T{x}_c)\\ & =\mathrm{tr}(\mathbf{I})-2\mathrm{E}\left(\mathrm{tr}\left({\mathbf{z}}_c{\mathbf{x}}_c^T\right)\right)+\mathrm{tr}(\mathbf{\Sigma })\\ & =d-2\mathrm{tr}(\mathbf{\Phi })+\mathrm{tr}(\mathbf{V})\end{array}$$
维度$d$与$tr(V)={\sum }_{i=1}^d{\sigma }_i^2$与矩阵$W$都无关，最小化平方距离等价于最大化$\mathrm{tr}(\mathbf{\Phi })$。
$$\mathrm{tr}(\mathbf{\Phi })=\sum _{i=1}^d\mathrm{cov}\left(z_i,x_i\right)=\mathrm{tr}\left(Q_1{\Sigma }^{1/2}\right)\equiv g_1\left(Q_1\right)$$
可以证明当且仅当$W=W^{ZCA}$的时候$g_1(Q_1)$取得最大值。

Proof：

先对$\Lambda$进行特征值分解，后半部分可视为对矩阵$B$的帖子分解，最后一个等号由$\Lambda$是对角阵所保证

$g_1\left({\mathbf{Q}}_1\right)=\mathrm{tr}\left({\mathbf{Q}}_1\mathbf{U}{\mathbf{\Lambda }}^{1/2}{\mathbf{U}}^T\right)=\mathrm{tr}\left({\mathbf{\Lambda }}^{1/2}{\mathbf{U}}^T{\mathbf{Q}}_1\mathbf{U}\right)\equiv \mathrm{tr}\left({\mathbf{\Lambda }}^{1/2}\mathbf{B}\right)={\sum }_i{\mathbf{\Lambda }}_{ii}^{1/2}{B}_{ii}$

而$Q_1,U$的正交性则保证了$B=U^T{Q}_{1}U$的正交性。

$\therefore B_{ii}\le 1$. （正交阵行列向量都为单位向量，则对角线数值必≤1）

$iff\mathbf{B}=\mathbf{I},i.e.{\mathbf{Q}}_1=\mathbf{I},i.e.\mathbf{W}={\mathbf{W}}^{\mathbf{Z}\mathbf{C}\mathbf{A}}$ , $g_1(Q_1)$可以取最大值

因而，我们发现ZCA白化是最大化了$z$和$x$每一个成分之间的平均互协方差。

5.2 ZCA-cor

同样的，对于标准化的向量$V^{-1/2}{x}_c$

${\mathbf{W}}^{\mathrm{P}\mathrm{C}\mathrm{A}}={\mathbf{\Lambda }}^{-1/2}{\mathbf{U}}^T$

4.3 Cholesky分解

cholesky分解即将一个方阵分解为一个上三角形和一个下三角形矩阵的乘积，两者互为转置。
$$\mathrm{E}\left({\left({\mathbf{z}}_c-{\mathbf{V}}^{-1/2}{\mathbf{x}}_c\right)}^T\left({\mathbf{z}}_c-{\mathbf{V}}^{-1/2}{\mathbf{x}}_c\right)\right)=2d-2\mathrm{tr}(\mathbf{\Psi })$$
可以证明当且仅当$W=W^{ZCA-cor}$的时候$g_2(Q_2)$取得最大值。

结果,我们确定ZCA-cor是确保$z$的分量保持最大限度地与相应的$x$的分量的相似性。

5.3 PCA

PCA目标是数据压缩，我们的目标是用重构的数据$z_1,\cdots ,z_d$中的前几个向量就能够尽可能的表示出原来的向量所$x_1,\cdots ,x_d$包含的信息。

PCA中我们的目标是找到一个白化向量$z$，使得$z$与$x$的组间向量（即行向量）互协方差${\sum }_{j=1}^d{\varphi }_{ij}^2={\sum }_{j=1}^d\mathrm{cov}{\left(z_i,x_j\right)}^2$

${\left({\varphi }_1,\dots ,{\varphi }_d\right)}^T=\mathrm{diag}\left(\Phi {\Phi }^T\right)=\mathrm{diag}\left(Q_1\Sigma Q_1^T\right)\equiv h_1\left(Q_1\right)$，满足于${\varphi }_i\ge {\varphi }_{i+1}$,使得${\varphi }_i$最大。

$h_1(Q_1)=\text{diag}(Q_1\Sigma Q_1^T)=\text{diag}(Q_{1}U\Lambda U^T{Q}_1^T)$

当且仅当$Q_1=U^{,}i.e.W=W^{PCA}$时，$h_1(Q_1)={\sum }_j{\Lambda }_{jj}{D}_{ij}^2,whereD=Q_{1}U=I$

可以取得最大值。

因此，PCA基于互协协方差作为基础度量，最大化了球形向量$z$的每个分量对应的原始向量$x$的所有分量的积分或压缩。因此，PCA在降维方面是最优排序的。

5.4 PCA-cor

同样地，PCA-cor中我们的目标是找到一个白化向量$z$，使得$z$与$x$的组间向量（即行向量）互自相关系数${\sum }_{j=1}^d{\psi }_{ij}^2={\sum }_{j=1}^d\mathrm{cov}{\left(z_i,x_j\right)}^2$

${\left({\psi }_1,\dots ,{\psi }_d\right)}^T=\mathrm{diag}\left(\Psi {\Psi }^T\right)=\mathrm{diag}\left(Q_{2}P{Q}_2^T\right)=h_2\left(Q_2\right)$

同理可证 $W=W^{PCA-cor}$时成立。

因此，PCA-cor白化变换是用互相关ψ为基础测度,将原始向量$x$的所有分量在白化向量$z$的每个分量中的积分或压缩最大化的唯一变换。

5.5 Cholesky

Cholesky能够由下三角正对角互协方差和互自相关矩阵这两个特征唯一确定。

使用Cholesky因子分解进行白化的一个结果是，我们隐含地假定了变量的排序。这对于解释自相关的时间过程分析是非常有用的。

6 实验与结论分析

6.1 实验结果

选用了鸢尾花$4\times 150$的数据进行实验，分别进行五种白化方法进行实验：

最好的结果加粗表示；第二好的结果斜体表示

• 可以看到$tr(\hat{\Psi }),tr(\hat{\Phi })$中表现最好和第二好的都是ZCA/ZCA-cor，数值越大表明$z,x$的平方距离越小，拟合能力越强；

• $\max \mathrm{diag}\left(\Psi {\Psi }^T\right),\max ,\mathrm{diag}\left(\Phi {\Phi }^T\right)$中表现最好和第二好的都是PCA/PCA-cor，数值越大说明最大的对角线元素打，在相同压缩力度下，第一个元素$z_1$所能代表$x$的能力强。

• 且根据上半部分可以看出，ZCA/ZCA-cor中$z_i,x_i$的相关系数整体都很高，而PCA/PCA-cor中$\widehat{cor}(x_1,z_1),\widehat{cor}(x_2,z_2)$自相关系数很高，但后面$\widehat{cor}(x_3,z_3),\widehat{cor}(x_4,z_4)$就很低了，说明前两维已经能很好地刻画$x$。

• 无论如何，Cholesky白化总是排名第三，要么排在ZCA/ZCA-cor白化之后，要么排在PCA/PCA-cor白化之后。此外，它是唯一能通过构造一对$(z_4,x_4)$在完美地相互关联。

6.2 结论

一个合适的白化变换取决于应用目的，特别是是否需要最小的调整，还是数据压缩。
在前者中，白化后的变量与原始变量保持高度相关，从而保持其原始的解释。这在如变量选择上市有利的。
相比之下，在压缩环境中，通过构造而变白的变量与原始数据没有可解释的关系，而是反映了它们内在的有效维度。

总的来说，作者提倡使用以互自相关系数矩阵$\mathbf{\Psi }$为优化基础的尺度不变的最优性函数。
如果目的是获得与原始变量最大相似的球形变量，我们建议采用ZCA-cor白化；相反，如果需要最大压缩，我们建议使用PCA-cor白化。

7 总结

作者通过理论推导对五种白化方法的进行了原理性的解释和证明，并结合经典数据集的实验给出了何种情况下适合何种白化方法。
但也存在着一定的问题：

1. 只用一个经典数据集是否缺乏共性（可能是数理统计学的普遍做法）
2. 并没有给出为何基于Self-cor mat $\Psi$ 是更合理做法的更详尽的解释。实验结果中ZCA、PCA白化也都在$tr(\hat{\Phi })$,$\max \text{diag}(\hat{\Phi }{\hat{\Phi }}^T)$取得了最好。可能是ZCA-cor,PCA-cor都是基于$P$的分解所以在$tr(\hat{\Psi })$,$\max \text{diag}(\hat{\Psi }{\hat{\Psi }}^T)$的尺度上更统一？
3. 在文章开头引出$V$的时候并没有指明$\Sigma =V^{1/2}P{V}^{1/2}$是哈达玛积，即element-wise按照位置相乘的。$V$其实是形如$\left|\begin{array}{lll}Var(X)& \sqrt{Var(X)\cdot Var(Y)}& \sqrt{Var(X)\cdot Var(Z)}\\ \sqrt{Var(Y)\cdot Var(X)}& Var(Y)& \sqrt{Var(Y)\cdot Var(Z)}\\ \sqrt{Var(Z)\cdot Var(X)}& \sqrt{Var(Z)\cdot Var(Y)}& Var(Z)\end{array}\right|$的方差-标准差矩阵。$diag(V)=({\sigma }_1,{\sigma }_2,\cdots ,{\sigma }_d)$。而第二部分Notation开始就说"the diagonal varience matirx V"，不严谨。