## 2. 準備：中心化と標本共分散行列 $\boldsymbol{x}_1, \ldots, \boldsymbol{x}_n \in \mathbb{R}^p$（$n \ge 2$）を観測データとする。標本平均（重心）を $$ \bar{\boldsymbol{x}} = \frac{1}{n}\sum_{i=1}^{n} \boldsymbol{x}_i $$ とおき、中心化ベクトルを $\boldsymbol{y}_i = \boldsymbol{x}_i - \bar{\boldsymbol{x}}$ と定める。$\boldsymbol{y}_i^{\mathsf{T}}$ を第 $i$ 行にもつ行列 $X \in \mathbb{R}^{n \times p}$ を**中心化データ行列**と呼ぶ。さらに $$ S = \frac{1}{n} X^{\mathsf{T}} X = \frac{1}{n}\sum_{i=1}^{n} \boldsymbol{y}_i \boldsymbol{y}_i^{\mathsf{T}} \in \mathbb{R}^{p \times p} $$ を**標本共分散行列**と呼ぶ。その $(j, l)$ 成分は $s_{jl} = \frac{1}{n}\sum_{i=1}^n y_{ij} y_{il}$ であり、$s_{jj}$ は第 $j$ 変数の標本分散、$s_{jl}$ は第 $j$ 変数と第 $l$ 変数の標本共分散である。定義の中で $\frac{1}{n} X^{\mathsf{T}} X$ と $\frac{1}{n}\sum_i \boldsymbol{y}_i \boldsymbol{y}_i^{\mathsf{T}}$ を同一視しました。これは行列積の定義から出ます。$X^{\mathsf{T}}$ の第 $i$ 列が $\boldsymbol{y}_i$ なので、$(X^{\mathsf{T}} X)_{jl} = \sum_{i=1}^{n} X_{ij} X_{il} = \sum_{i=1}^{n} y_{ij} y_{il} = \left(\sum_i \boldsymbol{y}_i \boldsymbol{y}_i^{\mathsf{T}}\right)_{jl}$ となります。不偏分散を使う流儀では $\frac{1}{n-1}$ で割ります。この行列を $S'$ と書くと $S' = \frac{n}{n-1} S$ のように両者は正の定数倍しか違わないので、固有ベクトル（＝主成分方向）は完全に一致し、固有値は同じ定数倍だけずれます。寄与率は固有値の比なので、これも一致します。以下では表記が軽くなる $\frac{1}{n}$ を採ります。 PCA が見るのはただ一つの量、「ある方向へ射影したときの分散」です。これを定義します。 $\boldsymbol{u} \in \mathbb{R}^p$ を $\|\boldsymbol{u}\| = 1$ を満たすベクトルとする。第 $i$ 個体の $\boldsymbol{u}$ 方向の**得点**を $$ z_i(\boldsymbol{u}) = \langle \boldsymbol{u}, \boldsymbol{y}_i \rangle = \boldsymbol{u}^{\mathsf{T}}(\boldsymbol{x}_i - \bar{\boldsymbol{x}}) $$ と定め、その標本分散を $\boldsymbol{u}$ 方向の**射影分散**と呼ぶ。得点の平均は $\frac{1}{n}\sum_i z_i(\boldsymbol{u}) = \boldsymbol{u}^{\mathsf{T}}\left(\frac{1}{n}\sum_i \boldsymbol{y}_i\right) = \boldsymbol{u}^{\mathsf{T}}\boldsymbol{0} = 0$ です（中心化しているので $\sum_i \boldsymbol{y}_i = \boldsymbol{0}$）。したがって射影分散は二乗平均そのもので、 $$ \frac{1}{n}\sum_{i=1}^{n} z_i(\boldsymbol{u})^2 = \frac{1}{n}\sum_{i=1}^{n} \boldsymbol{u}^{\mathsf{T}} \boldsymbol{y}_i \, \boldsymbol{y}_i^{\mathsf{T}} \boldsymbol{u} = \boldsymbol{u}^{\mathsf{T}} \left( \frac{1}{n}\sum_{i=1}^{n} \boldsymbol{y}_i \boldsymbol{y}_i^{\mathsf{T}} \right) \boldsymbol{u} = \boldsymbol{u}^{\mathsf{T}} S \boldsymbol{u} $$ と書けます。途中で使ったのは、スカラー $z_i = \boldsymbol{u}^{\mathsf{T}}\boldsymbol{y}_i$ について $z_i^2 = z_i z_i^{\mathsf{T}} = (\boldsymbol{u}^{\mathsf{T}}\boldsymbol{y}_i)(\boldsymbol{y}_i^{\mathsf{T}}\boldsymbol{u})$ という書き換えと、$\boldsymbol{u}$ が $i$ に依存しないので和の外に出せることだけです。**この式が記事全体の出発点です**。データの散らばりを方向 $\boldsymbol{u}$ ごとに測る関数が、行列 $S$ の二次形式として書けてしまう。これで問題は線形代数に移りました。 $S$ がどんな行列なのかを押さえておきます。の $S$ は $S^{\mathsf{T}} = S$ を満たし、さらに任意の $\boldsymbol{u} \in \mathbb{R}^p$ に対して $\boldsymbol{u}^{\mathsf{T}} S \boldsymbol{u} \ge 0$ である。対称性: 転置の性質 $(AB)^{\mathsf{T}} = B^{\mathsf{T}}A^{\mathsf{T}}$ と $(A^{\mathsf{T}})^{\mathsf{T}} = A$ から $$ S^{\mathsf{T}} = \left(\frac{1}{n}X^{\mathsf{T}}X\right)^{\mathsf{T}} = \frac{1}{n}X^{\mathsf{T}}(X^{\mathsf{T}})^{\mathsf{T}} = \frac{1}{n}X^{\mathsf{T}}X = S . $$ 半正定値性: 任意の $\boldsymbol{u}$ に対し $$ \boldsymbol{u}^{\mathsf{T}} S \boldsymbol{u} = \frac{1}{n}\boldsymbol{u}^{\mathsf{T}}X^{\mathsf{T}}X\boldsymbol{u} = \frac{1}{n}(X\boldsymbol{u})^{\mathsf{T}}(X\boldsymbol{u}) = \frac{1}{n}\|X\boldsymbol{u}\|^2 \ge 0 . $$ ノルムの二乗は非負なので結論を得ます。なお $\|\boldsymbol{u}\| = 1$ のときこの量は射影分散に等しく、「分散は負にならない」という当たり前の事実の言い換えになっています。対称行列については、次の定理が使えます。これが PCA の理論的な土台です。 $S \in \mathbb{R}^{p \times p}$ が $S^{\mathsf{T}} = S$ を満たすとする。このとき $\mathbb{R}^p$ の正規直交基底 $\boldsymbol{v}_1, \ldots, \boldsymbol{v}_p$ と実数 $\lambda_1 \ge \lambda_2 \ge \cdots \ge \lambda_p$ が存在して、すべての $j$ について $S\boldsymbol{v}_j = \lambda_j \boldsymbol{v}_j$ が成り立つ。行列で書けば、$V = [\boldsymbol{v}_1\ \cdots\ \boldsymbol{v}_p]$、$\Lambda = \operatorname{diag}(\lambda_1, \ldots, \lambda_p)$ として $V^{\mathsf{T}}V = I_p$ かつ $S = V \Lambda V^{\mathsf{T}}$ である。さらに $S$ が半正定値ならば $\lambda_p \ge 0$、すなわちすべての固有値は非負である。前半（正規直交固有基底の存在と固有値の実数性）は本記事では認めて使います。所在はに示します。ここでは後半だけ示します。$S$ が半正定値のとき、$\boldsymbol{v}_j$ は単位ベクトルなので $$ \lambda_j = \lambda_j \|\boldsymbol{v}_j\|^2 = \lambda_j \boldsymbol{v}_j^{\mathsf{T}}\boldsymbol{v}_j = \boldsymbol{v}_j^{\mathsf{T}}(\lambda_j \boldsymbol{v}_j) = \boldsymbol{v}_j^{\mathsf{T}} S \boldsymbol{v}_j \ge 0 $$ となります。最後の不等号が半正定値性です。とくにより標本共分散行列の固有値はすべて非負です。スペクトル定理の証明は [スペクトル定理](/mathematics/linear-algebra/spectral-theorem) にあります。ここで使う実対称行列の場合はがそのままの形です。固有値・固有ベクトルの定義と特性多項式による計算は [固有値と固有ベクトル](/mathematics/linear-algebra/eigenvalues)（）、直交射影と正規直交基底の扱いは [内積空間とグラム・シュミット直交化](/mathematics/linear-algebra/inner-product-spaces)（）を参照してください。以下では $\boldsymbol{v}_1, \ldots, \boldsymbol{v}_p$ と $\lambda_1 \ge \cdots \ge \lambda_p \ge 0$ を、この定理が与える標本共分散行列 $S$ の正規直交固有系と固有値（降順）として固定します。 ## 3. 第 1 主成分：分散最大化がレイリー商になる PCA が解きたい問題は、§1 の設計原理をそのまま式にしたものです。 $$ \text{最大化} \quad \boldsymbol{u}^{\mathsf{T}} S \boldsymbol{u} \quad \text{（条件 } \|\boldsymbol{u}\| = 1 \text{ のもとで）} $$ なぜ $\|\boldsymbol{u}\| = 1$ という条件を付けるのでしょうか。付けないと問題が壊れるからです。$\boldsymbol{u}$ を $c$ 倍すると得点も $c$ 倍になり、射影分散は $(c\boldsymbol{u})^{\mathsf{T}}S(c\boldsymbol{u}) = c^2 \boldsymbol{u}^{\mathsf{T}}S\boldsymbol{u}$ と $c^2$ 倍になります。$S \ne O$ なら $c \to \infty$ で発散し、最大値が存在しません。私たちが知りたいのは「どの**向き**にデータが伸びているか」であって、物差しの長さではありません。長さを $1$ に固定すれば、比較しているのは向きだけになります。 $S \in \mathbb{R}^{p\times p}$ を実対称行列とし、の記号を使う。このとき $$ \max_{\|\boldsymbol{u}\| = 1} \boldsymbol{u}^{\mathsf{T}} S \boldsymbol{u} = \lambda_1 $$ が成り立ち、最大値は $\boldsymbol{u} = \boldsymbol{v}_1$ で達成される。さらに、最大値を達成する単位ベクトルの全体は、固有値 $\lambda_1$ の固有空間 $\ker(S - \lambda_1 I)$ に含まれる単位ベクトルの全体と一致する。 **ステップ 1（成分表示）.** より $\boldsymbol{v}_1, \ldots, \boldsymbol{v}_p$ は $\mathbb{R}^p$ の正規直交基底なので、任意の $\boldsymbol{u}$ は $\boldsymbol{u} = \sum_{j=1}^p c_j \boldsymbol{v}_j$（$c_j = \langle \boldsymbol{v}_j, \boldsymbol{u}\rangle$）と一意に書けます。正規直交性から $$ \|\boldsymbol{u}\|^2 = \left\langle \sum_j c_j \boldsymbol{v}_j, \sum_l c_l \boldsymbol{v}_l \right\rangle = \sum_{j,l} c_j c_l \langle \boldsymbol{v}_j, \boldsymbol{v}_l\rangle = \sum_{j=1}^{p} c_j^2 $$ なので、条件 $\|\boldsymbol{u}\| = 1$ は $\sum_j c_j^2 = 1$ と同値です。 **ステップ 2（二次形式の計算）.** $S\boldsymbol{u} = \sum_j c_j S\boldsymbol{v}_j = \sum_j c_j \lambda_j \boldsymbol{v}_j$ なので、同じく正規直交性から $$ \boldsymbol{u}^{\mathsf{T}} S \boldsymbol{u} = \left\langle \sum_l c_l \boldsymbol{v}_l, \sum_j c_j \lambda_j \boldsymbol{v}_j \right\rangle = \sum_{j=1}^{p} \lambda_j c_j^2 . $$ つまり二次形式は、固有値を重み $c_j^2$ で平均したものです（重みの和はステップ 1 より $1$）。 **ステップ 3（上からの評価）.** すべての $j$ で $\lambda_j \le \lambda_1$ かつ $c_j^2 \ge 0$ なので $\lambda_j c_j^2 \le \lambda_1 c_j^2$、したがって $$ \boldsymbol{u}^{\mathsf{T}} S \boldsymbol{u} = \sum_j \lambda_j c_j^2 \le \lambda_1 \sum_j c_j^2 = \lambda_1 . $$ **ステップ 4（達成）.** $\boldsymbol{u} = \boldsymbol{v}_1$ は単位ベクトルであり、$\boldsymbol{v}_1^{\mathsf{T}} S \boldsymbol{v}_1 = \boldsymbol{v}_1^{\mathsf{T}}(\lambda_1 \boldsymbol{v}_1) = \lambda_1$。よって上限 $\lambda_1$ は実際に達成され、最大値です。 **ステップ 5（等号成立の条件）.** ステップ 3 の不等式の差を取ると $$ \lambda_1 - \boldsymbol{u}^{\mathsf{T}} S \boldsymbol{u} = \sum_{j=1}^{p} (\lambda_1 - \lambda_j) c_j^2 $$ であり、右辺は非負項の和です。したがって等号が成り立つのは、$\lambda_j < \lambda_1$ であるすべての $j$ について $c_j = 0$ となるとき、かつそのときに限ります。これは $\boldsymbol{u}$ が $\lambda_1$ に属する固有ベクトルたちの張る空間 $\ker(S - \lambda_1 I) = \operatorname{span}\{\boldsymbol{v}_j : \lambda_j = \lambda_1\}$ に入ることと同値です。この定理をデータの言葉に翻訳すれば、第 1 主成分が得られます。の正規直交固有系 $\boldsymbol{v}_1, \ldots, \boldsymbol{v}_p$（固有値の降順）に対し、$\boldsymbol{v}_k$ を**第 $k$ 主成分方向**と呼ぶ。第 $i$ 個体の**第 $k$ 主成分得点**を $$ z_{ik} = \langle \boldsymbol{v}_k, \boldsymbol{y}_i \rangle = \boldsymbol{v}_k^{\mathsf{T}}(\boldsymbol{x}_i - \bar{\boldsymbol{x}}) $$ と定める。第 $k$ 主成分得点を並べた $n$ 個の数の組を第 $k$ 主成分と呼ぶこともある。すべての単位ベクトル $\boldsymbol{u}$ のうち、射影分散 $\frac{1}{n}\sum_i z_i(\boldsymbol{u})^2$ を最大にするのは $S$ の最大固有値 $\lambda_1$ に属する単位固有ベクトルであり、そのときの射影分散の値は $\lambda_1$ に等しい。 §2 で示したとおり射影分散は $\boldsymbol{u}^{\mathsf{T}} S \boldsymbol{u}$ に等しく、より $S$ は実対称なのでがそのまま適用できます。最大値は $\lambda_1$、達成するのは $\lambda_1$ の固有空間内の単位ベクトルです。条件付き最大化ですから、[ラグランジュの未定乗数法](/mathematics/calculus/multivariable-differentiation)でも同じ結論に達します。$L(\boldsymbol{u}, \mu) = \boldsymbol{u}^{\mathsf{T}}S\boldsymbol{u} - \mu(\boldsymbol{u}^{\mathsf{T}}\boldsymbol{u} - 1)$ とおきます。$S$ が対称であることを使うと $$ \frac{\partial}{\partial u_k}\sum_{j,l} S_{jl}u_j u_l = \sum_{l} S_{kl}u_l + \sum_{j} S_{jk}u_j = 2 (S\boldsymbol{u})_k $$ なので、$\nabla_{\boldsymbol{u}} L = 2S\boldsymbol{u} - 2\mu\boldsymbol{u} = \boldsymbol{0}$、すなわち $S\boldsymbol{u} = \mu\boldsymbol{u}$ が停留条件です。**制約付き最適化の停留点は、そのまま固有方程式になります**。しかも停留点では $\boldsymbol{u}^{\mathsf{T}}S\boldsymbol{u} = \mu\boldsymbol{u}^{\mathsf{T}}\boldsymbol{u} = \mu$ なので、目的関数の値は乗数そのものです。よって最大を与えるのは最大固有値。ただしこの議論だけでは「最大値が存在すること」を別途言う必要があります（単位球面はコンパクト、二次形式は連続なので最大値は存在します）。の証明はこの点も込みで初等的に片付いています。主成分方向には二つの不定性があります。第一に符号です。$\boldsymbol{v}$ が単位固有ベクトルなら $-\boldsymbol{v}$ もそうで、得点の符号がすべて反転するだけですから、両者は同じ主成分を表します。数値ライブラリの出力の符号が実行環境によって変わることがあるのはこのためで、異常ではありません。第二に、固有値が重複すると固有空間が 2 次元以上になり、その中のどの正規直交基底を採っても構わなくなります（）。「第 1 主成分方向は $\boldsymbol{u}$ である」と一意に言えるのは $\lambda_1 > \lambda_2$ のときだけです。 $S = \sigma^2 I_p$（$\sigma > 0$）となるデータを考えます。どの方向を選んでも射影分散は $\boldsymbol{u}^{\mathsf{T}}(\sigma^2 I)\boldsymbol{u} = \sigma^2\|\boldsymbol{u}\|^2 = \sigma^2$ で同じです。固有値は $\sigma^2$ のみ（重複度 $p$）で、任意の正規直交基底が主成分系になります。$k$ 次元に落とすと、どの部分空間を選んでも $(p-k)\sigma^2$ の分散を失います。次元削減で得をする余地がまったくない状況です。実データでも上位の固有値が接近していると、標本のわずかな揺らぎで主成分方向が大きく回転します。「第 2 主成分は $\lambda_2 \approx \lambda_3$ のとき解釈しない」というのは実務上の鉄則です。 ## 4. 第 k 主成分と主成分得点の性質第 1 主成分だけでは足りないとき、次の軸をどう選ぶか。分散が大きい方向をもう一度探すと、また $\boldsymbol{v}_1$ の近くが選ばれてしまい、同じ情報を二度記録することになります。そこで「すでに採った方向と直交する」という条件を課します。直交する方向の得点は、後で見るように互いに無相関になり、情報の重複が起きません。の記号のもとで、$1 \le k \le p$ とする。集合 $$ W_k = \{\boldsymbol{u} \in \mathbb{R}^p : \|\boldsymbol{u}\| = 1,\ \langle \boldsymbol{u}, \boldsymbol{v}_1\rangle = \cdots = \langle \boldsymbol{u}, \boldsymbol{v}_{k-1}\rangle = 0 \} $$ （$k = 1$ のときは単位球面全体）の上で、$\boldsymbol{u}^{\mathsf{T}}S\boldsymbol{u}$ の最大値は $\lambda_k$ であり、$\boldsymbol{u} = \boldsymbol{v}_k$ で達成される。 $\boldsymbol{u} \in W_k$ を取り、の証明のステップ 1 と同様に $\boldsymbol{u} = \sum_{j=1}^p c_j \boldsymbol{v}_j$、$c_j = \langle \boldsymbol{v}_j, \boldsymbol{u}\rangle$ と展開します。$W_k$ の定義から $c_1 = \cdots = c_{k-1} = 0$ です。したがって $\sum_{j \ge k} c_j^2 = \|\boldsymbol{u}\|^2 = 1$ であり、ステップ 2 の計算から $$ \boldsymbol{u}^{\mathsf{T}} S \boldsymbol{u} = \sum_{j = k}^{p} \lambda_j c_j^2 \le \lambda_k \sum_{j=k}^{p} c_j^2 = \lambda_k $$ を得ます。不等号では、$j \ge k$ のとき固有値が降順に並んでいることから $\lambda_j \le \lambda_k$ であることを使いました。一方 $\boldsymbol{v}_k$ は単位ベクトルで、$j \le k-1$ について $\langle \boldsymbol{v}_k, \boldsymbol{v}_j \rangle = 0$（正規直交性）なので $\boldsymbol{v}_k \in W_k$ であり、$\boldsymbol{v}_k^{\mathsf{T}}S\boldsymbol{v}_k = \lambda_k$。よって上限は達成され、最大値です。つまり「分散が最大の方向を採る → それに直交する中で分散が最大の方向を採る → …」という素朴な逐次手続きが、固有値を大きい順に並べる操作とぴったり一致します。で固有値の降順に番号を付けたのは、この一致を保つためです。の主成分得点 $z_{ik}$ について、次が成り立つ。 1. 各 $k$ について $\frac{1}{n}\sum_{i=1}^{n} z_{ik} = 0$。 2. 各 $k, l$ について $\frac{1}{n}\sum_{i=1}^{n} z_{ik} z_{il} = \lambda_k \delta_{kl}$。すなわち第 $k$ 主成分得点の標本分散は $\lambda_k$ であり、異なる主成分得点どうしの標本共分散は $0$ である。 3. $\displaystyle\sum_{k=1}^{p} \lambda_k = \operatorname{tr} S = \sum_{j=1}^{p} s_{jj} = \frac{1}{n}\sum_{i=1}^{n}\|\boldsymbol{y}_i\|^2$。 **1.** 内積の線形性から $\frac{1}{n}\sum_i z_{ik} = \left\langle \boldsymbol{v}_k, \frac{1}{n}\sum_i \boldsymbol{y}_i \right\rangle = \langle \boldsymbol{v}_k, \boldsymbol{0}\rangle = 0$。ここでの中心化から $\sum_i \boldsymbol{y}_i = \boldsymbol{0}$ を使いました。 **2.** $z_{ik} z_{il} = (\boldsymbol{v}_k^{\mathsf{T}}\boldsymbol{y}_i)(\boldsymbol{y}_i^{\mathsf{T}}\boldsymbol{v}_l)$ と書けるので $$ \frac{1}{n}\sum_i z_{ik} z_{il} = \boldsymbol{v}_k^{\mathsf{T}}\left(\frac{1}{n}\sum_i \boldsymbol{y}_i\boldsymbol{y}_i^{\mathsf{T}}\right)\boldsymbol{v}_l = \boldsymbol{v}_k^{\mathsf{T}} S \boldsymbol{v}_l = \lambda_l \boldsymbol{v}_k^{\mathsf{T}}\boldsymbol{v}_l = \lambda_l \delta_{kl}, $$ 最後は正規直交性です。$k = l$ なら $\lambda_k$、$k \ne l$ なら $0$。1. より得点の平均が $0$ なので、この量はそのまま標本共分散です。 **3.** まず $\operatorname{tr}(\boldsymbol{y}\boldsymbol{y}^{\mathsf{T}}) = \sum_j y_j^2 = \|\boldsymbol{y}\|^2$ とトレースの線形性から $$ \operatorname{tr} S = \operatorname{tr}\left(\frac{1}{n}\sum_i \boldsymbol{y}_i\boldsymbol{y}_i^{\mathsf{T}}\right) = \frac{1}{n}\sum_i \|\boldsymbol{y}_i\|^2 . $$ また $\operatorname{tr}S = \sum_j s_{jj}$ はトレースの定義そのものです。最後にの $S = V\Lambda V^{\mathsf{T}}$ と $\operatorname{tr}(AB) = \operatorname{tr}(BA)$、$V^{\mathsf{T}}V = I$ から $$ \operatorname{tr}S = \operatorname{tr}(V\Lambda V^{\mathsf{T}}) = \operatorname{tr}(\Lambda V^{\mathsf{T}}V) = \operatorname{tr}\Lambda = \sum_k \lambda_k . $$ 2. は PCA のもう一つの顔です。もとの変数は互いに相関していますが、**主成分得点に移ると相関が消えます**。$V^{\mathsf{T}}$ による座標変換は、共分散行列を対角化する変換にほかなりません（[対角化とジョルダン標準形](/mathematics/linear-algebra/diagonalization-and-jordan-form)、）。3. は「全分散はどの座標系で測っても同じで、それが固有値の総和に等しい」という保存則で、次節の寄与率の分母になります。 ## 5. もう一つの顔：当てはめ誤差の最小化 Pearson の出発点は分散ではなく距離でした。ここでは、その定式化が同じ固有ベクトルに行き着くことを示します。$V \subseteq \mathbb{R}^p$ を線形部分空間、$P_V$ をそれへの直交射影とします。 $1 \le k \le p$ とする。$k$ 次元線形部分空間 $V \subseteq \mathbb{R}^p$ と点 $\boldsymbol{c} \in \mathbb{R}^p$ が定めるアフィン部分空間 $A = \boldsymbol{c} + V$ に対し、当てはめ誤差を $$ E(\boldsymbol{c}, V) = \frac{1}{n}\sum_{i=1}^{n} \operatorname{dist}(\boldsymbol{x}_i, A)^2 = \frac{1}{n}\sum_{i=1}^{n} \left\| (I - P_V)(\boldsymbol{x}_i - \boldsymbol{c}) \right\|^2 $$ と定める。このとき $$ \min_{\boldsymbol{c} \in \mathbb{R}^p,\ \dim V = k} E(\boldsymbol{c}, V) = \lambda_{k+1} + \lambda_{k+2} + \cdots + \lambda_p $$ であり、最小値は $\boldsymbol{c} = \bar{\boldsymbol{x}}$、$V = \operatorname{span}\{\boldsymbol{v}_1, \ldots, \boldsymbol{v}_k\}$ で達成される（$k = p$ のとき右辺は $0$ と読む）。 **ステップ 0（距離の式）.** $A = \boldsymbol{c} + V$ 上の点は $\boldsymbol{c} + \boldsymbol{w}$（$\boldsymbol{w} \in V$）と書けます。$\|\boldsymbol{x} - \boldsymbol{c} - \boldsymbol{w}\|$ を $\boldsymbol{w} \in V$ について最小にするのは直交射影 $\boldsymbol{w} = P_V(\boldsymbol{x} - \boldsymbol{c})$ であり（）、そのときの残差が $(I - P_V)(\boldsymbol{x}-\boldsymbol{c})$ です。よって定義式の二つの表現は一致します。 **ステップ 1（重心が最適であること）.** $\boldsymbol{x}_i - \boldsymbol{c} = \boldsymbol{y}_i + (\bar{\boldsymbol{x}} - \boldsymbol{c})$ と分解し、$\boldsymbol{d} = (I-P_V)(\bar{\boldsymbol{x}} - \boldsymbol{c})$ とおくと $(I-P_V)(\boldsymbol{x}_i - \boldsymbol{c}) = (I-P_V)\boldsymbol{y}_i + \boldsymbol{d}$。二乗して和を取ると $$ n E(\boldsymbol{c}, V) = \sum_i \|(I-P_V)\boldsymbol{y}_i\|^2 + 2\left\langle \sum_i (I-P_V)\boldsymbol{y}_i,\ \boldsymbol{d} \right\rangle + n\|\boldsymbol{d}\|^2 . $$ 中央の項は、$\sum_i \boldsymbol{y}_i = \boldsymbol{0}$ と $I - P_V$ の線形性から $\sum_i (I-P_V)\boldsymbol{y}_i = (I-P_V)\boldsymbol{0} = \boldsymbol{0}$ となって消えます。残りは $\|\boldsymbol{d}\|^2 \ge 0$ なので、$\boldsymbol{c} = \bar{\boldsymbol{x}}$（このとき $\boldsymbol{d} = \boldsymbol{0}$）が最適です。以下 $\boldsymbol{c} = \bar{\boldsymbol{x}}$ とし、$E(V) = \frac{1}{n}\sum_i \|(I-P_V)\boldsymbol{y}_i\|^2$ と書きます。**データを中心化してよい理由がここにあります。** **ステップ 2（ピタゴラス：誤差最小化＝分散最大化）.** 直交射影の性質から $P_V \boldsymbol{y} \perp (I - P_V)\boldsymbol{y}$ なので、 $$ \|\boldsymbol{y}_i\|^2 = \|P_V \boldsymbol{y}_i\|^2 + \|(I-P_V)\boldsymbol{y}_i\|^2 . $$ $i$ について平均を取りの 3. を使うと $$ E(V) = \operatorname{tr}S - \frac{1}{n}\sum_i \|P_V\boldsymbol{y}_i\|^2 . $$ $\operatorname{tr}S$ は $V$ に依らない定数なので、**$E(V)$ の最小化は射影後の分散の和 $\frac{1}{n}\sum_i\|P_V\boldsymbol{y}_i\|^2$ の最大化と同値**です。Pearson の問題と Hotelling の問題が一致するのは、このピタゴラスの定理 1 本によります。 **ステップ 3（最大化）.** $\boldsymbol{u}_1, \ldots, \boldsymbol{u}_k$ を $V$ の正規直交基底とすると $P_V\boldsymbol{y} = \sum_{m=1}^k \langle \boldsymbol{u}_m, \boldsymbol{y}\rangle \boldsymbol{u}_m$、したがって $\|P_V\boldsymbol{y}\|^2 = \sum_{m=1}^k \langle \boldsymbol{u}_m,\boldsymbol{y}\rangle^2$ です。§2 の射影分散の計算を各 $\boldsymbol{u}_m$ に適用すると $$ \frac{1}{n}\sum_i \|P_V\boldsymbol{y}_i\|^2 = \sum_{m=1}^{k} \boldsymbol{u}_m^{\mathsf{T}} S \boldsymbol{u}_m . $$ これを $\lambda_1 + \cdots + \lambda_k$ で上から抑えます。$w_j = \sum_{m=1}^{k}\langle \boldsymbol{v}_j, \boldsymbol{u}_m\rangle^2$ とおくと、の証明ステップ 2 と同じ展開から $\sum_m \boldsymbol{u}_m^{\mathsf{T}}S\boldsymbol{u}_m = \sum_{j=1}^{p}\lambda_j w_j$ です。重み $w_j$ は二つの制約を満たします。 - $0 \le w_j \le 1$：$w_j = \|P_V \boldsymbol{v}_j\|^2$ であり、$\|P_V\boldsymbol{v}_j\| \le \|\boldsymbol{v}_j\| = 1$（、あるいはステップ 2 のピタゴラス）。 - $\sum_{j=1}^p w_j = k$：和の順序を変えると $\sum_j w_j = \sum_{m=1}^{k}\sum_{j=1}^{p}\langle \boldsymbol{v}_j,\boldsymbol{u}_m\rangle^2 = \sum_{m=1}^k \|\boldsymbol{u}_m\|^2 = k$（$\boldsymbol{v}_j$ が正規直交基底であることによるパーセバルの等式）。このとき $$ \sum_{j=1}^{p}\lambda_j w_j - \sum_{j=1}^{k}\lambda_j = \sum_{j\le k}\lambda_j (w_j - 1) + \sum_{j > k}\lambda_j w_j \le \sum_{j\le k}\lambda_k (w_j - 1) + \sum_{j>k}\lambda_k w_j = \lambda_k\left(\sum_{j=1}^p w_j - k\right) = 0 . $$ 不等号の根拠は次のとおりです。$j \le k$ では $w_j - 1 \le 0$ かつ $\lambda_j \ge \lambda_k$ なので $\lambda_j(w_j-1) \le \lambda_k(w_j-1)$。$j > k$ では $w_j \ge 0$ かつ $\lambda_j \le \lambda_k$ なので $\lambda_j w_j \le \lambda_k w_j$。したがって $\sum_m \boldsymbol{u}_m^{\mathsf{T}}S\boldsymbol{u}_m \le \lambda_1 + \cdots + \lambda_k$ です。 **ステップ 4（達成と結論）.** $V = \operatorname{span}\{\boldsymbol{v}_1,\ldots,\boldsymbol{v}_k\}$、$\boldsymbol{u}_m = \boldsymbol{v}_m$ と取れば $\sum_m \boldsymbol{v}_m^{\mathsf{T}}S\boldsymbol{v}_m = \sum_{m\le k}\lambda_m$ で等号が成立します。ステップ 2 に戻して $$ \min_V E(V) = \operatorname{tr}S - \sum_{j=1}^{k}\lambda_j = \sum_{j=1}^{p}\lambda_j - \sum_{j=1}^{k}\lambda_j = \sum_{j=k+1}^{p}\lambda_j $$ を得ます（の 3. を使いました）。この定理と[線形回帰](/computer-science/math-for-ml/linear-regression)の違いに注意してください。回帰は目的変数の方向（縦方向）に測った残差の二乗和を最小にします（）。PCA は直線・平面への**垂直距離**を最小にします。したがって $y$ を $x$ で回帰した直線と、$x$ を $y$ で回帰した直線と、PCA の第 1 主成分軸は、一般に三本とも別の直線です。どれが正しいかは「何を誤差とみなすか」という問題設定で決まります。説明変数にも測定誤差がある場合に垂直距離を採る立場は、全最小二乗法（total least squares）と呼ばれます。 ## 6. 寄与率：どこまで採ればよいかは、$k$ 次元まで採ったときに捨てる誤差が $\lambda_{k+1} + \cdots + \lambda_p$ ちょうどであることを教えてくれます。全分散が $\sum_j \lambda_j$（）ですから、比を取れば「何割を残したか」が固有値だけで読めます。 $\operatorname{tr}S = \sum_{j=1}^p \lambda_j > 0$ とする。第 $k$ 主成分の**寄与率**を $$ r_k = \frac{\lambda_k}{\lambda_1 + \cdots + \lambda_p} $$ と定め、上位 $k$ 本までの**累積寄与率**を $$ R_k = \sum_{m=1}^{k} r_m = \frac{\lambda_1 + \cdots + \lambda_k}{\lambda_1 + \cdots + \lambda_p} $$ と定める。 $0 \le r_k \le 1$ かつ $\sum_k r_k = 1$ は定義と $\lambda_k \ge 0$ から直ちに従います。意味は次の等式に集約されます。の最小誤差を $E_k^{*}$ と書けば $$ \frac{E_k^{*}}{\operatorname{tr}S} = \frac{\lambda_{k+1}+\cdots+\lambda_p}{\lambda_1+\cdots+\lambda_p} = 1 - R_k . $$ つまり**累積寄与率 $R_k$ は、$k$ 次元に落としたときに保たれる分散の割合であり、同時に $1 - R_k$ が相対的な平均二乗再構成誤差**です。「累積寄与率 90% で打ち切る」という操作は、「もとの散らばりの 10% ぶんの二乗誤差を許容して圧縮する」と言い換えられます。 $\mathbb{R}^2$ の 5 点 $(1,2), (2,3), (3,5), (4,4), (5,6)$ を PCA にかけます。 **重心と中心化.** $\bar{\boldsymbol{x}} = \frac{1}{5}(1+2+3+4+5,\ 2+3+5+4+6) = (3, 4)$。よって中心化ベクトルは $$ \boldsymbol{y}_1 = (-2,-2),\quad \boldsymbol{y}_2 = (-1,-1),\quad \boldsymbol{y}_3 = (0,1),\quad \boldsymbol{y}_4 = (1,0),\quad \boldsymbol{y}_5 = (2,2). $$ 検算として $\sum_i \boldsymbol{y}_i = (-2-1+0+1+2,\ -2-1+1+0+2) = (0,0)$。 **共分散行列.** $\sum_i y_{i1}^2 = 4+1+0+1+4 = 10$、$\sum_i y_{i2}^2 = 4+1+1+0+4 = 10$、$\sum_i y_{i1}y_{i2} = 4+1+0+0+4 = 9$ なので $$ S = \frac{1}{5}\begin{pmatrix} 10 & 9 \\ 9 & 10 \end{pmatrix} = \begin{pmatrix} 2 & 1.8 \\ 1.8 & 2 \end{pmatrix}. $$ **固有値と固有ベクトル.** $\det(S - \lambda I) = (2-\lambda)^2 - 1.8^2 = 0$ より $2 - \lambda = \pm 1.8$、すなわち $\lambda_1 = 3.8$、$\lambda_2 = 0.2$。$\lambda_1$ に対しては $(S - 3.8 I)\boldsymbol{u} = \begin{pmatrix} -1.8 & 1.8 \\ 1.8 & -1.8\end{pmatrix}\boldsymbol{u} = \boldsymbol{0}$ から $u_1 = u_2$、よって $\boldsymbol{v}_1 = \frac{1}{\sqrt{2}}(1,1)^{\mathsf{T}}$。同様に $\lambda_2$ に対しては $u_1 = -u_2$ から $\boldsymbol{v}_2 = \frac{1}{\sqrt{2}}(1,-1)^{\mathsf{T}}$。内積は $\frac{1}{2}(1 - 1) = 0$ で直交しており、と整合します。 **主成分得点.** $z_{i1} = \frac{y_{i1}+y_{i2}}{\sqrt{2}}$、$z_{i2} = \frac{y_{i1}-y_{i2}}{\sqrt{2}}$ なので $$ (z_{i1})_{i=1}^{5} = \tfrac{1}{\sqrt2}(-4,\, -2,\, 1,\, 1,\, 4), \qquad (z_{i2})_{i=1}^{5} = \tfrac{1}{\sqrt2}(0,\, 0,\, -1,\, 1,\, 0). $$ 分散を確かめます。$\frac{1}{5}\sum_i z_{i1}^2 = \frac{1}{5}\cdot\frac{16+4+1+1+16}{2} = \frac{38}{10} = 3.8 = \lambda_1$、$\frac{1}{5}\sum_i z_{i2}^2 = \frac{1}{5}\cdot\frac{0+0+1+1+0}{2} = \frac{2}{10} = 0.2 = \lambda_2$。共分散は $\frac{1}{5}\cdot\frac{1}{2}\big(0 + 0 + (1)(-1) + (1)(1) + 0\big) = 0$。の 2. のとおりです。 **寄与率.** $\operatorname{tr}S = 4 = 3.8 + 0.2$ なので $r_1 = 3.8/4 = 0.95$、$r_2 = 0.05$。第 1 主成分だけで全分散の 95% を説明します。 **再構成誤差.** $k = 1$ に落とすと、$\boldsymbol{y}_i$ の残差は $z_{i2}\boldsymbol{v}_2$ で、その二乗ノルムは $z_{i2}^2 = 0, 0, 0.5, 0.5, 0$。平均は $\frac{1}{5}(0+0+0.5+0.5+0) = 0.2 = \lambda_2$ で、と一致します。復元値は $\hat{\boldsymbol{x}}_i = \bar{\boldsymbol{x}} + z_{i1}\boldsymbol{v}_1$、たとえば $i=3$ なら $z_{31} = 1/\sqrt2$ より $\hat{\boldsymbol{x}}_3 = (3,4) + \frac{1}{\sqrt2}\cdot\frac{1}{\sqrt2}(1,1) = (3.5, 4.5)$ です。

