Independent Components Analysis

接下来我们要讲的主体是独立成分分析（Independent Components Analysis，缩写为 ICA）。这个方法和主成分分析（PCA）类似，也是要找到一组新的基向量（basis）来表征（represent）样本数据。然而，这两个方法的目的是非常不同的。

还是先用“鸡尾酒会问题（cocktail party problem）”为例。在一个聚会场合中，有 个人同时说话，而屋子里的任意一个话筒录制到底都只是叠加在一起的这 个人的声音。但如果假设我们也有 个不同的话筒安装在屋子里，并且这些话筒与每个说话人的距离都各自不同，那么录下来的也就是不同的组合形式的所有人的声音叠加。使用这样布置的 个话筒来录音，能不能区分开原始的 个说话者每个人的声音信号呢？

把这个问题用方程的形式来表示，我们需要先假设有某个样本数据 ，这个数据是由 个独立的来源（independent sources）生成的。我们观察到的则为：

上面式子中的 是一个未知的正方形矩阵（square matrix），叫做混合矩阵（mixing matrix）。 通过重复的观察，我们就得到了训练集 %7D%20%3B%20i%20%3D%201%2C%20.%20.%20.%20%2C%20m%5C%7D#card=math&code=%5C%7Bx%5E%7B%28i%29%7D%20%3B%20i%20%3D%201%2C%20.%20.%20.%20%2C%20m%5C%7D&height=19&width=111)，然后我们的目的是恢复出生成这些样本 %7D%20%3D%20As%5E%7B(i)%7D#card=math&code=x%5E%7B%28i%29%7D%20%3D%20As%5E%7B%28i%29%7D&height=16&width=66) 的原始声音源 %7D#card=math&code=s%5E%7B%28i%29%7D&height=16&width=18) 。

在咱们的“鸡尾酒会问题”中，%7D#card=math&code=s%5E%7B%28i%29%7D&height=16&width=18) 就是一个 维度向量，而 %7D#card=math&code=s_j%5E%7B%28i%29%7D&height=23&width=18) 是第 个说话者在第 次录音时候发出的声音。%7D#card=math&code=x%5E%7B%28i%29%7D&height=16&width=20) 同样也是一个 维度向量，而 %7D#card=math&code=x_j%5E%7B%28i%29%7D&height=23&width=20)是第 个话筒在第 次录制到的声音。

设混合矩阵 的逆矩阵 是混合的逆向过程，称之为还原矩阵（unmixing matrix）。 那么咱们的目标就是找出这个 ，这样针对给定的话筒录音 %7D#card=math&code=x%5E%7B%28i%29%7D&height=16&width=20)，我们就可以通过计算 %7D%20%3D%20Wx%5E%7B(i)%7D#card=math&code=s%5E%7B%28i%29%7D%20%3D%20Wx%5E%7B%28i%29%7D&height=16&width=70) 来还原出来声音源。为了方便起见，我们就用 来表示 的第 行，这样就有：

这样就有 ，通过计算 %7D%20%3D%20w_j%5ET%20x%5E%7B(i)%7D#card=math&code=s_j%5E%7B%28i%29%7D%20%3D%20w_j%5ET%20x%5E%7B%28i%29%7D&height=23&width=74) 就可以恢复出第 个声源了。

13.1 独立成分分析（ICA）的模糊性（ambiguities）

能恢复到怎样的程度呢？如果我们对声源和混合矩阵都有预先的了解（prior knowledge），那就不难看出，混合矩阵 当中存在的某些固有的模糊性，仅仅给定了 %7D#card=math&code=x%5E%7B%28i%29%7D&height=16&width=20) 可能无法恢复出来。

例如，设 是一个 的置换矩阵（permutation matrix）。这就意味着矩阵 的每一行和每一列都只有一个 。下面就是几个置换矩阵的样例：

如果 是一个向量，那么 就是另外一个向量，这个向量包含了 坐标的置换版本（permuted version）。如果只给出 %7D#card=math&code=x%5E%7B%28i%29%7D&height=16&width=20)，是没有办法区分出 和 的。具体来说，原始声源的排列（permutation）是模糊的（ambiguous），这一点也不奇怪。好在大多数情况下，这个问题都并不重要。

