独立成分分析

Independent Component Analysis (ICA) と呼ばれる多変量解析の手法です．だいたい90年代頃に確立されました．まずどういったものかと言うと，観測された信号を独立な複数の信号の線形な重ね合わせとして再表現するものです．

といってもよくわからないと思うので，我々も実は日頃からICAをやっていますよという話からしましょう．

• お気持ち
• 求め方

お気持ち

皆さん毎日夜には駅前の居酒屋やバーで楽しくお酒を嗜んでいる事かと思いますが，実はこの時我々の脳は独立成分分析をしているのです．

飲み屋では多くの人間が同時に声を発し，食器の音やどこかのコールの音，厨房の音と無数の音が同時に我々の耳=脳に送られてきます．しかしどういうわけか，その音がどこの誰の声か，厨房の音なのか食器の音なのか判断できますよね．つまり，与えられた音の時系列データ(これは一つのデータに重ねられて聴こえてる) を複数の信号源に分解しているわけです．

我々は特に意識せずとも当たり前のように出来るこの作業ですが，実はPCにやらせるととんでもなく難しかったのです．それをどうにか実現しましょうという事で出来たのが現在の ICA です．

脳波の研究に何の役に立つのかですが，端的に言えばノイズの除去や脳活動の特徴抽出です．脳波はだいたい 32 とか64ch の電極を使って計測するわけですが，そこには残念ながら眼球運動由来の電位や筋肉の電位など，脳波ではない成分も乗ってしまっています．これらは脳波に比べて電位も大きいので，脳波の解析の際に邪魔以外のなにものでもなく，親の仇のように憎むべき存在です．

もう分かるかと思いますが， "脳波" 信号を分解し，真の脳波信号と眼電，筋電，などに分解してあげるのに使えるのが独立成分分析です．実用では，こうして分解したデータのうち眼電や筋電由来と思われる成分だけ取り除いてあげれば綺麗に脳波だけ見れるわけですね．すくなくとも目的はそんなとこです．

求め方

長くなりましたが日本語はこれくらいにして，早速数理に入りましょう．まず元の時系列信号を $\mathbf{x}(t)$ とします．これが，初めに仮定より線形な $(N \geq 2)$ 個の信号源 ($\mathbf{s}(t)=[s_i(t), s_2(t),...s_N(t)]$) からの信号の重ね合わせによってN個のセンサに観測されているとします．ここで信号源とセンサーの個数が一緒なのに違和感を覚える人もいると思いますが，とりあえず置いてください．実際脳波のICAも基本的にチャンネル数と同じ数に分解します．するとこの関係は

$$\mathbf{x}(t) = A \mathbf{s}(t)^T$$

と表せます．ここで A は $N\times N$ の係数行列で，$a_{ij}$ は $i$ 番目のセンサで観測される $j$ 番目の信号源からの信号の係数です．この値が大きいのは信号源の影響を強く受け，小さいのはあまり受けていない事を意味します．

ではここで，今回の問題で求めたいのは $\mathbf{s}(t)$ ですよね．なので $\hat{\mathbf{s}}(t)$ 的な行列として $\mathbf{y}(t)$ を定義し，

$$\mathbf{y}(t) = \mathbf{W} \mathbf{x}(t)$$

という式をたてます．実数の分離行列$\mathbf{W}$を元信号にかけて $\mathbf{y}(t) = [y_1(t), y_2(t),...,y_N(t)]$ で表される分離信号に分解する事を考えます．分離行列も $N\times N$ の成分をもっていて，それぞれのセンサの値をどれだけ反映させて分離信号を生成するのか，を分離信信号の数だけ行います．この $\mathbf{y}(t)$ が綺麗に独立するように $\mathbf{W}$ を更新していくのが独立成分分析の手順です．

ではその独立性をどう表すか，というところですが，$\mathbb{D}_{KL}$ を使います．正確には相互情報量です．相互情報量は，複数の変数があったときにその同時分布と周辺分布の積との間のKL距離を測るものでしたね．両者が全く同じ分布である時，0 になるものだったので，相互情報量が0になるときは変数が全て独立である，という事を表せるのでした．これを使います．

$$p (\mathbf{y}) = p(y_1, y_2,...,y_N) \qquad \text{同時分布}\\ p (\mathbf{y}) = \Pi_{i=1}^N p(y_i) \qquad \text{周辺分布}$$

こいつらのKL距離が相互情報量だから

$$\mathbb{D}_{KL}(\mathbf{y}) = \int p(\mathbf{y}) \log \frac{p(\mathbf{y})}{\Pi_{i=1}^N p(y_i)} d\mathbf{y} \tag{1}$$

と表せます．この時，式 (1) が0になるのは $p(y_1, y_2,...,y_N) =\Pi_{i=1}^N p(y_i)$ が成り立つ，つまり独立な信号源として仮定した$\mathbf{y}$がちゃんと独立である時になります．よって最小化しましょう．

$$\mathbf{\hat W} = \argmin_{\mathbf{W}} \mathbb{D}_{KL}(\mathbf{y}|\mathbf{W}) \nonumber \\ = \argmin_{\mathbf{W}} \int p(\mathbf{y}|\mathbf{W}) \log \frac{p(\mathbf{y}|\mathbf{W})}{\Pi_{i=1}^N p(y_i|\mathbf{W})} d\mathbf{y} \tag{2}$$

式(1)が0になるように分離行列$\mathbf{W}$を逐次更新していけばいいわけですね．最小化する逐次更新の仕方はいくつかあって，最有力?なのが勾配法を用いるものです．勾配法についてまだちゃんと勉強できていないので今は式だけ．

$$\mathbf{W} \leftarrow \mathbf{W} - \eta \mathbf{E}[\phi (\mathbf{y})\mathbf{y}^T - I]\mathbf{W}^{-T} \tag{3}$$

この更新式(2)に従って評価式 (1) を更新していけば相互情報量を最小化できて，無事に元の$\mathbf{x}(t)$を独立成分$\mathbf{y}(t)$を取り出すことが出来るわけですね．正直最後だけまだ分からん．でも最小化したいってのは分かるのでとりあえずヨシ！

Math 基礎数学編 掛け算 三角関数 複素数 微分・積分 基礎統計 信号処理編 フーリエ変換 wavelet変換 ヒルベルト変換 基本の脳波解析 位相同期解析 統計編 確率分布 大数の法則と中心極限定理 統計量と標本分布 統計的検定 分散分析	Analysis EEGLAB 環境構築 データのインポート 基本的な下処理 発展的な下処理 単被験者での解析 被験者群での解析	Experiment	Simulations 環境構築 Python環境構築 pythonでのGPUセットアップ Jupyterセットアップ Juliaセットアップ 非線形力学 非線形力学系とは 線形安定性解析 非線形系の安定性解析