すべての点が 1 本の直線上にある場合、すなわち単位ベクトル $\boldsymbol{a}$ と実数 $t_1, \ldots, t_n$（$\sum_i t_i = 0$ としてよい）によって $\boldsymbol{y}_i = t_i \boldsymbol{a}$ と書ける場合を考えます。このとき $$ S = \frac{1}{n}\sum_i (t_i\boldsymbol{a})(t_i\boldsymbol{a})^{\mathsf{T}} = \left(\frac{1}{n}\sum_i t_i^2\right)\boldsymbol{a}\boldsymbol{a}^{\mathsf{T}} = s_t^2\, \boldsymbol{a}\boldsymbol{a}^{\mathsf{T}} $$ です。この行列は $\boldsymbol{a}$ を $s_t^2\boldsymbol{a}$ に写し（$\boldsymbol{a}^{\mathsf{T}}\boldsymbol{a} = 1$ より）、$\boldsymbol{a}$ に直交するベクトル $\boldsymbol{w}$ を $s_t^2 \boldsymbol{a}(\boldsymbol{a}^{\mathsf{T}}\boldsymbol{w}) = \boldsymbol{0}$ に写します。したがって固有値は $\lambda_1 = s_t^2$（固有ベクトル $\boldsymbol{a}$）と $0$（重複度 $p-1$）。寄与率は $r_1 = 1$、$k=1$ での再構成誤差は $\lambda_2 + \cdots + \lambda_p = 0$。PCA は、データが実際に乗っている直線をそのまま見つけ出します。逆に言えば、**固有値 $0$ はデータに厳密な線形従属関係があることの検出器**です。 $S\boldsymbol{u} = \frac{1}{n}\sum_i \boldsymbol{y}_i \langle \boldsymbol{y}_i, \boldsymbol{u}\rangle$ は常に $\operatorname{span}\{\boldsymbol{y}_1,\ldots,\boldsymbol{y}_n\}$ に属するので $\operatorname{rank}S \le \dim\operatorname{span}\{\boldsymbol{y}_i\}$ です。さらに $\sum_i \boldsymbol{y}_i = \boldsymbol{0}$ という線形関係があるため、$n$ 本のベクトル $\boldsymbol{y}_i$ が張る空間の次元は高々 $n-1$。したがって $$ \operatorname{rank}S \le \min(p,\ n-1) $$ です。標本数が次元以下（$n \le p$）のデータでは必ず $\lambda_p = 0$ となり、共分散行列は正則になりません。遺伝子発現データのように $p \gg n$ の場面では、非零の固有値は高々 $n-1$ 個しかなく、PCA が返せる主成分の本数も $n-1$ 本までです。 ## 7. 使う前に確認すべきこと

