PyTorch : GPyTorch tutorials : GPyTorch 回帰

PyTorch : GPyTorch tutorials : GPyTorch 回帰チュートリアル (翻訳)

翻訳 : (株)クラスキャット セールスインフォメーション
作成日時 : 11/22/2018 (0.1.0.rc5)

* 本ページは、GPyTorch のドキュメント tutorials : GPyTorch Regression Tutorial を
翻訳した上で適宜、補足説明したものです:

* サンプルコードの動作確認はしておりますが、必要な場合には適宜、追加改変しています。
* ご自由にリンクを張って頂いてかまいませんが、sales-info@classcat.com までご一報いただけると嬉しいです。

 

GPyTorch tutorials : GPyTorch 回帰チュートリアル

イントロダクション

この notebook では、最も単純な例を使用し、単純な関数上で RBF カーネル・ガウス過程を訓練し、GPyTorch の多くのデザインの特徴を示します。次の関数を 11 訓練サンプルでモデル化して
\[
\begin{align*}
y &= \sin(2\pi x) + \epsilon \\
\epsilon &\sim \mathcal{N}(0, 0.2)
\end{align*}
\]

そして 51 テストサンプル上でテストしていきます。

Note: この notebook はガウス過程の数学的背景を教えることを必ずしも意図していません、むしろ GPyTorch でどのように単純なものを訓練して予測を行なうかです。数学的扱いについては、Gaussian Processes for Machine Learning の 2 章が GP 回帰への非常に完全なイントロダクションを提供しています (このテキスト全体は高く推奨されます) :

http://www.gaussianprocess.org/gpml/chapters/RW2.pdf

import math
import torch
import gpytorch
from matplotlib import pyplot as plt

%matplotlib inline
%load_ext autoreload
%autoreload 2

 

訓練データのセットアップ

次のセルでは、この例のために訓練データをセットアップします。[0, 1] 上の 11 の規則的な間隔のポイントを使用します、その上で関数を評価して訓練ラベルを得るためにガウシアンノイズを追加します。

# Training data is 11 points in [0,1] inclusive regularly spaced
train_x = torch.linspace(0, 1, 100)
# True function is sin(2*pi*x) with Gaussian noise
train_y = torch.sin(train_x * (2 * math.pi)) + torch.randn(train_x.size()) * 0.2

 

モデルをセットアップする

次のセルは GPyTorch のユーザ定義ガウス過程モデルの最も重要な特徴を示します。GPyTorch の GP モデルの構築は様々な方法で異なります。

最初に多くの既存の GP パッケージとは対照的にユーザのために full GP モデルは提供しません。むしろ、迅速に一つを構築するために必要なツールを提供します。これは何故ならば、標準的な PyTorch でニューラルネットワークを構築することに類似して、どのような必要なコンポーネントでも含む柔軟性を持つことが重要だと私達は信じるからです。深層学習とガウス過程を結合する CIFAR10 深層カーネル学習 例のような、より複雑な例で見られるように、これはカスタムモデルの設計においてユーザに素晴らしい柔軟性を許容します。

殆どの GP 回帰モデルについて、次の GPyTorch オブジェクトを構築することが必要です :

  1. GP モデル (gpytorch.models.ExactGP) – これは殆どの推論を扱います。
  2. 尤度 (= Likelihood) (gpytorch.likelihoods.GaussianLikelihood) – これは GP 回帰のために使用される最も一般的な尤度です。
  3. 平均 (= Mean) – これは GP の事前平均を定義します。
    • どの平均を使用するべきか知らない場合には、gpytorch.means.ConstantMean() で始めましょう。
  1. Kernel – これは GP の事前共分散を定義します。
    • どのカーネルを使用するべきか知らない場合には、gpytorch.kernels.ScaleKernel(gpytorch.kernels.RBFKernel()) で始めましょう。
  1. 多変量正規分布 (= MultivariateNormal Distribution) (gpytorch.distributions.MultivariateNormal) – これは多変量正規分布を表わすために使用されるオブジェクトです。

 

