fMRI解析の基礎 (10)：主成分分析

最終更新日: 2020年7月19日

fMRI解析の基礎シリーズ

主成分分析（PCA）は、多変量統計手法の１つで、他の手法の前処理としてノイズや次元を削減するために使われています。

多変量正規分布に従うデータ $[x_i,y_i]$ を考えます。

[x_i,y_i] \sim N \left([\mu_x,\mu_y], \begin{bmatrix} \sigma_x^2 & cov_{xy} \\ cov_{xy} & \sigma_y^2 \end{bmatrix} \right)

$[x_i,y_i]$ に相関がある場合は、合わせて１つの指標として用いることが次元削減に繋がります。まずは、分散共分散行列を固有値分解します。

\begin{bmatrix} \sigma_x^2 & cov_{xy} \\ cov_{xy} & \sigma_y^2 \end{bmatrix} = \begin{bmatrix} \bold{v}_1 & \bold{v}_2\end{bmatrix} \begin{bmatrix} d_1 & 0 \\ 0 & d_2 \end{bmatrix} \begin{bmatrix} \bold{v}_1^T \\ \bold{v}_2^T \end{bmatrix}

分散共分散行列は半正定値なので、固有値 $d$ は正となります。固有ベクトル $v$ は、データがつくる楕円形の等高線の軸に直交し、それら同士も直交します。固有ベクトル（主成分）の方向に新しい軸をとれば、多数の値を１つの軸で説明できるわけです。

固有値は、固有値ベクトルの方向における分散と一致します。ということは、各軸における情報量（寄与率）は、固有値に比例することになります。

第 r 主成分までの累積寄与率は、以下の式で求められます。これは、第 r 主成分までに説明しうる情報量の割合を示しています。

\frac{d_1+d_2+\ldots+d_r}{d_1+d_2+\ldots+d_n}

$y_{j}^{(i)}$ が $i$ 番目のTRの $j$ 番目のボクセルにおける値だとします。PCAを適用するときは、各ボクセルの平均を0にした行列 $Y_C$ を用います。

Y_C=\begin{bmatrix} y_{1}^{(1)} & \ldots & y_{N_v}^{(1)} \\ \vdots & \ddots & \vdots \\ y_{1}^{(N_{TR})} & \ldots & y_{N_v}^{(N_{TR})} \end{bmatrix}- \begin{bmatrix} 1 \\ \vdots \\ 1 \end{bmatrix} \begin{bmatrix} \bar{y}_{1} & \cdots & \bar{y}_{N_v}\end{bmatrix}

これを用いて、分散共分散行列 $\Sigma \in \R^{N_v \times N_v}$ を推定できます。

\hat{\Sigma} = \frac{1}{N_{TR}-1}Y_C^TY_C

$\hat{\Sigma}$ を固有値分解して固有ベクトルを算出し、その方向に射影します。 $\hat{\Sigma}$ のランクは $Y_C$ と一致するので、 $0$ ではない固有値 $d$ は最大でTRの数 $N_{TR}$ だけあります。

\hat{\Sigma} = \begin{bmatrix} \bold{v}_1 & \cdots & \bold{v}_{N_v}\end{bmatrix} \begin{bmatrix} d_1 & \cdots & 0 \\ \vdots & \ddots & \vdots \\ 0 & \cdots & d_{N_v} \end{bmatrix} \begin{bmatrix} \bold{v}_1^{\top} \\ \vdots \\ \bold{v}_{N_v}^{\top} \end{bmatrix}

z_{j}^{(i)}=\bar{y}_j+d_i v_{ij}

タスクに敏感な信号は固有値が大きく、ノイズは固有値が小さくなると期待するわけです。

固有値が小さい要素を除くことで、ノイズを削減できます。

$Y$ を正規化して $Y_{c_s}$ とし、分散共分散行列を計算し、固有値分解します。しかし、分散共分散行列 $\Sigma_v=Y_{c_s}^\top Y_{c_s} \in \R^{N_v \times N_v}$ は大きいので、固有値を出すためなら

T=Y_{c_s}Y_{c_s}^\top=UDU^\top \in \R^{N_{TR} \times N_{TR}}

と近似できます。最初の r 個を選択して $U_r \in \R^{N_{TR} \times r}$ とし、これを用いて $Y_C$ を射影します。

Y_{C(r)}=U_r^\top Y_C \in \R^{r \times N_{TR}}

この $Y_{C(r)}$ が、ノイズを削減した成分となります。

ノイズ削減後に元の空間に戻して、始めに引いた平均を加算してあげます。

\hat{Y}_r=U_rU_r^\top Y_C+ \bold{1} \bold{\bar{y}}^\top

固有値を大きい順に並べて、Kinkと呼ばれる歪みの位置を基準にします（スクリープロット）。視覚的に見つからなければ、累積寄与率が90％となるようにします。

特異値分解 (SVD) がを使うと、分散共分散行列を直接必要としません。特異値分解は、次のように表されます。

Y_C=UGV^\top

ただし、 $U,V$ は直交行列で、 $G$ は対角行列と零行列を並べたものです。

このとき、 $T$ は以下のように簡単に計算できます。

\begin{array}{rcl} T &=& Y_{C_s}Y_{C_s}^\top \\ \ \\ &=& \frac{1}{N_{TR}-1}Y_CY_C^\top \\ \ \\ &=& \frac{1}{N_{TR}-1}(UGV^\top )(VG^\top U) \\ \ \\ &=& \frac{1}{N_{TR}-1}UGG^\top U \\ \ \\ &=& U \begin{bmatrix} \frac{1}{N_{TR}-1}b^2 & 0\\ 0 & 0 \end{bmatrix}U^\top \end{array}

SPMでは、SVDを適用して $U$ から相関を取ってきているようです。