B["重心を引いて中心化"] B --> C["標本共分散行列 S を作る"] C --> D["固有値分解: λ_1 ≥ ... ≥ λ_p ≥ 0"] D --> E["累積寄与率などで次元 k を決める"] E --> F["上位 k 本の固有ベクトルに射影"] F --> G["k 次元表現 z_i を得る"] D --> H["捨てる誤差 = λ_(k+1) + ... + λ_p"]`} /> **中心化を忘れない.** 平均を引かずに $\frac{1}{n}\sum_i \boldsymbol{x}_i\boldsymbol{x}_i^{\mathsf{T}}$ を固有値分解すると何が起きるでしょうか。$\boldsymbol{x}_i = \boldsymbol{y}_i + \bar{\boldsymbol{x}}$ を代入して展開すると、交差項が $\sum_i \boldsymbol{y}_i = \boldsymbol{0}$ で消えるので $$ \frac{1}{n}\sum_i \boldsymbol{x}_i\boldsymbol{x}_i^{\mathsf{T}} = S + \bar{\boldsymbol{x}}\,\bar{\boldsymbol{x}}^{\mathsf{T}} $$ となります。余計に足された $\bar{\boldsymbol{x}}\bar{\boldsymbol{x}}^{\mathsf{T}}$ は $\bar{\boldsymbol{x}}$ 方向に固有値 $\|\bar{\boldsymbol{x}}\|^2$ をもつ階数 1 の行列です。原点から遠いデータでは $\|\bar{\boldsymbol{x}}\|^2$ が $\lambda_1$ より大きくなりがちで、そのとき「第 1 主成分」は散らばりの方向ではなく重心の方向を指してしまいます。のステップ 1 が示すとおり、重心を通らせるのが当てはめの意味でも最適です。 **スケールに依存する.** 変数ごとに単位を取り替えると、共分散行列は $\tilde{S} = D S D$（$D$ は正の対角行列）に変わりますが、固有ベクトルは $D\boldsymbol{v}_j$ にはなりません。 $S = \begin{pmatrix} 4 & 2 \\ 2 & 4\end{pmatrix}$ とします。固有値は $4 \pm 2$、すなわち $\lambda_1 = 6$（固有ベクトル $\frac{1}{\sqrt2}(1,1)$）と $\lambda_2 = 2$（$\frac{1}{\sqrt2}(1,-1)$）で、第 1 主成分軸は $45^\circ$ の方向です。ここで第 2 変数の単位を $\frac{1}{10}$ に細かくする（たとえば m を dm にする）と、$D = \operatorname{diag}(1,10)$ として $$ \tilde{S} = D S D = \begin{pmatrix} 4 & 20 \\ 20 & 400 \end{pmatrix}. $$ $\operatorname{tr}\tilde S = 404$、$\det \tilde S = 1600 - 400 = 1200$ なので、固有値は $\lambda = \frac{404 \pm \sqrt{404^2 - 4\cdot 1200}}{2} = \frac{404\pm\sqrt{158416}}{2}$、数値では $\lambda_1 \approx 401.01$、$\lambda_2 \approx 2.99$。第 1 主成分方向は $(4 - 401.01)u_1 + 20u_2 = 0$ から $u_2 \approx 19.85\,u_1$、すなわち第 2 座標軸から約 $2.9^\circ$ しか離れていません。もとの $45^\circ$ とはまったく別の方向です。寄与率も $0.75$ から $0.993$ へ跳ね上がります。これは PCA の欠陥ではなく、分散という量が単位に依存する（次元をもつ）ことの当然の帰結です。身長 [cm] と年収 [円] のように単位が比較不能な変数が混ざるときは、各変数を標準偏差で割ってから PCA にかけます。これは相関行列の固有値問題を解くことと同じで、相関行列 PCA と呼ばれます。 **教師なしであることを忘れない.** PCA はラベルを見ません。分散が大きい方向が、知りたい区別に対応しているとは限らないのです。 $\mathbb{R}^2$ に 2 つのクラスがあり、クラス A は $(0, 1)$ を中心に、クラス B は $(0,-1)$ を中心に、どちらも第 1 座標方向に標準偏差 $10$、第 2 座標方向に標準偏差 $0.1$ でばらついているとします。全データをまとめた共分散行列は、対称性から近似的に $$ S \approx \begin{pmatrix} 100 & 0 \\ 0 & 1 + 0.01 \end{pmatrix} $$ です（第 2 成分の分散は、クラス中心の差による $1$ とクラス内のばらつき $0.01$ の和）。よって $\lambda_1 = 100$、第 1 主成分方向は第 1 座標軸で、寄与率は $100/101.01 \approx 0.99$。ところが 2 つのクラスは第 1 座標では完全に重なっており、区別しているのは第 2 座標だけです。第 1 主成分だけを残すと、寄与率 99% を確保しながらクラスの情報を 100% 失います。ラベルが使える場合は、クラス間分散とクラス内分散の比を最大化する線形判別分析（LDA）のような教師あり手法が適切です。PCA を前処理に使うなら、「分散の大きさ」と「タスクにとっての有用性」は別物だと意識してください。 **線形であることを忘れない.** PCA が探すのは部分空間、すなわち平らな集合です。データが曲がった曲面（渦巻き状の 2 次元曲面など）に沿って分布している場合、それを平面で近似することになり、うまく展開できません。非線形に拡張する道具としてはカーネル PCA や、[ニューラルネットワーク](/computer-science/math-for-ml/backpropagation)を使ったオートエンコーダがあります。逆に、隠れ層の活性化が恒等写像で損失が二乗誤差の場合、オートエンコーダの大域最適解が張る部分空間は PCA のそれと一致することが知られています（Baldi–Hornik, 1989）。線形の世界では PCA が最適解であり、非線形性を入れて初めてその先へ進めるわけです。 **確率モデルとして見ることもできる.** ここまでは標本共分散行列という記述統計の量だけを扱いましたが、$\boldsymbol{x}$ が期待値 $\boldsymbol{\mu}$、共分散行列 $\Sigma$ の確率ベクトルであるとして、$\operatorname{Var}(\boldsymbol{u}^{\mathsf{T}}\boldsymbol{x}) = \boldsymbol{u}^{\mathsf{T}}\Sigma\boldsymbol{u}$ を最大化する問題としても同じ理論が成り立ちます（[確率変数と期待値](/mathematics/probability/random-variables)、）。この視点を進めると、潜在変数モデルとしての確率的 PCA が得られます（[確率論とベイズ統計の役割](/computer-science/math-for-ml/bayesian-statistics)）。 ## 8. 演習 2 変数データの標本共分散行列が $S = \begin{pmatrix} 5 & 2 \\ 2 & 2\end{pmatrix}$ であるとする。(1) 固有値を求めよ。(2) 正規化した主成分方向を求め、直交していることを確かめよ。(3) 第 1 主成分の寄与率を求めよ。(4) 第 1 主成分だけを残したときの平均二乗再構成誤差を求めよ。 **(1)** 特性多項式は $$ \det(S - \lambda I) = (5-\lambda)(2-\lambda) - 2\cdot 2 = \lambda^2 - 7\lambda + 6 = (\lambda - 6)(\lambda - 1) $$ なので $\lambda_1 = 6$、$\lambda_2 = 1$。検算：$\lambda_1 + \lambda_2 = 7 = \operatorname{tr}S$、$\lambda_1\lambda_2 = 6 = \det S$。 **(2)** $\lambda_1 = 6$ について $S - 6I = \begin{pmatrix} -1 & 2 \\ 2 & -4 \end{pmatrix}$。第 1 行から $-u_1 + 2u_2 = 0$、すなわち $u_1 = 2u_2$（第 2 行 $2u_1 - 4u_2 = 0$ は同じ条件です）。長さを $1$ にして $\boldsymbol{v}_1 = \frac{1}{\sqrt5}(2,1)^{\mathsf{T}}$。 $\lambda_2 = 1$ について $S - I = \begin{pmatrix} 4 & 2 \\ 2 & 1\end{pmatrix}$。第 1 行から $4u_1 + 2u_2 = 0$、すなわち $u_2 = -2u_1$。よって $\boldsymbol{v}_2 = \frac{1}{\sqrt5}(1,-2)^{\mathsf{T}}$。内積は $\frac{1}{5}(2\cdot 1 + 1\cdot(-2)) = 0$ で直交しており、と整合します。 **(3)** $\operatorname{tr}S = 7$ なので $r_1 = 6/7 \approx 0.857$。 **(4)** より、捨てられる誤差は残りの固有値の和 $\lambda_2 = 1$。全分散 $7$ の $1/7 \approx 14.3\%$ で、$1 - r_1$ に一致します。 2 変数を標準化した後の相関行列 $R = \begin{pmatrix} 1 & \rho \\ \rho & 1 \end{pmatrix}$（$-1 < \rho < 1$）を考える。固有値と正規化した固有ベクトルを $\rho$ の式で求め、第 1 主成分の寄与率を書け。さらに $\rho \to 0$ と $\rho \to 1$ の極限で何が起きるかを述べよ。 $\det(R - \lambda I) = (1-\lambda)^2 - \rho^2 = 0$ から $1 - \lambda = \pm\rho$、よって固有値は $1+\rho$ と $1-\rho$。 $\rho \ne 0$ とすると、固有値 $1+\rho$ については $R - (1+\rho)I = \begin{pmatrix} -\rho & \rho \\ \rho & -\rho\end{pmatrix}$ より $u_1 = u_2$、固有ベクトルは $\frac{1}{\sqrt2}(1,1)^{\mathsf{T}}$。固有値 $1-\rho$ については同様に $u_1 = -u_2$ で $\frac{1}{\sqrt2}(1,-1)^{\mathsf{T}}$。大小関係は $\rho$ の符号で決まります。$\rho > 0$ なら $\lambda_1 = 1+\rho$（方向 $\frac{1}{\sqrt2}(1,1)$）、$\rho<0$ なら $\lambda_1 = 1-\rho$（方向 $\frac{1}{\sqrt2}(1,-1)$）。いずれの場合も $\lambda_1 = 1 + |\rho|$、$\operatorname{tr}R = 2$ なので $$ r_1 = \frac{1 + |\rho|}{2}. $$ $\rho \to 0$：$r_1 \to 1/2$。実際 $\rho = 0$ では $R = I_2$ となりの状況で、主成分方向が定まりません。2 変数が無相関なら、1 次元に落とすと必ず分散の半分を失います。次元削減の利得がない、というのが正しい結論です。 $\rho \to 1$：$r_1 \to 1$。2 変数がほぼ同じ情報をもち、$\frac{1}{\sqrt2}(1,1)$ 方向の 1 次元でほぼ完全に表せます。捨てる誤差は $\lambda_2 = 1 - \rho \to 0$。なお、相関行列は各変数の分散を $1$ に揃えた後の共分散行列なので、$\operatorname{tr}R = p$（ここでは $2$）が常に成り立ちます。カイザー基準「$\lambda > 1$ の主成分だけ採る」は、「平均的な 1 変数分（$\operatorname{tr}R/p = 1$）より多く説明する主成分を残す」という意味です。 $Q \in \mathbb{R}^{p\times p}$ を直交行列（$Q^{\mathsf{T}}Q = I$）、$\boldsymbol{b} \in \mathbb{R}^p$ を定ベクトルとし、$\tilde{\boldsymbol{x}}_i = Q\boldsymbol{x}_i + \boldsymbol{b}$ とおく。(1) 新しいデータの標本共分散行列が $\tilde{S} = QSQ^{\mathsf{T}}$ となることを示せ。(2) $\tilde S$ の固有値は $S$ の固有値と一致し、対応する正規直交固有ベクトルが $Q\boldsymbol{v}_j$ で与えられることを示せ。(3) PCA の寄与率と再構成誤差が回転・鏡映・平行移動で不変であることを結論せよ。(4) 一方 $D = \operatorname{diag}(1,10)$、$\tilde{\boldsymbol{x}}_i = D\boldsymbol{x}_i$ の場合には、$S = \begin{pmatrix}4&2\\2&4\end{pmatrix}$ に対して $D\boldsymbol{v}_1$ が $\tilde S = DSD$ の固有ベクトルにならないことを確かめよ。 **(1)** 新しい重心は平均の線形性から $\frac{1}{n}\sum_i (Q\boldsymbol{x}_i + \boldsymbol{b}) = Q\bar{\boldsymbol{x}} + \boldsymbol{b}$。よって中心化ベクトルは $$ \tilde{\boldsymbol{y}}_i = (Q\boldsymbol{x}_i + \boldsymbol{b}) - (Q\bar{\boldsymbol{x}} + \boldsymbol{b}) = Q(\boldsymbol{x}_i - \bar{\boldsymbol{x}}) = Q\boldsymbol{y}_i $$ で、平行移動 $\boldsymbol{b}$ は消えます。したがって $$ \tilde S = \frac{1}{n}\sum_i (Q\boldsymbol{y}_i)(Q\boldsymbol{y}_i)^{\mathsf{T}} = Q\left(\frac{1}{n}\sum_i \boldsymbol{y}_i\boldsymbol{y}_i^{\mathsf{T}}\right)Q^{\mathsf{T}} = QSQ^{\mathsf{T}}. $$ **(2)** $S\boldsymbol{v}_j = \lambda_j\boldsymbol{v}_j$ の両辺に左から $Q$ を掛け、途中に $I = Q^{\mathsf{T}}Q$ を挿入します。 $$ \tilde S (Q\boldsymbol{v}_j) = QSQ^{\mathsf{T}}Q\boldsymbol{v}_j = QS\boldsymbol{v}_j = \lambda_j (Q\boldsymbol{v}_j). $$ また $\langle Q\boldsymbol{v}_j, Q\boldsymbol{v}_l\rangle = \boldsymbol{v}_j^{\mathsf{T}}Q^{\mathsf{T}}Q\boldsymbol{v}_l = \boldsymbol{v}_j^{\mathsf{T}}\boldsymbol{v}_l = \delta_{jl}$ なので $\{Q\boldsymbol{v}_j\}$ は正規直交系、とくに $p$ 個の一次独立なベクトルです。よってこれが $\tilde S$ の正規直交固有基底を与え、固有値の組（重複度込み）は $S$ のそれと一致します。 **(3)** 寄与率は固有値の比だけで決まるので不変です。再構成誤差もより $\lambda_{k+1}+\cdots+\lambda_p$ で、これも不変。主成分方向だけがデータと一緒に $Q$ で回ります。幾何的には、PCA が座標系ではなく点の配置そのものを見ていることを意味します。 **(4)** $S = \begin{pmatrix}4&2\\2&4\end{pmatrix}$ の固有ベクトルは $\boldsymbol{v}_1 = \frac{1}{\sqrt2}(1,1)^{\mathsf{T}}$（$\lambda_1 = 6$）です。$D\boldsymbol{v}_1$ は $(1,10)^{\mathsf{T}}$ の定数倍。一方 $\tilde S = DSD = \begin{pmatrix} 4 & 20 \\ 20 & 400\end{pmatrix}$ なので $$ \tilde S \begin{pmatrix} 1 \\ 10\end{pmatrix} = \begin{pmatrix} 4 + 200 \\ 20 + 4000 \end{pmatrix} = \begin{pmatrix} 204 \\ 4020 \end{pmatrix}. $$ もし $(1,10)^{\mathsf{T}}$ が固有ベクトルなら、像も第 2 成分が第 1 成分の $10$ 倍、すなわち $(204, 2040)^{\mathsf{T}}$ でなければなりません。$4020 \ne 2040$ なので固有ベクトルではありません。で見たとおり、実際の第 1 主成分方向は $u_2/u_1 \approx 19.85$ の方向です。結論として、PCA は直交変換（等長変換）に対しては同変ですが、一般の可逆線形変換に対してはそうではありません。だからこそ変数のスケーリングが本質的な前処理になります。 ## 参考文献 - K. Pearson, "On Lines and Planes of Closest Fit to Systems of Points in Space", *Philosophical Magazine* 2 (1901), 559–572. PCA を「点集合への最良当てはめ直線・平面」として導入した原論文。 - H. Hotelling, "Analysis of a Complex of Statistical Variables into Principal Components", *Journal of Educational Psychology* 24 (1933), 417–441, 498–520. 分散最大化による定式化と principal component の命名。 - I. T. Jolliffe, *Principal Component Analysis*, 2nd ed., Springer, 2002. PCA の標準的なモノグラフ。第 1〜3 章に本記事の内容と、寄与率・主成分の解釈に関する詳細な議論があります。 - C. M. Bishop, *Pattern Recognition and Machine Learning*, Springer, 2006 — 第 12 章「Continuous Latent Variables」。分散最大化と誤差最小化の両方の定式化、および確率的 PCA。 - G. H. Golub and C. F. Van Loan, *Matrix Computations*, 4th ed., Johns Hopkins University Press, 2013 — 特異値分解とその数値計算を扱う章。 - P. Baldi and K. Hornik, "Neural networks and principal component analysis: Learning from examples without local minima", *Neural Networks* 2 (1989), 53–58. 線形オートエンコーダの大域最適解が PCA の部分空間と一致することを示した論文。 ## Appendix: 特異値分解で計算する **共分散行列を作らずに済ませる.** 実務の PCA は、$S$ を作ってから固有値分解するのではなく、中心化データ行列 $X$ の特異値分解（SVD）を直接計算するのが標準です。理由は次の命題と、その後の数値的な注意にあります。 $X \in \mathbb{R}^{n\times p}$ をの中心化データ行列、$r = \operatorname{rank}X$ とし、その簡約特異値分解を $$ X = U \Sigma V^{\mathsf{T}}, \qquad U \in \mathbb{R}^{n\times r},\ \ \Sigma = \operatorname{diag}(\sigma_1 \ge \cdots \ge \sigma_r > 0),\ \ V \in \mathbb{R}^{p \times r} $$ （$U^{\mathsf{T}}U = V^{\mathsf{T}}V = I_r$）とする。このとき $V$ の第 $j$ 列は $S = \frac{1}{n}X^{\mathsf{T}}X$ の固有ベクトルであり、対応する固有値は $\lambda_j = \sigma_j^2/n$（$1 \le j \le r$）である。残る $p - r$ 個の固有値は $0$ である。さらに第 $j$ 主成分得点は $z_{ij} = \sigma_j u_{ij}$ で与えられる。 $U^{\mathsf{T}}U = I_r$ を使って $$ S = \frac{1}{n}X^{\mathsf{T}}X = \frac{1}{n}V\Sigma U^{\mathsf{T}} U \Sigma V^{\mathsf{T}} = \frac{1}{n} V\Sigma^2 V^{\mathsf{T}} = \sum_{j=1}^{r}\frac{\sigma_j^2}{n}\, \boldsymbol{v}_j\boldsymbol{v}_j^{\mathsf{T}} $$ （$\boldsymbol{v}_j$ は $V$ の第 $j$ 列）。両辺に $\boldsymbol{v}_l$ を掛けると、$V$ の列の正規直交性から $S\boldsymbol{v}_l = \frac{\sigma_l^2}{n}\boldsymbol{v}_l$ を得ます。一方 $\boldsymbol{w} \in \mathbb{R}^p$ が $V$ の全列と直交すれば $S\boldsymbol{w} = \boldsymbol{0}$ なので、$V$ の列直交補空間（次元 $p-r$）はすべて固有値 $0$ の固有空間に含まれます。以上で $p$ 個の固有値が出そろいました。得点については $z_{ij} = \langle \boldsymbol{v}_j, \boldsymbol{y}_i\rangle$ が行列 $XV$ の $(i,j)$ 成分であり、$XV = U\Sigma V^{\mathsf{T}}V = U\Sigma$ なので $z_{ij} = \sigma_j u_{ij}$ です。 **なぜ SVD のほうがよいのか.** $S = \frac{1}{n}X^{\mathsf{T}}X$ を明示的に作ると、$X$ の特異値が二乗されるため条件数も二乗されます。$\sigma_r/\sigma_1$ が計算機の丸め誤差の水準に近いとき、$X^{\mathsf{T}}X$ の段階で小さい固有値の情報が失われます。$X$ に直接 SVD をかければこの二乗を経由しません。計算量の面でも、$p \gg n$ なら $p\times p$ 行列を作らずに済みます。 **NumPy での確認.** の数値を両方の方法で再現します。 ```python import numpy as np X_raw = np.array([[1., 2.], [2., 3.], [3., 5.], [4., 4.], [5., 6.]]) X = X_raw - X_raw.mean(axis=0) # 中心化 n = X.shape[0] # 方法 1: 共分散行列の固有値分解（eigh は実対称行列用、固有値は昇順） S = X.T @ X / n lam, V = np.linalg.eigh(S) lam, V = lam[::-1], V[:, ::-1] # 降順に並べ替える # 方法 2: 中心化データ行列の特異値分解 U, sigma, Vt = np.linalg.svd(X, full_matrices=False) lam_svd = sigma ** 2 / n print(lam) # [3.8 0.2] print(lam_svd) # [3.8 0.2] print(lam / lam.sum()) # 寄与率 [0.95 0.05] print(V[:, 0]) # 第 1 主成分方向（符号は不定） print(X @ V[:, :1]) # 第 1 主成分得点 ``` `eigh` が固有値を昇順で返すので並べ替えが要ります。出力される固有ベクトルの符号は環境によって反転しうる点はのとおりで、得点の符号がまとめて変わるだけなので寄与率も再構成誤差も影響を受けません。