GP モデル

GPyTorch のユーザ構築 (Exact, i.e. 非変分) GP モデルのコンポーネントは、大まかに言えば :

  1. __init__ メソッド、これは訓練データと尤度を取り、そしてどのようなオブジェクトであれモデルの forward メソッドのために必要なものを構築します。
  2. forward メソッドは、ある \(n \times d\) データ x を取りそして x で評価された事前平均と共分散を持つ MultivariateNormal を返します。換言すれば、GP の事前平均と共分散行列を表わすベクトル \(\mu(x)\) と \(n \times n\) 行列 \(K_{xx}\) を返します。

この仕様はモデルを定義するとき巨大な総量の柔軟性を残します。例えば、2 つの kernel を加算を通じて組み立てるために、kernel モジュールを直接加算できます :

self.covar_module = ScaleKernel(RBFKernel() + WhiteNoiseKernel())

あるいは forward メソッドで kernel の出力を加算できます :

covar_x = self.rbf_kernel_module(x) + self.white_noise_module(x)
# We will use the simplest form of GP model, exact inference
class ExactGPModel(gpytorch.models.ExactGP):
    def __init__(self, train_x, train_y, likelihood):
        super(ExactGPModel, self).__init__(train_x, train_y, likelihood)
        self.mean_module = gpytorch.means.ConstantMean()
        self.covar_module = gpytorch.kernels.ScaleKernel(gpytorch.kernels.RBFKernel())

    def forward(self, x):
        mean_x = self.mean_module(x)
        covar_x = self.covar_module(x)
        return gpytorch.distributions.MultivariateNormal(mean_x, covar_x)

# initialize likelihood and model
likelihood = gpytorch.likelihoods.GaussianLikelihood()
model = ExactGPModel(train_x, train_y, likelihood)

 

モデル・モード

殆どの PyTorch モジュールのように、ExactGP は .train() と .eval() モードを持ちます。- .train() モードはモデルのハイパーパラメータを最適化するためです。- .eval() モードはモデル事後を通して予測を計算するためです。

 

モデルを訓練する

次のセルでは、ガウス過程のハイパーパラメータを訓練するために Type-II MLE を使用して処理します。

多くの他の GP 実装と比較してここでの最も明白な違いは、標準的な PyTorch でのように、コア訓練ループはユーザにより書かれます。GPyTorch では、torch.optim からの標準的な PyTorch optimizer を利用し、そしてモデルの総ての訓練可能なパラメータはタイプ torch.nn.Parameter です。GP モデルは torch.nn.Module を直接拡張していますので、model.parameters() や model.named_parameters() 関数のようなメソッドへの呼び出しは期待するように PyTorch 由来です。

殆どの場合、下のボイラープレートなコードは上手く動作します。それは標準的な PyTorch 訓練ループと同じ基本的なコンポーネントを持ちます :

  1. 総てのパラメータ勾配をゼロにします。
  2. モデルを呼び出して損失を計算します。
  3. 勾配を書き込むために loss 上で backward を呼び出します。
  4. optimizer 上で step します。

けれども、カスタム訓練ループを定義すればより素晴らしい柔軟性を可能にします。例えば、訓練の各ステップでパラメータをセーブしたり、異なるパラメータのために異なる学習率を使用することは容易です (これは例えば深層 kernel 学習において有用です)。

# Find optimal model hyperparameters
model.train()
likelihood.train()

# Use the adam optimizer
optimizer = torch.optim.Adam([
    {'params': model.parameters()},  # Includes GaussianLikelihood parameters
], lr=0.1)

# "Loss" for GPs - the marginal log likelihood
mll = gpytorch.mlls.ExactMarginalLogLikelihood(likelihood, model)

