【详细版】独立成分分析：算法和应用（一）

独立成分分析：算法和应用 翻译原文+自己的理解 原文：Independent component analysis: algorithms and applications 作者：A. Hyvarinen, E. Oja 时间：2000年 注：详细版原文太长，分两篇来写；另外可以也看简洁版

  假如在一个房间里有两个人同时讲话，我们使用放置在不同的位置两个麦克风来记录这两个人的讲话（信号）可以用 x 1 ( t ) x_1(t) x1​(t)和 x 2 ( t ) x_2(t) x2​(t)表示，其中 x 1 x_1 x1​和 x 2 x_2 x2​是幅度， t t t是时间。这些记录的信号中的每一个都是两个人发出的语音信号的加权和，我们用 s 1 ( t ) s_1(t) s1​(t)和 s 2 ( t ) s_2(t) s2​(t)表示。我们可以将其表达为线性方程： { x 1 ( t ) = a 11 s 1 + a 12 s 2 x 2 ( t ) = a 21 s 1 + a 22 s 2 (1) \left\{ \begin{aligned} x_1(t) & = & a_{11}s_{1}+a_{12}s_2 \\ x_2(t) & = & a_{21}s_{1}+a_{22}s_2 \\ \end{aligned} \right.\tag{1} {x1​(t)x2​(t)​==​a11​s1​+a12​s2​a21​s1​+a22​s2​​(1)   式中， a 11 , a 12 , a 21 a_{11},a_{12},a_{21} a11​,a12​,a21​和 a 22 a_{22} a22​是与麦克风和说话者之间距离相关参数。我们可以仅使用 x 1 ( t ) x_1(t) x1​(t)和 x 2 ( t ) x_2(t) x2​(t)来有效地估计 s 1 s_{1} s1​和 s 2 s_{2} s2​，这就是鸡尾酒问题（the cocktail party problem）。目前，我们忽略了任何时间延迟或其他额外因素以简化混合模型。

上述插图来源于博客（侵删）：https://blog.csdn.net/edogawachia/article/details/78886626

  作为例子，我们考虑图1和图2中的波形。当然，这些不是逼真的语音信号，但足以满足此说明的需要。原始语音信号看起来像图1的信号，混合信号看起来像图2中的信号。现在的问题就是，仅利用图2的数据恢复图1的数据。   实际上，如果我们知道参数 a i j a_{ij} aij​的值，就可以使用传统的方法求解公式1的线性方程。但是，重点是如果我们不知道 a i j a_{ij} aij​，问题就会变得很困难。   解决此类问题的一种方法是使用一些有关信号 s i ( t ) s_i(t) si​(t)的统计特性来估计 a i i a_{ii} aii​。实际上，也许令人惊讶的是，事实证明，假设在每个 t t t时刻 s 1 ( t ) s_1(t) s1​(t)和 s 2 ( t ) s_2(t) s2​(t)在统计上都是独立的。在许多情况下，这并非不切实际的假设，在实践中也不必完全正确。可以使用最近提出的基于其独立性信息来估计 a i j a_{ij} aij​的独立成分分析算法，这使我们能够从两个原始源信号 s 1 ( t ) s_1(t) s1​(t)和 s 2 ( t ) s_2(t) s2​(t)的混合信号 x 1 ( t ) x_1(t) x1​(t)和 x 2 ( t ) x_2(t) x2​(t)。图3给出了通过ICA方法估计的两个信号。可以看出，它们与原始源信号非常接近（它们的符号相反，但这没有意义）。   独立成分分析最初是为处理与鸡尾酒会问题密切相关的问题而提出的。自从最近人们对ICA的兴趣增加以来，该原理在许多其他方面的应用也越来越明显。例如，考虑由脑电图(EEG)给出的大脑活动的电记录。脑电图数据包括头皮上许多不同位置的电势记录。这些电势大概是通过混合一些大脑活动的潜在成分而产生的。这种情况与鸡尾酒会问题非常相似：我们想找到大脑活动的原始成分，但是我们只能观察到这些成分的混合。 ICA可以访问其独立成分，从而揭示有关大脑活动的有趣信息。   ICA的另一个非常不同的应用是特征提取。数字信号处理中的一个基本问题是为图像、音频或其他类型的数据（例如压缩和降噪）找到合适的表示形式。数据表示通常基于（离散）线性变换。在图像处理中广泛使用的标准线性变换是傅立叶变换、Haar变换（哈尔变换），余弦变换等。它们每个都有自己的优势（Gonzalez＆Wintz，1987）。   从数据本身估计线性变换将是最有用的，在这种情况下，该变换可以理想地适合于正在处理的数据类型。图4显示了ICA从自然图像块（窗口）中获得的ICA基函数（基向量）。训练图像集中的每个图像窗口都是这些窗口的叠加，因此叠加中的系数是独立的。稍后将详细解释ICA进行的特征提取。   实际上，上述所有应用都可以用统一的数学框架ICA来表述。这是一种非常通用的信号处理和数据分析方法。   在这篇综述中，我们将在第2节和第3节中介绍ICA的定义和基本原理。然后，从第4节，ICA问题是通过最小化或最大化某些对比函数（contrast functions）来解决的。这将ICA问题转换为数值优化问题（numerical optimization problem）。给出了许多对比函数，并阐明了它们之间的关系。第5节介绍了有用的预处理方法，可极大地帮助解决ICA问题，第6节介绍了解决问题的最有效的实践学习规则之一，即FastICA算法。然后，在第7节中，将介绍ICA的典型应用：从大脑信号记录中去除伪像，在金融时间序列中找到隐藏的因素，并减少自然图像中的噪声。第八节总结全文。

  为了严格定义ICA（Comon，1994年；Jutten＆Herault，1991年），我们可以使用统计“潜在变量”模型。假设我们观察到 n n n个独立分量的 n n n个线性混合 x 1 , ⋯   , x n x_1,\cdots,x_n x1​,⋯,xn​： x j = a j 1 s 1 + a j 2 s 2 + ⋯ + a j n s n    f o r   a l l   j (2) \begin{aligned} x_{j} & = & a_{j1}s_{1}+a_{j2}s_{2}+\cdots +a_{jn}s_{n}\ \ for \ all \ j \end{aligned}\tag{2} xj​​=​aj1​s1​+aj2​s2​+⋯+ajn​sn​  for all j​(2)   现在我们删除时间 t t t，在ICA模型中，我们假设每一个变量 x j x_j xj​和每一个独立分量 s k s_k sk​都是随机变量，而不是确定性信号（proper time signal，deterministic signal）。比如，在鸡尾酒会问题上观察到的麦克风信号 x j ( t ) x_j(t) xj​(t)就是一个随机变量样本。不失一般性（Without loss of generality），我们假设所有混合变量和独立分量都是零均值。如果不是的话，我们要对观测变量 x i x_i xi​进行中心化（减去均值），使其变成零均值模型。   为了方便，我们使用矢量–矩阵表示法（vector–matrix notation）代替前面公式中的求和形式。我们将使用 x \mathbf x x和 s \mathbf s s表示随机向量，其元素分别是： x 1 , ⋯   , x n x_1,\cdots,x_n x1​,⋯,xn​和 s 1 , ⋯   , s n s_1,\cdots,s_n s1​,⋯,sn​；用 A \mathbf A A表示矩阵，其元素是 a i j a_{ij} aij​。通常，粗体小写字母（bold lower case letters）表示向量，粗体大写字母（bold upper-case letters）表示矩阵。所有的向量都是列向量（column vector）， x T x^T xT（the transpose of x \mathbf x x）就是行向量（row vector）。使用矢量–矩阵表示法，上面的混合模型可以表示为（is written as）： x = A s (3) \begin{aligned} \mathbf {x}=\mathbf{As} \end{aligned}\tag{3} x=As​(3)   有时，我们需要矩阵 A \mathbf A A的列向量，用 a j \mathbf a_j aj​表示列向量，那么模型可以表示为： x = A s = [ a 1 , a 2 , ⋯   , a n ] s = [ a 1 , a 2 , ⋯   , a n ] [ s 1 s 2 ⋮ s n ] = ∑ i = 1 n a i s i (4) \mathbf x = \mathbf {As} =\mathbf {[ a _ { 1 } , a _ { 2 } , \cdots , a _ { n } ] s} =\mathbf {[ a _ { 1 } , a _ { 2 } , \cdots , a _ { n } ]} \left[ \begin{matrix} s_1 \\ s_2 \\ \vdots \\ s_n \end{matrix} \right] = \sum _ { i = 1 } ^ { n } \mathbf a _ { i } s _ { i }\tag{4} x=As=[a1​,a2​,⋯,an​]s=[a1​,a2​,⋯,an​]⎣⎢⎢⎢⎡​s1​s2​⋮sn​​⎦⎥⎥⎥⎤​=i=1∑n​ai​si​(4)   等式（3）中的统计模型称为独立成分分析或ICA模型。ICA模型是一个生成模型(generative model)，这意味着它描述了如何通过混合分量 的过程来生成观测数据。独立分量是潜在变量，这意味着它们不能直接观察。混合矩阵(the mixing matrix)也被认为是未知的。我们观察到的只是随机向量 ，我们必须使用它来估计 和 。这必须在尽可能一般的假设下进行。   ICA的出发点(The starting point for ICA)是非常简单的假设，即分量 s i s_i si​是统计独立的（statistically independent）。统计独立性（statistical independence）将在第3节中严格定义。从下面可以看出，我们还必须假设独立分量必须具有非高斯分布(non-Gaussian distributions)。但是，在基本模型中，我们假定这些分布是未知的（如果已知则可以大大简化问题）。为简单起见(For simplicity)，我们还假设未知混合矩阵为方阵，但有时可以放宽此假设，如第4.5节中所述。然后，在估计矩阵 A \mathbf A A之后，我们可以计算其逆，即 W \mathbf W W，并通过以下方法简单地获得独立分量： s = W x (5) \mathbf{s}=\mathbf{Wx}\tag{5} s=Wx(5)   ICA与盲源分离（blind source separation, BSS）或盲信号分离(blind signal separation)的方法密切相关。“源”在这里是指原始信号，即独立成分，例如鸡尾酒会问题中的说话者。“盲”是指我们对混合矩阵了解甚少，甚至对源信号也没有任何假设。ICA是执行盲源分离的一种方法，也许是使用最广泛的方法。 在许多应用中，假设测量中存在一些噪声（Hyva-rinen，1998a；Hyva-rinen，1999c）会更真实，这意味着在模型中添加噪声项。为简单起见，我们忽略了任何噪声项，因为无噪声模型的估计本身很难，而且对于许多应用来说似乎足够了。

  在公式（3）的ICA模型中，我们可以很容易的看出ICA具有以下不足： 1、无法确定独立分量的方差（能量）。 原因是，在 s \mathbf s s和 A \mathbf A A都不知道的情况下，总是可以通过将 A \mathbf A A的相应列 a i \mathbf a_i ai​除以相同的标量来抵消其中一个源中的任何标量乘积： x = ∑ i ( 1 α i a i ) ( s i α i ) (6) \mathbf{x}=\sum\limits_{i}{(\frac{1}{{{\alpha }_{i}}}{\mathbf{a}_{i}})({{s}_{i}}{{\alpha }_{i}})}\tag{6} x=i∑​(αi​1​ai​)(si​αi​)(6)   作为推论，我们也可以用同样的方式对独立分量的幅值进行修正。因为它们是随机变量，所以（进行修正）最自然的方法是假定每个变量都有单位方差： E { s i 2 } = 1 E\{{s_i^2}\}=1 E{si2​}=1。然后在ICA求解方法中，可以通过调整矩阵 A \mathbf A A来实现该归一化约束。但是，这仍然留下符号的不确定性问题：我们可以将某一独立分量乘以 − 1 -1 −1，而不会影响模型。幸运的是，这种歧义在大多数应用中都是微不足道的。 2、无法确定独立分量的顺序。   原因是，再次不知道 s \mathbf s s和 A \mathbf A A，我们可以自由地更改公式（4）中求和项的顺序，并将任何独立分量称为第一个分量。形式上，置换矩阵 P \mathbf P P及其逆可以在模型中替换为 x = A P − 1 P s \mathbf x=\mathbf{AP^{-1}Ps} x=AP−1Ps。 P s \mathbf {Ps} Ps的元素是原来的独立变量 s j s_j sj​，但顺序发生了变化。矩阵 A P − 1 \mathbf {AP^{-1}} AP−1只是一个新的未知混合矩阵，将由ICA算法求解。

为了以统计方式说明ICA模型，请考虑具有以下均匀分布的两个独立分量： p ( s i ) = { 1 2 3     i f   ∣ s i ∣ ≤ 3 0           o t h e r w i s e (7) p({{s}_{i}})=\left\{ \begin{aligned} & \frac{1}{2\sqrt{3}}\ \ \ if\ \left| {{s}_{i}} \right|\le \sqrt{3} \\ & 0\ \ \ \ \ \ \ \ \ otherwise \end{aligned} \right.\tag{7} p(si​)=⎩⎪⎨⎪⎧​​23 ​1​   if ∣si​∣≤3 ​0         otherwise​(7)   如上一节所述，选择这种均匀分布的取值范围，以使均值为零而方差为1。 和 分布的联合密度在一个正方形区域上是均匀的。这是从独立性的基本定义得出的，即两个独立变量的联合密度仅是其边际密度的乘积（参见等式10)：我们需要简单地计算两个概率密度函数的乘积即可。这样的联合概率密度的直观图如图5所示，该图是通过随机绘制该分布中的数据点而得到的。 现在将这两个独立的成分混合在一起，让我们采用以下混合矩阵： A 0 = [ 5 10 10 2 ] (8) \mathbf{{A}_{0}}=\left[ \begin{matrix} 5 & 10 \\ 10 & 2 \\ \end{matrix} \right]\tag{8} A0​=[510​102​](8)   这给了我们两个混合变量 x 1 x_{1} x1​和 x 2 x_2 x2​。容易计算出混合数据在平行四边形上具有均匀分布，如图6所示。 注意到这时随机变量 x 1 x_{1} x1​和 x 2 x_2 x2​不再相互独立：一个简单的验证方法是考虑是否可能从其中一个变量取值预测另外一个的取值。显然，如果 x 1 x_{1} x1​取其最大值或最小值，则 x 2 x_{2} x2​的值将完全确定。因此它们不是独立的。对于变量 s 1 s_1 s1​和 s 2 s_2 s2​情况不同：从图5可以看出，知道 s 1 s_1 s1​的值对猜测 s 2 s_2 s2​的值没有任何帮助。   现在，估计ICA的数据模型的问题是仅使用混合变量 x 1 x_{1} x1​和 x 2 x_2 x2​中包含的信息来估计混合矩阵 A 0 \mathbf {A_0} A0​。实际上，我们可以从图6中看到一种估算 A \mathbf {A} A的直观方法：平行四边形的边沿 A \mathbf {A} A列的方向。从原则上讲，这意味着我们可以通过首先估算 x 1 x_{1} x1​和 x 2 x_2 x2​的联合概率密度，然后找到边缘概率密度来估算ICA模型。因此问题似乎已解决。   然而，实际上这将是一种非常糟糕的方法，因为它仅适用于具有完全一致分布的变量。而且，它在计算上非常复杂。我们需要的是一种适用于独立分量的任何分布的方法，并且可以快速可靠地工作。   接下来，在开始提出用于评估ICA模型的方法之前，我们将考虑独立性的确切定义。

  要定义独立性的概念，请考虑两个标量值随机变量 y 1 y_1 y1​和 y 2 y_2 y2​。基本上，如果关于 y 1 y_1 y1​的值的信息没有给出关于 y 2 y_2 y2​的值的任何信息，则变量 y 1 y_1 y1​和 y 2 y_2 y2​可以说是独立的，反之亦然。 上面我们注意到，变量 s 1 、 s 2 s_1、s_2 s1​、s2​就是这种情况，而混合变量 x 1 、 x 2 x_1、x_2 x1​、x2​则不是。.   从技术上讲，独立性可以通过概率密度来定义。让我们用 p ( y 1 , y 2 ) p(y1,y2) p(y1,y2)表示 y 1 y_1 y1​和 y 2 y_2 y2​的联合概率密度函数（pdf）。让我们进一步用 p 1 ( y 1 ) p_1(y_1) p1​(y1​)表示 y 1 y_1 y1​的边际pdf，即单独考虑时 y 1 y_1 y1​的pdf： p 1 ( y 1 ) = ∫ p ( y 1 , y 2 ) d y 2 (9) p _ { 1 } ( y _ { 1 } ) = \int p ( y _ { 1 } , y _ { 2 } ) d y _ { 2 }\tag{9} p1​(y1​)=∫p(y1​,y2​)dy2​(9)   对于 y 2 y_2 y2​同样。然后，当且仅当联合pdf可通过以下方式分解时，我们定义 y 1 y_1 y1​和 y 2 y_2 y2​是独立的： p ( y 1 , y 2 ) = p 1 ( y 1 ) p 2 ( y 2 ) (10) p ( y _ { 1 } , y _ { 2 } ) = p _ { 1 } ( y _ { 1 } ) p _ { 2 } ( y _ { 2 } )\tag{10} p(y1​,y2​)=p1​(y1​)p2​(y2​)(10)  该定义自然适用于任意数量的 n n n个随机变量，在这种情况下，联合密度必须是 n n n个项的乘积。  该定义可用于推导独立随机变量的最重要属性。给定两个函数 h 1 h_1 h1​和 h 2 h_2 h2​，我们总是有： E { h 1 ( y 1 ) h 2 ( y 2 ) } = E { h 1 ( y 1 ) } E { h 2 ( y 2 ) } (11) E \{ h _ { 1 } ( y _ { 1 } ) h _ { 2 } ( y _ { 2 } ) \} = E \{ h _ { 1 } ( y _ { 1 } ) \} E \{ h _ { 2 } ( y _ { 2 } ) \}\tag{11} E{h1​(y1​)h2​(y2​)}=E{h1​(y1​)}E{h2​(y2​)}(11)  可以证明如下： E { h ( y 1 ) h 2 ( y 2 ) } = ∫ ∫ h 1 ( y 1 ) h 2 ( y 2 ) p ( y 1 , y 2 ) d y 1 d y 2 = ∫ ∫ h 1 ( y 1 ) p 1 ( y 1 ) h 2 ( y 2 ) p 2 ( y 2 ) d y 1 d y 2 = E { h 1 ( y 1 ) } E { h 2 ( y 2 ) } (12) \begin{aligned} E\{h( y _ { 1 } ) h _ { 2 } ( y _ { 2 } ) \} &= \int \int h _ { 1 } ( y _ { 1 } ) h _ { 2 } ( y _ { 2 } ) p ( y _ { 1 } , y _ { 2 } ) d y _ { 1 }dy_{2} \\ &= \int \int h _ { 1 } ( y _ { 1 } ) p _ { 1 } ( y _ { 1 } ) h _ { 2 } ( y _ { 2 } ) p _ { 2 } ( y _ { 2 } ) d y _ { 1 } d y _ { 2 } \\ &= E \{ h _ { 1 } ( y _ { 1 } ) \} E \{ h _ { 2 } ( y _ { 2 } ) \}\tag{12} \end{aligned} E{h(y1​)h2​(y2​)}​=∫∫h1​(y1​)h2​(y2​)p(y1​,y2​)dy1​dy2​=∫∫h1​(y1​)p1​(y1​)h2​(y2​)p2​(y2​)dy1​dy2​=E{h1​(y1​)}E{h2​(y2​)}​(12)

 独立性的较弱形式是不相关。如果两个随机变量 y 1 y_1 y1​和 y 2 y_2 y2​的协方差为零，则它们是不相关的： E { y 1 y 2 } − E { y 1 } E { y 2 } = 0 (13) E \{ y _ { 1 } y _ { 2 } \} - E \{ y _ { 1 } \} E \{ y _ { 2 } \} = 0\tag{13} E{y1​y2​}−E{y1​}E{y2​}=0(13) 如果变量是独立的，则它们是不相关的。取 h 1 ( y 1 ) = y 1 h_1(y_1)=y_1 h1​(y1​)=y1​和 h 2 ( y 2 ) = y 2 h_2(y_2)=y_2 h2​(y2​)=y2​则直接从等式(11)可以得出。  另一方面，不相关并不意味着独立。例如，假设 ( y 1 , y 2 ) (y_1,y_2) (y1​,y2​)是离散值，并且遵循这样的分布，即该对的概率1/4等于以下任意值：(0,1),(0,-1),(1,0),(-1,0)。可以简单地计算出 y 1 y_1 y1​和 y 2 y_2 y2​不相关。另一方面， E { y 1 2 y 2 2 } = 0 ≠ 1 4 = E { y 1 2 } E { y 2 2 } (14) E \{ y _ { 1 } ^ { 2 } y _ { 2 } ^ { 2 } \} = 0 \neq \frac { 1 } { 4 } = E \{ y _ { 1 } ^ { 2 } \} E \{ y _ { 2 } ^ { 2 } \}\tag{14} E{y12​y22​}=0​=41​=E{y12​}E{y22​}(14) 所以不符合等式(11)中的条件，变量不是独立的。  由于独立性包含不相关性，许多ICA方法都限制了估计过程，因此它总是给出独立分量的不相关性估计。这减少了可用参数的数量，并简化了问题。

  要了解为什么高斯变量使ICA不可能，我们假设混合矩阵是正交的，并且 是高斯的。那么 x 1 x_{1} x1​和 x 2 x_2 x2​是高斯的，不相关的并且具有单位方差。它们的联合密度由下式给出： p ( x 1 , x 2 ) = 1 2 π e x p ( − x 1 2 + x 2 2 2 ) (15) p({{x}_{1}},{{x}_{2}})=\frac{1}{2\pi }exp\left( -\frac{x_{1}^{2}+x_{2}^{2}}{2} \right)\tag{15} p(x1​,x2​)=2π1​exp(−2x12​+x22​​)(15)  此分布如图7所示。   该图显示密度是完全对称的。因此，它不包含有关混合矩阵A的列方向的任何信息。这就是为什么无法估计A的原因。  更严格地说，可以证明高斯变量 ( x 1 , x 2 ) (x_1,x_2) (x1​,x2​)的任何正交变换的分布具有与 ( x 1 , x 2 ) (x_1,x_2) (x1​,x2​)完全相同的分布，并且 x 1 x_1 x1​和 x 2 x_2 x2​是独立的。因此，在高斯变量的情况下，我们只能估计ICA模型直到正交变换。换句话说，对于高斯独立分量，矩阵A是不可识别的。实际上，如果仅独立分量之一是高斯分布的，则仍可以估计ICA模型。