进一步来说，就是没有什么办法能恢复出 的正确的缩放规模。例如，如果把 替换成了 ，那么每个 %7D#card=math&code=s%5E%7B%28i%29%7D&height=16&width=18) 都替换成了 s%5E%7B(i)%7D#card=math&code=%280.5%29s%5E%7B%28i%29%7D&height=19&width=46)，那么观测到的 %7D%20%3D%202A%20%C2%B7%20(0.5)s%5E%7B(i)%7D#card=math&code=x%5E%7B%28i%29%7D%20%3D%202A%20%C2%B7%20%280.5%29s%5E%7B%28i%29%7D&height=19&width=111) 还是跟原来一样的。再进一步说，如果 当中的某一列，都用一个参数 来进行缩放，那么对应的音源就被缩放到了 ，这也表明，仅仅给出 %7D#card=math&code=x%5E%7B%28i%29%7D&height=16&width=20)，是没办法判断这种情况是否发生的。因此，我们并不能还原出音源的“正确”缩放规模。然而，在我们应用的场景中，例如本文提到的这个“鸡尾酒会问题”中，这种不确定性并没有关系。具体来说，对于一个说话者的声音信号 %7D#card=math&code=s%5E%7B%28i%29%7D&height=16&width=18) 的缩放参数 只影响说话者声音的大小而已。另外，符号变换也没有影响，因为%7D#card=math&code=s_j%5E%7B%28i%29%7D&height=23&width=18) 和 %7D#card=math&code=-s_j%5E%7B%28i%29%7D&height=23&width=28) 都表示了扬声器中同样的声音大小。所以，如果算法找到的 被乘以任意一个非零数进行了缩放，那么对应的恢复出来的音源 也进行了同样的缩放；这通常都不要紧。（这些考量也适用于课堂上讨论的对 Brain/MEG 数据使用的 ICA 算法。）

上面这些是 ICA 算法模糊性的唯一来源么？还真是这样，只要声源 是非高斯分布（non-Gaussian）的即可。如果是高斯分布的数据（Gaussian data），例如一个样本中，，而 #card=math&code=s%5Csim%20N%280%2CI%29&height=16&width=66) 。（译者注：即 是一个以 和 为参数的正态分布，正态分布属于高斯分布） 其中的 是一个 的单位矩阵（identity matrix）。要注意，这是一个标准正态分布，其密度（density）轮廓图（contour）是以圆点为中心的圆，其密度是旋转对称的（rotationally symmetric）。

接下来，假如我们观测到了某个 ，其中的 就是混合矩阵（mixing matrix）。这样得到的 也是一个高斯分布的，均值为 ，协方差 。 然后设 为任意的正交矩阵（不太正式地说，也可以说成是旋转（rotation）矩阵或者是反射（reflection）矩阵），这样则有 ，然后设 。如果使用 而不是 作为混合矩阵，那么观测到的数据就应该是 。这个 也还是个高斯分布，依然是均值为 ，协方差为 %5ET%20%5D%20%3D%20E%5BA’ss%5ET%20(A’)%5ET%20%5D%20%3D%20E%5BARss%5ET%20(AR)%5ET%20%5D%20%3D%20ARR%5ET%20A%5ET%20%3D%20AA%5ET#card=math&code=E%5Bx%27%28x%27%29%5ET%20%5D%20%3D%20E%5BA%27ss%5ET%20%28A%27%29%5ET%20%5D%20%3D%20E%5BARss%5ET%20%28AR%29%5ET%20%5D%20%3D%20ARR%5ET%20A%5ET%20%3D%20AA%5ET&height=18&width=402)。看到没，无论混合矩阵使用 还是 ，得到的数据都是一个正态分布 #card=math&code=N%20%280%2C%20AA%5ET%20%29&height=18&width=63)（以 0 为均值，协方差为 ） 。这样就根本不能区分出来混合矩阵使用的是 还是 。所以，只要混合矩阵中有一个任意的旋转分量（arbitrary rotational component），并且不能从数据中获得，那么就不能恢复出原始数据源了。

上面这些论证，是基于多元标准正态分布（multivariate standard normal distribution）是旋转对称（rotationally symmetric）的这个定理。这些情况使得 ICA 面对高斯分布的数据（Gaussian data）的时候很无力，但是只要数据不是高斯分布的，然后再有充足的数据，那就还是能恢复出 个独立的声源的。

