Sun Dec 01 2019

ガウス過程 (Gaussian Process)

ガウス過程とは何か

ガウス過程でググると上の方に www.jstage.jst.go.jp/article/isciesci/62/10/62_390/_pdf という神記事が見つかる。

ガウス過程の説明は色々あるが、私は結局この記事が一番理解できた。 そして結論だけが欲しかったら式 (21)-(23) と読むと分かる。

しかしこれはいわば学習データが \(N\) 個あるとすると \(N \times N\) 行列の逆行列を求めたりする必要があり、計算量的に辛いので、上手くデータを取捨選択して \(N\) を小さく抑えよう、みたいな話がそれより後ろに書いてある。

それは今回は飛ばす。

事前分布

import math
import numpy
import numpy.random
from numpy.linalg import inv


# 今回は RBF カーネルを共分散に使う
def k(x1, x2, beta):
    return math.exp(-beta * (x1 - x2) ** 2)

def gp(x, beta=0.1):
    mu = numpy.zeros(x.shape)  # 平均
    K = numpy.array([[k(xi, xj, beta) for xj in x] for xi in x])  # 共分散
    y = numpy.random.multivariate_normal(mu, K)
    return y

x = numpy.array(list(range(-5, 6)), dtype='f')
y = gp(x)
plot(x, y)

そうするとこんなものが出力される.

カーネルの中に出現した beta というパラメータが共分散の強さを表している. 数字が小さいほど共分散の影響が強く出るので, 近い \(x\) には近い値が出て, 関数が滑らかに見える. 逆に大きいほど影響は小さくなり, バラバラの値に見える.

それを表したのが次図. 特に 0 のときは完全に定数関数になってしまう。

事後分布

では \(N\) 個の点組 \(\{(x_i, y_i)\}\) にフィッティングしたあとのガウス過程はどうなっているかが実際の興味の対象である. これはさっきも言ったように資料の (21) を見ると計算式そのものが書いてある.

例えば観測データ（学習データ）として

# [-5, 6)
x = numpy.array(list(range(-5, 6)), dtype='f')
n = len(x)
y = 0.1 * x * x - x

が与えられたときの未知データ

# [6, 10)
x_unk = numpy.array(list(range(6, 10)), dtype='f')
m = len(x_unk)

への予測は次のように求まる.

# 観測データ内の共分散
Knn = numpy.array([[k(x[i], x[j]) for j in range(n)] for i in range(n)])
# 未知データとの間の共分散
Kmn = numpy.array([[k(x_unk[i], x[j]) for j in range(n)] for i in range(m)])
# 未知データ内の共分散
Kmm = numpy.array([[k(x_unk[i], x_unk[j]) for j in range(m)] for i in range(m)])

# 事後分布
mu = Kmn @ inv(Knn) @ y
vr = Kmm - Kmn @ inv(Knn) @ Kmn.T

# 予測値
y_pred = numpy.random.multivariate_normal(mu, vr)

結果を見ると,

かなり賢い.