はじめに
前回の記事(ノルム最小解の射影による導出) で、連立線形方程式の
- 最小二乗解
- ノルム最小解
についてまとめました。ここでは機械学習の分野でこれらが登場する場面として、線形回帰問題、およびその双対系であるカーネル回帰をとりあげます。
具体的には PRML の §3.1, §6.1 で登場する話題です。
線形回帰問題の定式化
基本的にPRMLと同じ表記を使います。
学習データとして、$N$ 個の入力データ $\{x_i\}$、出力データ $\{t_i\}$ の組がある。
入力データ $x$ をある関数セット $\phi_j\ (j=1,\dots,M)$ に通し、これを並べたもの
$$\phi(x) = \begin{pmatrix}
\phi_1(x) \\
\vdots \\
\phi_M(x)
\end{pmatrix} $$ を特徴量ベクトル1と呼ぶことにします。関数セットのことは基底関数と呼びます。
線形回帰モデルでは、特徴量の線形和を出力とします。つまり、
$$\sum_{j=1}^M \phi_j(x_i) w_j \simeq t_i~~~(i=1, \cdots , N)$$ であるような係数 $w_j$ を学習します。係数をベクトルにしてまとめて書くと、
$$\Phi w \simeq t$$ と書けます。この行列 $\Phi$ は計画行列と呼ばれます。重要なので、きちんと書き下しておきます。
$$\Phi := \begin{pmatrix}
\phi(x_1)^T \\
\vdots \\
\phi(x_N)^T \\
\end{pmatrix}
=\begin{pmatrix}
\phi_1(x_1) &\cdots & \phi_M(x_1) \\
\vdots & & \vdots \\
\phi_1(x_N) &\cdots & \phi_M(x_N) \\
\end{pmatrix}$$ 計画行列は「データ数(N)×特徴量の数(M)」の行列であることをしっかり覚えておくことが重要です。
結局、線形回帰問題とは、$w$ についての「解けない」連立線形方程式
$$\Phi w = t$$ について(何らかの意味で)最良な解を求める問題であることが分かります。
以前の記事でまとめた最小二乗解、ノルム最小解(もしくは一般逆行列)がここで登場します。
最小二乗回帰
通常の回帰問題では、データ数 $N$ が多いため、特徴量ベクトル $\phi$ の次元 $M$ よりも大きいです。
このような場合は、最小二乗解を使えば二乗誤差を最小化する近似となります。
$$w = (\Phi^T\Phi)^{-1}\Phi^T t$$
最小二乗回帰では、右辺に出てくる $\Phi^T\Phi $ について解くことが計算の中心になります。
$$ \Phi^T\Phi = \sum_{i=1}^N \phi(x_i)\phi(x_i)^T $$ これは次数 $M$ の正方行列です。$M$ が小さい、すなわち、基底関数の数が少ないほど、少ない計算量で解けることが分かります。したがって、問題に応じて、そのような上手い基底関数を見つけることが重要になってきます。
多項式関数 $\phi_j(x) = x^j\ (j=0,1,\dots)$ を基底関数として最小二乗回帰を行う Python のコードは以下のようになります2。
# 多項式最小2乗回帰
class LeastSquare:
# N: 多項式次数
def __init__(self, N):
self.N = N
# データを多項式基底関数による表示に直す
def __data_vec(self, x):
return np.array([x**n for n in range(self.N + 1)])
# 計画行列 X を返す
def __design_mat(self, xv):
return np.array([self.__data_vec(x) for x in xv])
# 学習データを入力して係数wを決定する
def fit(self, x, y):
X = self.__design_mat(x)
# solve (X^t X) w = X^t y
self.w = np.linalg.solve(np.dot(X.T, X), np.dot(X.T, y))
# 予測値を返す
def predict(self, x):
# 配列かどうかで分けている(ダサい)
try:
X = self.__design_mat(x)
except TypeError:
X = self.__data_vec(x)
return np.dot(X, self.w)
入力データとして正弦波にガウス乱数を足したものを与えて、5次多項式で回帰した例を示します。青線がデータ生成に用いた正弦波、水色の点がそれに乱数を足した学習データ、赤線が回帰曲線です。
カーネル回帰
カーネル回帰では、特徴量ベクトル $\phi(x)$ の次元 $M$ がデータ数 $N$ よりも大きいことを仮定します。
$M>N$ のときは、方程式 $\Phi w = t$ は不定で、ノルム最小解
$$ w = \Phi^T(\Phi\Phi^T)^{-1} t $$ が使えます。ただし、ノルム最小解は二乗誤差を0とするため、回帰問題にそのまま適用すると過学習になってしまいます。そこで、(L2)正則化項を入れて、
$$ w = \Phi^T(\Phi\Phi^T + \lambda I)^{-1} t $$ を使います。こうすることで、誤差を一定以下にしつつ、$||w||^2$ を最小化する解が得られます3。
カーネル回帰では右辺に現れる行列 $\Phi\Phi^T$ について解くことが計算の中心になります。この行列はグラム行列と呼ばれ、その行列要素は特徴量ベクトルの内積です。
$$ K_{ij} := (\Phi\Phi^T)_{ij} = \phi(x_i)^T\phi(x_j) $$ これは次数 $N$ の正方行列です。$N$ が小さい、すなわち、データ数が少ないほど、少ない計算量で解けます、したがって、データ数が少ない場合にカーネル回帰は適用しやすくなります。
カーネルトリック
この記事の主題からは逸れるのですが、カーネル回帰では基底関数 $\phi(x)$ を明示的に作らずにすべての計算を済ませることができます。その概要は以下の通りです。
カーネル回帰では基底関数 $\phi(x)$ を与える代わりに、その内積 $\phi(x_i)^T\phi(x_j)$ を関数として与えます。
$$ k(x, x') := \phi(x)^T\phi(x') $$
以下に示すように、(驚くべきことに)基底関数そのものではなく、その内積さえ与えれば、カーネル回帰の解は決まってしまいます。
まず、グラム行列は、カーネル関数を使って
$$ K_{ij} = k(x_i, x_j) $$ と求められます。次に、$ w = \Phi^T a$ とおいて、$w$ の代わりに $a$ を求めるようにすると、
$$ a = (\Phi\Phi^T + \lambda I)^{-1} t = (K + \lambda I)^{-1}t $$ となります。$a$ を使うと、$x$ における回帰予測値は、
$$ \phi(x)^T w = \phi(x)^T \Phi^T a = \sum_{j=1}^N \phi(x)^T \phi(x_j)a_j = \sum_{j=1}^N k(x, x_j) a_j $$ となるので、$\phi(x)$ を使わずにカーネル関数だけで得られます。ただし、学習に使った入力データセット $\{x_i\}$ を保持しておく必要があります。
以上に述べたカーネル回帰を実装すると、以下のようになります。
カーネル関数として、ガウスカーネル $\exp(-\beta(x-x')^2)$ を使っています。
# ガウスカーネル回帰
class KernelRegression:
# beta: ガウスカーネルパラメータ
# lam : 正則化係数
def __init__(self, beta, lam):
self.beta = beta
self.lam = lam
# gauss kernel
def __kernel(self, x1, x2):
return np.exp(-self.beta * (x1 - x2)**2)
# グラム行列を返す
def __gram_mat(self, x):
N = x.size
K = np.fromfunction(lambda i, j: self.__kernel(x[i], x[j]),
(N, N), dtype=int)
return K
# 学習データを入力して係数aを決定する
def fit(self, x, y):
# 学習データを保持する
self.x = x
# solve (K + lambda I) a = y
K = self.__gram_mat(x)
self.a = np.linalg.solve(K + self.lam * np.eye(x.size), y)
# 予測値を返す
def predict(self, x):
# \sum_j k(x, x_j) a_j
K = np.fromfunction(lambda i, j: self.__kernel(x[i], self.x[j]),
(x.size, self.x.size), dtype=int)
return K.dot(self.a)
先ほどと同じ問題で、カーネル回帰を実行した例を示します。ガウスカーネルのパラメータを $\beta=1.0$, 正則化係数を $\lambda=0.1$ としました。
まとめ
数式やコードが増えてしまいましたが、この記事で強調したかったのは、具体的な計算方法ではなく、線形回帰問題を線形方程式の一般的な問題として捉えたときの、最小二乗解とカーネル法の解の間の関係を対比的に見ることです。
- 最小二乗回帰 : $N>M$ が前提、最小二乗解 $(\Phi^T\Phi)^{-1}\Phi^Tt$、$M$ 次正方行列 $\Phi^T\Phi$ が主役
- カーネル回帰 : $N<M$ が前提、ノルム最小解 $\Phi^T(\Phi\Phi^T)^{-1}t$、$N$ 次正方行列 $\Phi\Phi^T$ が主役
実は主成分分析とカーネル主成分分析にも似たような関係があるのですが、それについてはまた別の記事で書ければと思います。