13.2 密度（Densities）和线性变换（linear transformations）

在继续去推导 ICA 算法之前，我们先来简要讲一讲对密度函数进行线性变换的效果（effect）。

加入我们有某个随机变量 ，可以根据某个密度函数 #card=math&code=p_s%28s%29&height=16&width=29) 来绘制。简单起见，咱们现在就把 当做是一个实数，即 。然后设 为某个随机变量，定义方式为 （其中 ）。然后设 是 的密度函数。那么这个 是多少呢？

设 。要计算 取某一个特定值的“概率（probability）”，可以先计算对于 ，在这一点上的 ，然后推导出%20%3D%20p_s(Wx)%E2%80%9D#card=math&code=%E2%80%9Cp_x%28x%29%20%3D%20p_s%28Wx%29%E2%80%9D&height=16&width=108)。然而，这是错误的。例如，假设 ，即其密度函数 %20%3D%201%5C%7B0%20%E2%89%A4%20s%20%E2%89%A4%201%5C%7D#card=math&code=p_s%28s%29%20%3D%201%5C%7B0%20%E2%89%A4%20s%20%E2%89%A4%201%5C%7D&height=16&width=124)。然后设 ，这样 。很明显， 在 这个区间均匀分布（distributed uniformly）。所以其密度函数也就是 %20%3D%20(0.5)1%5C%7B0%20%E2%89%A4%20x%20%E2%89%A4%202%5C%7D#card=math&code=p_x%28x%29%20%3D%20%280.5%291%5C%7B0%20%E2%89%A4%20x%20%E2%89%A4%202%5C%7D&height=16&width=156)。这并不等于 #card=math&code=p_s%28W%20x%29&height=16&width=45)，其中的 。所以正确的推导公式应该是 %20%3D%20p_s(Wx)%7CW%7C#card=math&code=p_x%28x%29%20%3D%20p_s%28Wx%29%7CW%7C&height=16&width=117)。

推广一下，若 是一个向量值的分布，密度函数为 ，而 ，其中的 是一个可逆的正方形矩阵，那么 的密度函数则为：

%20%3D%20p_s(Wx)%20%C2%B7%20%7CW%7C%0A#card=math&code=px%28x%29%20%3D%20p_s%28Wx%29%20%C2%B7%20%7CW%7C%0A&height=16&width=127)

上式中 。

Remark. If you’ve seen the result that A maps [0, 1]n to a set of volume |A|, then here’s another way to remember the formula for px given above, that also generalizes our previous 1-dimensional example.

备注。 可能你已经看到了用 映射 得到的就是一个由 组成的集合（译者注：这里的 volume 我不确定该怎么翻译），然后就又有了一个办法可以记住上面给出的关于 的公式了，这也是对之前讨论过的 维样例的一个泛化扩展。具体来说，设给定了 ，然后还按照惯例设 。接着设 是一个 维度超立方体，然后设 为由 给定的映射下的 的投影图像。这就是线性代数里面，用 来表示 的体积的标准结果，另外也是定义行列式（determinants）的一种方式。接下来，设 在 上均匀分布（uniformly distributed），这样其密度函数为 %20%3D%201%5C%7Bs%20%5Cin%20C_1%5C%7D#card=math&code=p_s%28s%29%20%3D%201%5C%7Bs%20%5Cin%20C_1%5C%7D&height=16&width=107)。然后很明显， 也是在 内均匀分布（uniformly distributed）。因此可以知道其密度函数为 %20%3D%201%5C%7Bx%20%5Cin%20C_2%5C%7D%2Fvol(C_2)#card=math&code=p_x%28x%29%20%3D%201%5C%7Bx%20%5Cin%20C_2%5C%7D%2Fvol%28C_2%29&height=16&width=161)，必须在整个 累积为（integrate to ，这是概率的性质）。但利用逆矩阵的行列式等于行列式的倒数这个定理，就有了 %20%3D%201%2F%7CA%7C%20%3D%20%7CA%5E%7B-1%7D%7C%20%3D%20%7CW%7C#card=math&code=1%2Fvol%28C_2%29%20%3D%201%2F%7CA%7C%20%3D%20%7CA%5E%7B-1%7D%7C%20%3D%20%7CW%7C&height=18&width=195)。所以则有 %20%3D%201%5C%7Bx%20%5Cin%20C_2%5C%7D%7CW%7C%20%3D%201%5C%7BWx%20%5Cin%20C_1%5C%7D%7CW%20%7C%20%3D%20p_s(W%20x)%7CW%20%7C#card=math&code=p_x%28x%29%20%3D%201%5C%7Bx%20%5Cin%20C_2%5C%7D%7CW%7C%20%3D%201%5C%7BWx%20%5Cin%20C_1%5C%7D%7CW%20%7C%20%3D%20p_s%28W%20x%29%7CW%20%7C&height=16&width=332)。