training_iter = 50
for i in range(training_iter):
    # Zero gradients from previous iteration
    optimizer.zero_grad()
    # Output from model
    output = model(train_x)
    # Calc loss and backprop gradients
    loss = -mll(output, train_y)
    loss.backward()
    print('Iter %d/%d - Loss: %.3f   log_lengthscale: %.3f   log_noise: %.3f' % (
        i + 1, training_iter, loss.item(),
        model.covar_module.base_kernel.log_lengthscale.item(),
        model.likelihood.log_noise.item()
    ))
    optimizer.step()
Iter 1/50 - Loss: 1.084   log_lengthscale: 0.000   log_noise: 0.000
Iter 2/50 - Loss: 1.043   log_lengthscale: -0.100   log_noise: -0.100
Iter 3/50 - Loss: 1.004   log_lengthscale: -0.196   log_noise: -0.200
Iter 4/50 - Loss: 0.964   log_lengthscale: -0.293   log_noise: -0.300
Iter 5/50 - Loss: 0.922   log_lengthscale: -0.387   log_noise: -0.399
Iter 6/50 - Loss: 0.877   log_lengthscale: -0.479   log_noise: -0.499
Iter 7/50 - Loss: 0.825   log_lengthscale: -0.572   log_noise: -0.598
Iter 8/50 - Loss: 0.767   log_lengthscale: -0.667   log_noise: -0.698
Iter 9/50 - Loss: 0.705   log_lengthscale: -0.762   log_noise: -0.799
Iter 10/50 - Loss: 0.644   log_lengthscale: -0.860   log_noise: -0.899
Iter 11/50 - Loss: 0.590   log_lengthscale: -0.960   log_noise: -1.001
Iter 12/50 - Loss: 0.543   log_lengthscale: -1.058   log_noise: -1.102
Iter 13/50 - Loss: 0.502   log_lengthscale: -1.150   log_noise: -1.204
Iter 14/50 - Loss: 0.462   log_lengthscale: -1.234   log_noise: -1.306
Iter 15/50 - Loss: 0.426   log_lengthscale: -1.303   log_noise: -1.408
Iter 16/50 - Loss: 0.389   log_lengthscale: -1.360   log_noise: -1.509
Iter 17/50 - Loss: 0.360   log_lengthscale: -1.404   log_noise: -1.611
Iter 18/50 - Loss: 0.321   log_lengthscale: -1.432   log_noise: -1.712
Iter 19/50 - Loss: 0.280   log_lengthscale: -1.454   log_noise: -1.812
Iter 20/50 - Loss: 0.250   log_lengthscale: -1.465   log_noise: -1.911
Iter 21/50 - Loss: 0.227   log_lengthscale: -1.469   log_noise: -2.010
Iter 22/50 - Loss: 0.188   log_lengthscale: -1.461   log_noise: -2.108
Iter 23/50 - Loss: 0.158   log_lengthscale: -1.442   log_noise: -2.204
Iter 24/50 - Loss: 0.125   log_lengthscale: -1.411   log_noise: -2.300
Iter 25/50 - Loss: 0.095   log_lengthscale: -1.377   log_noise: -2.393
Iter 26/50 - Loss: 0.070   log_lengthscale: -1.340   log_noise: -2.485
Iter 27/50 - Loss: 0.050   log_lengthscale: -1.298   log_noise: -2.574
Iter 28/50 - Loss: 0.032   log_lengthscale: -1.256   log_noise: -2.662
Iter 29/50 - Loss: 0.014   log_lengthscale: -1.218   log_noise: -2.746
Iter 30/50 - Loss: 0.003   log_lengthscale: -1.182   log_noise: -2.828
Iter 31/50 - Loss: -0.001   log_lengthscale: -1.148   log_noise: -2.906
Iter 32/50 - Loss: -0.008   log_lengthscale: -1.121   log_noise: -2.980
Iter 33/50 - Loss: -0.012   log_lengthscale: -1.102   log_noise: -3.049
Iter 34/50 - Loss: -0.011   log_lengthscale: -1.103   log_noise: -3.114
Iter 35/50 - Loss: -0.014   log_lengthscale: -1.114   log_noise: -3.174
Iter 36/50 - Loss: -0.014   log_lengthscale: -1.138   log_noise: -3.228
Iter 37/50 - Loss: -0.010   log_lengthscale: -1.169   log_noise: -3.275
Iter 38/50 - Loss: -0.011   log_lengthscale: -1.204   log_noise: -3.317
Iter 39/50 - Loss: -0.008   log_lengthscale: -1.239   log_noise: -3.352
Iter 40/50 - Loss: -0.001   log_lengthscale: -1.270   log_noise: -3.380
Iter 41/50 - Loss: -0.005   log_lengthscale: -1.296   log_noise: -3.401
Iter 42/50 - Loss: 0.008   log_lengthscale: -1.317   log_noise: -3.415
Iter 43/50 - Loss: 0.001   log_lengthscale: -1.331   log_noise: -3.422
Iter 44/50 - Loss: 0.009   log_lengthscale: -1.343   log_noise: -3.423
Iter 45/50 - Loss: 0.001   log_lengthscale: -1.350   log_noise: -3.419
Iter 46/50 - Loss: -0.001   log_lengthscale: -1.360   log_noise: -3.410
Iter 47/50 - Loss: 0.007   log_lengthscale: -1.374   log_noise: -3.397
Iter 48/50 - Loss: 0.000   log_lengthscale: -1.388   log_noise: -3.380
Iter 49/50 - Loss: -0.010   log_lengthscale: -1.396   log_noise: -3.359
Iter 50/50 - Loss: -0.008   log_lengthscale: -1.404   log_noise: -3.337

 

