最近、特異値分解と主成分分析について復習する機会があったので、メモ代わりに記事にしておく。
特異値分解と主成分分析の概要
特異値分解（SVD）と主成分分析（PCA）は推薦システムではもちろんデータ分析でも重要な手法である。
特異値分解（SVD）
特異値分解 (Singular Value Decomposition, SVD) は行列 $A$ を以下のように分解する手法である。
$$ A = U \Sigma V^T $$
ここで $U$ は直交行列、$\Sigma$ は特異値の対角行列、$V$ は直交行列である。SVD は行列のランクを下げるために利用され、推薦システムでは評価行列の低ランク近似により評価予測を行う。
$$ \displaystyle R_k = U_k \Sigma_k V_k^T $$
主成分分析（PCA）
主成分分析 (Principal Component Analysis, PCA) はデータの次元削減手法である。共分散行列 $C$ の固有値分解を行い、固有ベクトルを利用して次元削減を行う。PCA も SVD を用いて実装可能である。
ソースコード
本記事で使用するソースコードは以下の通りである。
github
- jupyter notebook形式のファイルはこちら
google colaboratory
- google colaboratory で実行する場合はこちら
実行環境
OSはmacOSです。LinuxやUnixのコマンドとはオプションが異なりるので注意。
!sw_vers
ProductName: macOS
ProductVersion: 13.5.1
BuildVersion: 22G90
!python -V
Python 3.9.17
基本的なライブラリをインポートし watermark を利用してそのバージョンを確認しておきます。 ついでに乱数のseedの設定をします。
%matplotlib inline
%config InlineBackend.figure_format = 'svg'
import random
import numpy as np
from pprint import pprint
from watermark import watermark
seed = 123
random_state = 123
random.seed(seed)
np.random.seed(seed)
# 小数点を丸めたときに-0.0が出現するが、0.0 と -0.0 は等価であり、0.0として表示する関数
def ppprint(A):
"""
A: np.array
表示する行列
"""
pprint(np.where(A == -0.0, 0, A))
print(watermark(python=True, watermark=True, iversions=True, globals_=globals()))
Python implementation: CPython
Python version : 3.9.17
IPython version : 8.17.2
numpy: 1.25.2
Watermark: 2.4.3
Pythonによる実装と具体的な説明
特異値分解(SVD)
特異値分解は、任意の$m \times n$行列$A$を以下のように分解する方法である：
$$ A = U \Sigma V^T $$
ここで、$U$は$m \times m$の直交行列、$\Sigma$は$m \times n$の対角行列、$V$は$n \times n$の直交行列である。$\Sigma$の対角成分は$A$の特異値であり、非負の実数である。
SVDの計算例
例えば、$A$が$5 \times 2$の行列であるとする：
$$ A = \begin{pmatrix} 1 & 2 \\ 3 & 4 \\ 7 & 6 \\ 2 & 0 \\ 3 & 1 \\ \end{pmatrix} $$
この行列のSVDを求める。SVD自体はnumpyのlinalgモジュールで計算できる。以下に適当な行列$A$に対してSVDを行う実装を示す。
import numpy as np
A = np.array([[1, 2], [3, 4], [7, 6], [2, 0], [3, 1]])
U, S, V_T = np.linalg.svd(A, full_matrices=True)
print("=" * 40)
print("U : ")
pprint(U.round(2))
print("=" * 40)
print("S : ")
pprint(S.round(2))
# 特異値ベクトルSから特異値行列を作成
Sigma = np.zeros(A.shape)
for i in range(len(S)):
Sigma[i, i] = S[i]
print("=" * 40)
print("Sigma : ")
pprint(Sigma.round(2))
print("=" * 40)
print("V_T : ")
pprint(V_T.round(2))
========================================
U :
array([[-0.19, -0.37, -0.82, -0.21, -0.33],
[-0.44, -0.45, -0.05, 0.57, 0.52],
[-0.83, 0.06, 0.37, -0.25, -0.34],
[-0.13, 0.59, -0.27, 0.66, -0.36],
[-0.26, 0.55, -0.35, -0.35, 0.62]])
========================================
S :
array([11.13, 2.24])
========================================
Sigma :
array([[11.13, 0. ],
[ 0. , 2.24],
[ 0. , 0. ],
[ 0. , 0. ],
[ 0. , 0. ]])
========================================
V_T :
array([[-0.75, -0.66],
[ 0.66, -0.75]])
この結果、$U$、$S$、$V^T$が求まる。$S$は対角行列の対角成分を持つベクトルとして出力されるため、対角行列に変換する必要があることに注意する。
SVDの利用例：推薦システム
推薦システムにおいて、特異値分解を用いた方法は非常に有効である。例えば、ユーザーとアイテムの評価行列$R$をSVDにより分解することで、次元削減や潜在特徴の抽出が可能になる。これにより、ユーザーの評価傾向やアイテムの特徴を抽出し、精度の高い推薦が実現できる。
たとえば、特異値の大きさに基づいて特異値を降順に並べ、上位 $ k $ 個の特異値を選択する。 $ k $ は新しい次元数を表すパラメタである。選択した上位 $ k $ 個の特異値に対応する $ U $、$ \Sigma $、$ V^T $ の列ベクトルや行ベクトルを用いて、元の行列 $ A $ を近似する。具体的には、
$$ A_k = U[:, :k] \Sigma [:k, :k] V^T[:k, :] $$
となる。ここで、$ [:, :k] $ は行列の最初から $ k $ 列目までを抽出することを意味する。
この手順により、元の行列 $ A $ を特異値が大きい部分の情報のみを用いて近似することが可能である。これによって、次元削減やノイズの低減などの目的でデータをより効率的に表現することが可能になる。
簡単にメリットとデメリットを書いておく。
メリット
- 高次元データの次元削減が可能。
- 潜在特徴の抽出により、データの本質的な構造を把握できる。
デメリット
- 特に大規模データセットの場合計算コストが高い。
- データの欠損値に対する耐性が低い。
主成分分析(PCA)
主成分分析は、データセットの分散を最大化するようにデータを直交座標系に変換する方法である。
共分散行列の固有値分解を用いたPCA
PCAは、データの行列$\mathbf{X}$に対して、次のステップで実行される。
- データを示す行列$\mathbf{X}$を列方向に対して平均が0になるように変形し、共分散行列$\mathbf{C}$を計算する。
$\mathbf{X}_i$は$\mathbf{X}$の$i$列目を表す。
$$ C_{i j}=\mathrm{E}\left[\left(\mathbf{X}_i - \mu_i \cdot \mathbf{1} \right)\left(\mathbf{X}_j-\mu_j \cdot \mathbf{1} \right)\right]=\mathrm{E}\left(\mathbf{X}_i \mathbf{X}_j\right)-\mathrm{E}\left(\mathbf{X}_i\right) \mathrm{E}\left(\mathbf{X}_j\right) $$
ここで$\mathbf{1}$は要素が1のベクトルである。
$$ \mu_i=\mathrm{E}\left(\mathbf{X}_i\right) $$
$\mu_i$は$\mathbf{X}_i$の平均値である。
一般化すると、共分散行列は以下のように表される。
$$ \mathbf{C}=\mathrm{E}\left[(\mathbf{X}-\mathrm{E}[\mathbf{X}])^T(\mathbf{X}-\mathrm{E}[\mathbf{X}])\right] $$
共分散行列の固有値$\lambda_i$と固有ベクトル$\mathbf{v}_i$を計算する。
固有値が大きい順に固有ベクトルを選び、新しい直交基底とする。
レイリー商を用いたPCA
上記の手法だと、分散を最大化するように直交座標系に変換するというイメージと結びつかないので、レイリー商を利用して、分散を最大化するように直交座標系に変換することを考える。
具体的には、以下のような$\mathbf{X}^T\mathbf{X}$を対象としたレイリー商を最大化するようなベクトル$\mathbf{v}$を求める。
$$ \mathbf{v}_1=\underset{\mathbf{v} \neq \mathbf{0}}{\arg \max } \frac{||\mathbf{X} \mathbf{v}||^2}{||\mathbf{v}||^2} = \underset{\mathbf{v} \neq \mathbf{0}}{\arg \max } \frac{\mathbf{v}^T \mathbf{X}^T \mathbf{X} \mathbf{v}}{\mathbf{v}^T \mathbf{v}} $$
このように計算された$\mathbf{v}_1$は、レイリー商の特性により$\mathbf{X}^T\mathbf{X}$の最大固有値に対応する固有ベクトルとなる。 このことは以下のように確認することができる。
$$ \left(\mathbf{X} \mathbf{v}_1\right)^{T} \cdot\left(\mathbf{X} \mathbf{v}_1\right) =\mathbf{v}_1^{T} \mathbf{X}^{T} \mathbf{X} \mathbf{v}_1 = \mathbf{v}_1^{T} \lambda_1 \mathbf{v}_1 =\lambda_1||\mathbf{v}_1||^2 $$
上記の式は、$\mathbf{X}$を$\mathbf{v}_1$に射影した後の分散を表している。この値は$\mathbf{X}$の分散を最大化するような$\mathbf{v}_1$を求めることに対応している。
次に以下のように、$\mathbf{X}$から$\mathbf{v}_1$に射影したデータを引いた$\mathbf{X}_2$を計算する。
$$ \mathbf{X}_2=\mathbf{X}-\mathbf{X} \mathbf{v}_1 \mathbf{v}_1^{T} $$
次に、$\mathbf{v}_2$を求めるために、$\mathbf{v}_1$に直交するようなベクトル$\mathbf{v}_2$を求める。これは以下のように計算される。
$$ \mathbf{v}_2=\underset{\mathbf{v} \neq \mathbf{0}}{\arg \max } \frac{||\mathbf{X}_2 \mathbf{v}||^2}{||\mathbf{v}||^2} = \underset{\mathbf{v} \neq \mathbf{0}}{\arg \max } \frac{\mathbf{v}^T \mathbf{X}_2^T \mathbf{X}_2 \mathbf{v}}{\mathbf{v}^T \mathbf{v}} $$
このように計算された$\mathbf{v}_2$は、$\mathbf{X}^T\mathbf{X}$の2番目に大きい固有値に対応する固有ベクトルとなる。
このようにして、上位$k$個の固有ベクトルを求めることで、データの分散を最大化するような直交座標系に変換することができる。
PCAの計算例
次に、PythonでPCAを実装する。 適当なデータセットを用意し、PCAを実行する。主成分数は2とする。
最初にscikit-learnのPCAを使って、手書き数字データセットをPCAで次元削減する例を示す。
from sklearn.decomposition import PCA
# サンプルデータ行列
X = np.array([[1, 2, 1, 4, 2], [2, 5, 4, 2, 1], [2, 1, -1, -2, 3], [4, 8, 1, 2, -1]])
# 主成分数
k = 2
pca = PCA(n_components=k)
result = pca.fit_transform(X)
print("PCA後のデータ:")
pprint(result.round(2))
PCA後のデータ:
array([[ 1.26, 3.03],
[-1.85, 1.19],
[ 4.98, -2.15],
[-4.38, -2.08]])
次にsikit-learnのPCAを使わずに、上記の手法でPCAを実装する。
# サンプルデータ行列
X = np.array([[1, 2, 1, 4, 2], [2, 5, 4, 2, 1], [2, 1, -1, -2, 3], [4, 8, 1, 2, -1]])
# 観測値から平均を引いたデータを作成
C = X - np.mean(X, axis=0)
# 共分散行列の計算
cov = np.cov(C.T)
# 共分散行列の固有値と固有ベクトルの計算
vals, vecs = np.linalg.eig(cov)
# 固有値を降順にソートし、そのインデックスを取得
idx = np.argsort(vals)[::-1]
sorted_vals = vals[idx]
# 上位2つの固有値に対応する固有ベクトルを取得
top_vecs = vecs[:, idx[:2]]
# データを主成分に射影
PC = C.dot(top_vecs)
print("PCA後のデータ:")
pprint(PC.real.round(2))
PCA後のデータ:
array([[ 1.26, -3.03],
[-1.85, -1.19],
[ 4.98, 2.15],
[-4.38, 2.08]])
上記の二つの結果は第2列の符号を除き一致する事が確認できた。
PCAの利用例：推薦システム
PCAはデータの次元削減に利用され、推薦システムにおいても有効である。例えば、ユーザー評価行列の次元を削減し、計算コストを削減しつつも高い精度の推薦を実現できる。
メリット
- データの次元削減により、計算コストを削減できる。
- データの分散を最大化するため、情報の損失が少ない。
デメリット
- 線形変換に基づくため、非線形な関係を捉えるのが難しい。
- 次元削減の過程で重要な情報が失われる可能性がある。
結論
本記事では、特異値分解と主成分分析についてPythonの実装を交えつつ解説した。
特異値分解と主成分分析は、データの次元削減や特徴抽出において非常に有効な技術である。特異値分解は推薦システムにおいて、ユーザーとアイテムの潜在特徴を抽出し、高精度な推薦を実現する。一方、主成分分析はデータの分散を最大化することで、重要な特徴を抽出し、計算コストを削減する。両者ともにメリットとデメリットが存在するため、利用目的に応じて適切な手法を選択することが重要である。