13.3 独立成分分析算法（ICA algorithm）

现在就可以推导 ICA 算法了。我们这里描述的算法来自于 Bell 和 Sejnowski，然后我们对算法的解释也是基于他们的算法，作为一种最大似然估计（maximum likelihood estimation）的方法。（这和他们最初的解释不一样，那个解释里面要涉及到一个叫做最大信息原则（infomax principal） 的复杂概念，考虑到对 ICA 的现代理解，推导过程已经不需要那么复杂了。）

我们假设每个声源的分布 都是通过密度函数 给出，然后联合分布 则为：

%3D%5Cprod%7Bi%3D1%7D%5En%20p_s(s_i)%0A#card=math&code=p%28s%29%3D%5Cprod%7Bi%3D1%7D%5En%20p_s%28s_i%29%0A&height=40&width=96)

这里要注意，通过在建模中将联合分布（joint distribution）拆解为边界分布（marginal）的乘积（product），就能得出每个声源都是独立的假设（assumption）。利用上一节推导的共识，这就表明对 的密度函数为：

%3D%5Cprod%7Bi%3D1%7D%5En%20p_s(w_i%5ET%20x)%5Ccdot%20%7Cw%7C%0A#card=math&code=p%28s%29%3D%5Cprod%7Bi%3D1%7D%5En%20p_s%28w_i%5ET%20x%29%5Ccdot%20%7Cw%7C%0A&height=40&width=137)

剩下的就只需要去确定每个独立的声源的密度函数 了。

回忆一下，给定某个实数值的随机变量 ，其累积分布函数（cumulative distribution function，cdf） 的定义为%3DP(z%5Cle%20z0)%3D%5Cint%7B-%5Cinfty%7D%5E%7Bz0%7Dp_z(z)dz#card=math&code=F%28z_0%29%3DP%28z%5Cle%20z_0%29%3D%5Cint%7B-%5Cinfty%7D%5E%7Bz_0%7Dp_z%28z%29dz&height=35&width=197)。然后，对这个累积分布函数求导数，就能得到 z 的密度函数：%20%3D%20F’(z)#card=math&code=p_z%28z%29%20%3D%20F%27%28z%29&height=17&width=79)。

因此，要确定 的密度函数，首先要做的就是确定其累积分布函数（cdf）。这个 函数必然是一个从 到 的单调递增函数。根据我们之前的讨论，这里不能选用高斯分布的 ，因为 ICA 不适用于高斯分布的数据。所以这里我们选择一个能够保证从 增长到 的合理的“默认（default）” 函数就可以了，比如 形函数（sigmoid function） %20%3D%201%2F(1%20%2B%20e%5E%7B-s%7D)#card=math&code=g%28s%29%20%3D%201%2F%281%20%2B%20e%5E%7B-s%7D%29&height=17&width=108)。这样就有，%20%3D%20g’(s)#card=math&code=p_s%28s%29%20%3D%20g%27%28s%29&height=17&width=75)。

1 如果你对声源的密度函数的形式有了事先的了解，那么在这个位置替换过来就是个很好的办法。不过如果没有这种了解，就可以用 形函数（sigmoid function），可以把这个函数当做是一个比较合理的默认函数，在很多问题中，这个函数用起来效果都不错。另外这里讲述的是假设要么所有的数据 %7D#card=math&code=x%5E%7B%28i%29%7D&height=16&width=20) 已经被证明均值为 ，或者可以自然预期具有 均值，比如声音信号就是如此。这很有必要，因为我们的假设 %20%3D%20g’(s)#card=math&code=p_s%28s%29%20%3D%20g%27%28s%29&height=17&width=75) 就意味着期望 （这个逻辑函数（logistic function）的导数是一个对称函数，因此给出的就是均值为 的随机变量对应的密度函数），这也意味着 。