モデルで予測する

次のセルでは、モデルで予測を行ないます。これを行なうには、単純にモデルと尤度を eval モードにして、テストデータ上で両者のモジュールを呼び出します。

ちょうどユーザ定義の GP モデルが forward から事前平均と共分散を含む MultivariateNormal を返すように、eval モードの訓練された GP モデルは 事後平均と共分散を含む MultivariateNormal を返します。こうして、予測平均と分散を得て、それから与えられたテストポイントにおける GP からのサンプリング関数は次のような呼び出しで達成されるでしょう :

f_preds = model(test_x)
y_preds = likelihood(model(test_x))

f_mean = f_preds.mean
f_var = f_preds.variance
f_covar = f_preds.covariance_matrix
f_samples = f_preds.sample(sample_shape=torch.Size(1000,))

gpytorch.fast_pred_var コンテキストは必要ではありませんが、ここでは LOVE を使用してより高速な予測分布を得て、クールな特徴の一つを使用するプレビューを与えています :

# Get into evaluation (predictive posterior) mode
model.eval()
likelihood.eval()

# Test points are regularly spaced along [0,1]
# Make predictions by feeding model through likelihood
with torch.no_grad(), gpytorch.fast_pred_var():
    test_x = torch.linspace(0, 1, 51)
    observed_pred = likelihood(model(test_x))

 

モデル fit をプロットする

次のセルでは、ガウス過程モデルの平均と信頼 (= confidence) 領域をプロットします。confidence_region は平均の上と下の 2 標準偏差を返すヘルパー・メソッドです。

with torch.no_grad():
    # Initialize plot
    f, ax = plt.subplots(1, 1, figsize=(4, 3))

    # Get upper and lower confidence bounds
    lower, upper = observed_pred.confidence_region()
    # Plot training data as black stars
    ax.plot(train_x.numpy(), train_y.numpy(), 'k*')
    # Plot predictive means as blue line
    ax.plot(test_x.numpy(), observed_pred.mean.numpy(), 'b')
    # Shade between the lower and upper confidence bounds
    ax.fill_between(test_x.numpy(), lower.numpy(), upper.numpy(), alpha=0.5)
    ax.set_ylim([-3, 3])
    ax.legend(['Observed Data', 'Mean', 'Confidence'])


 

以上