是一个正方形矩阵，是模型中的参数。给定一个训练集合 %7D%3Bi%20%3D%201%2C…%2Cm%5C%7D#card=math&code=%5C%7Bx%5E%7B%28i%29%7D%3Bi%20%3D%201%2C…%2Cm%5C%7D&height=19&width=111)，然后对数似然函数（log likelihood）则为：

%3D%5Csum%7Bi%3D1%7D%5Em(%5Csum%7Bj%3D1%7D%5En%20log%20g’(wj%5ETx%5E%7B(i)%7D)%2Blog%7CW%7C))%0A#card=math&code=l%28W%29%3D%5Csum%7Bi%3D1%7D%5Em%28%5Csum_%7Bj%3D1%7D%5En%20log%20g%27%28w_j%5ETx%5E%7B%28i%29%7D%29%2Blog%7CW%7C%29%29%0A&height=41&width=235)

我们要做的就是上面这个函数找出关于 的最大值。通过求导，然后利用前面讲义中给出的定理 %5ET#card=math&code=%5Cnabla_W%7CW%7C%20%3D%20%7CW%7C%28W%5E%7B-1%7D%29%5ET&height=18&width=131)，就可以很容易推导出随机梯度上升（stochastic gradient ascent）学习规则（leaR^ning rule）。对于一个给定的训练样本 %7D#card=math&code=x%5E%7B%28i%29%7D&height=16&width=20)，这个更新规则为：

%7D)%20%5C%5C%0A1-2g(w_2%5ET%20x%5E%7B(i)%7D)%20%5C%5C%0A%5Cvdots%20%5C%5C%0A1-2g(w_n%5ET%20x%5E%7B(i)%7D)%0A%5Cend%7Bbmatrix%7Dx%5E%7B(i)T%7D%20%2B%20(W%5ET)%5E%7B-1%7D%0A%5Cend%7Bpmatrix%7D%0A#card=math&code=W%3A%3DW%2B%5Calpha%5Cbegin%7Bpmatrix%7D%0A%5Cbegin%7Bbmatrix%7D%0A1-2g%28w_1%5ET%20x%5E%7B%28i%29%7D%29%20%5C%5C%0A1-2g%28w_2%5ET%20x%5E%7B%28i%29%7D%29%20%5C%5C%0A%5Cvdots%20%5C%5C%0A1-2g%28w_n%5ET%20x%5E%7B%28i%29%7D%29%0A%5Cend%7Bbmatrix%7Dx%5E%7B%28i%29T%7D%20%2B%20%28W%5ET%29%5E%7B-1%7D%0A%5Cend%7Bpmatrix%7D%0A&height=88&width=305)

上式中的 是学习速率（leaR^ning rate）。

在算法收敛（converges）之后，就能计算出 %7D%20%3D%20Wx%5E%7B(i)%7D#card=math&code=s%5E%7B%28i%29%7D%20%3D%20Wx%5E%7B%28i%29%7D&height=16&width=70)，这样就能恢复出原始的音源了。

备注。 在写下数据的似然函数的时候，我们隐含地假设了这些 %7D#card=math&code=x%5E%7B%28i%29%7D&height=16&width=20) 都是彼此独立的（这里指的是对于不同的 值来说彼此独立；注意这个问题并不是说 %7D#card=math&code=x%5E%7B%28i%29%7D&height=16&width=20) 的不同坐标是独立的），这样对训练集的似然函数则为%7D%3BW)#card=math&code=%5Cprod_i%20p%28x%5E%7B%28i%29%7D%3BW%29&height=32&width=77)。很显然，对于语音数据和其他 %7D#card=math&code=x%5E%7B%28i%29%7D&height=16&width=20) 有相关性的时间序列数据来说，这个假设是不对的，但是这可以用来表明，只要有充足的数据，那么有相关性的训练样本并不会影响算法的性能。但是，对于成功训练的样本具有相关性的问题，如果我们把训练样本当做一个随机序列来进行访问，使用随机梯度上升（stochastic gradient ascent）的时候，有时候也能帮助加速收敛。（也就是说，在训练集的一个随机打乱的副本中运行随机梯度上升。）