PyTorch教程-18.3。高斯過程推理

在本節(jié)中，我們將展示如何使用我們?cè)谏弦还?jié)中介紹的 GP 先驗(yàn)執(zhí)行后驗(yàn)推理和進(jìn)行預(yù)測(cè)。我們將從回歸開始，我們可以在其中以封閉形式執(zhí)行推理。這是一個(gè)“簡(jiǎn)而言之 GP”部分，可在實(shí)踐中快速啟動(dòng)和運(yùn)行高斯過程。我們將從頭開始編寫所有基本操作的代碼，然后介紹 GPyTorch，這將使使用最先進(jìn)的高斯過程以及與深度神經(jīng)網(wǎng)絡(luò)的集成更加方便。我們將在下一節(jié)中深入探討這些更高級(jí)的主題。在該部分中，我們還將考慮需要近似推理的設(shè)置——分類、點(diǎn)過程或任何非高斯似然。

18.3.1?；貧w的后驗(yàn)推理

觀察模型與我們想要學(xué)習(xí)的功能相關(guān)聯(lián)， f(x), 根據(jù)我們的觀察y(x), 都由一些輸入索引x. 在分類上，x可以是圖像的像素，并且y可能是關(guān)聯(lián)的類標(biāo)簽。在回歸中， y通常表示連續(xù)輸出，例如地表溫度、海平面、CO2濃度等

在回歸中，我們通常假設(shè)輸出是由潛在的無噪聲函數(shù)給出的f(x)加上 iid 高斯噪聲 ?(x):

(18.3.1)y(x)=f(x)+?(x),

和?(x)～N(0,σ2). 讓 y=y(X)=(y(x1),…,y(xn))?是我們訓(xùn)練觀察的向量，并且 f=(f(x1),…,f(xn))?是潛在無噪聲函數(shù)值的向量，在訓(xùn)練輸入中查詢 X=x1,…,xn.

我們假設(shè)f(x)～GP(m,k)，這意味著任何函數(shù)值的集合f具有聯(lián)合多元高斯分布，均值向量 μi=m(xi)和協(xié)方差矩陣 Kij=k(xi,xj). RBF核 k(xi,xj)=a2exp?(?12?2||xi?xj||2) 將是協(xié)方差函數(shù)的標(biāo)準(zhǔn)選擇。為了符號(hào)簡(jiǎn)單起見，我們將假設(shè)均值函數(shù)m(x)=0; 我們的推導(dǎo)可以在以后很容易地推廣。

假設(shè)我們想對(duì)一組輸入進(jìn)行預(yù)測(cè)

(18.3.2)X?=x?1,x?2,…,x?m.

那么我們想要找到x2和 p(f?|y,X). 在回歸設(shè)置中，我們可以在找到聯(lián)合分布后，使用高斯恒等式方便地找到這個(gè)分布f?=f(X?) 和y.

如果我們?cè)谟?xùn)練輸入處評(píng)估等式(18.3.1)X， 我們有 y=f+?. 根據(jù)高斯過程的定義（見上一節(jié)）， f～N(0,K(X,X))在哪里K(X,X)是一個(gè) n×n通過在所有可能的輸入對(duì)上評(píng)估我們的協(xié)方差函數(shù)（又名kernel ）形成的矩陣xi,xj∈X. ?只是一個(gè)由 iid 樣本組成的向量N(0,σ2)從而有分布 N(0,σ2I).y因此是兩個(gè)獨(dú)立的多元高斯變量的總和，因此具有分布N(0,K(X,X)+σ2I). 也可以證明 cov(f?,y)=cov(y,f?)?=K(X?,X) 在哪里K(X?,X)是一個(gè)m×n通過評(píng)估所有對(duì)測(cè)試和訓(xùn)練輸入的內(nèi)核形成的矩陣。

(18.3.3)[yf?]～N(0,A=[K(X,X)+σ2IK(X,X?)K(X?,X)K(X?,X?)])

然后我們可以使用標(biāo)準(zhǔn)高斯恒等式從聯(lián)合分布中找到條件分布（例如，參見 Bishop 第 2 章）， f?|y,X,X?～N(m?,S?)， 在哪里m?=K(X?,X)[K(X,X)+σ2I]?1y， 和 S=K(X?,X?)?K(X?,X)[K(X,X)+σ2I]?1K(X,X?).

通常，我們不需要使用完整的預(yù)測(cè)協(xié)方差矩陣S, 而不是使用對(duì)角線S對(duì)于每個(gè)預(yù)測(cè)的不確定性。通常出于這個(gè)原因，我們?yōu)閱蝹€(gè)測(cè)試點(diǎn)編寫預(yù)測(cè)分布x?，而不是測(cè)試點(diǎn)的集合。

核矩陣有參數(shù)θ我們也希望估計(jì)，例如振幅a和長(zhǎng)度尺度?上面的 RBF 內(nèi)核。為了這些目的，我們使用邊際似然，p(y|θ,X)，我們已經(jīng)在計(jì)算邊際分布以找到聯(lián)合分布時(shí)得出y,f?. 正如我們將看到的，邊際似然劃分為模型擬合和模型復(fù)雜性項(xiàng)，并自動(dòng)編碼用于學(xué)習(xí)超參數(shù)的奧卡姆剃刀概念。有關(guān)完整討論，請(qǐng)參閱 MacKay Ch。28 （MacKay 和 Mac Kay，2003 年），以及 Rasmussen 和 Williams Ch。5 （拉斯穆森和威廉姆斯，2006 年）。

import math
import os
import gpytorch
import matplotlib.pyplot as plt
import numpy as np
import torch
from scipy import optimize
from scipy.spatial import distance_matrix
from d2l import torch as d2l

d2l.set_figsize()

18.3.2。在 GP 回歸中進(jìn)行預(yù)測(cè)和學(xué)習(xí)內(nèi)核超參數(shù)的方程式

我們?cè)谶@里列出了您將用于學(xué)習(xí)超參數(shù)和在高斯過程回歸中進(jìn)行預(yù)測(cè)的方程式。同樣，我們假設(shè)一個(gè)回歸目標(biāo)向量y, 由輸入索引 X={x1,…,xn}，我們希望在測(cè)試輸入中做出預(yù)測(cè)x?. 我們假設(shè)具有方差的 iid 加性零均值高斯噪聲σ2. 我們先使用高斯過程 f(x)～GP(m,k)對(duì)于潛在的無噪聲函數(shù)，具有均值函數(shù)m和核函數(shù)k. 內(nèi)核本身有參數(shù)θ我們想學(xué)習(xí)。例如，如果我們使用 RBF 內(nèi)核， k(xi,xj)=a2exp?(?12?2||x?x′||2), 我們想學(xué)習(xí)θ={a2,?2}. 讓K(X,X) 代表一個(gè)n×n對(duì)應(yīng)于為所有可能的對(duì)評(píng)估內(nèi)核的矩陣n培訓(xùn)投入。讓 K(x?,X)代表一個(gè)1×n通過評(píng)估形成的向量k(x?,xi),i=1,…,n. 讓?duì)淌峭ㄟ^計(jì)算均值函數(shù)形成的均值向量m(x)在每個(gè)訓(xùn)練點(diǎn)x.

通常在處理高斯過程時(shí)，我們遵循一個(gè)兩步過程。1.學(xué)習(xí)內(nèi)核超參數(shù)θ^通過最大化這些超參數(shù)的邊際似然。2.使用預(yù)測(cè)均值作為點(diǎn)預(yù)測(cè)器，2倍預(yù)測(cè)標(biāo)準(zhǔn)差形成95%可信集，以這些學(xué)習(xí)到的超參數(shù)為條件θ^.

對(duì)數(shù)邊際似然只是一個(gè)對(duì)數(shù)高斯密度，其形式為：

(18.3.4)log?p(y|θ,X)=?12y?[Kθ(X,X)+σ2I]?1y?12log?|Kθ(X,X)|+c

預(yù)測(cè)分布具有以下形式：

(18.3.5)p(y?|x?,y,θ)=N(a?,v?)

(18.3.6)a?=kθ(x?,X)[Kθ(X,X)+σ2I]?1(y?μ)+μ

(18.3.7)v?=kθ(x?,x?)?Kθ(x?,X)[Kθ(X,X)+σ2I]?1kθ(X,x?)

18.3.3。解釋學(xué)習(xí)和預(yù)測(cè)方程

關(guān)于高斯過程的預(yù)測(cè)分布，有一些關(guān)鍵點(diǎn)需要注意：

盡管模型類具有靈活性，但可以 對(duì)封閉形式的 GP 回歸進(jìn)行精確的貝葉斯推理。除了學(xué)習(xí)內(nèi)核超參數(shù)外，沒有任何訓(xùn)練。我們可以準(zhǔn)確地寫下我們想要用來進(jìn)行預(yù)測(cè)的方程式。高斯過程在這方面相對(duì)特殊，這極大地促進(jìn)了它們的便利性、多功能性和持續(xù)流行。

預(yù)測(cè)均值a?是訓(xùn)練目標(biāo)的線性組合y, 由內(nèi)核加權(quán) kθ(x?,X)[Kθ(X,X)+σ2I]?1. 正如我們將看到的，內(nèi)核（及其超參數(shù)）因此在模型的泛化特性中起著至關(guān)重要的作用。

預(yù)測(cè)均值明確取決于目標(biāo)值 y但預(yù)測(cè)方差沒有。預(yù)測(cè)不確定性反而隨著測(cè)試輸入而增長(zhǎng)x? 遠(yuǎn)離目標(biāo)位置X，由內(nèi)核函數(shù)控制。然而，不確定性將隱含地取決于目標(biāo)的價(jià)值y通過內(nèi)核超參數(shù)θ，這是從數(shù)據(jù)中學(xué)到的。

邊際似然劃分為模型擬合和模型復(fù)雜性（對(duì)數(shù)行列式）項(xiàng)。邊際似然傾向于選擇提供與數(shù)據(jù)仍然一致的最簡(jiǎn)單擬合的超參數(shù)。

關(guān)鍵的計(jì)算瓶頸來自解決線性系統(tǒng)和計(jì)算對(duì)數(shù)行列式n×n對(duì)稱正定矩陣K(X,X)為了n訓(xùn)練點(diǎn)。天真地，這些操作每個(gè)都會(huì)招致O(n3) 計(jì)算，以及O(n2)核（協(xié)方差）矩陣的每個(gè)條目的存儲(chǔ)，通常從 Cholesky 分解開始。從歷史上看，這些瓶頸將 GP 限制在訓(xùn)練點(diǎn)少于 10,000 的問題上，并給 GP 帶來了近十年來不準(zhǔn)確的“速度慢”的名聲。在高級(jí)主題中，我們將討論如何將 GP 擴(kuò)展到具有數(shù)百萬(wàn)個(gè)點(diǎn)的問題。

對(duì)于核函數(shù)的流行選擇，K(X,X)通常接近奇異，這可能會(huì)在執(zhí)行 Cholesky 分解或其他旨在求解線性系統(tǒng)的操作時(shí)導(dǎo)致數(shù)值問題。幸運(yùn)的是，在回歸中我們經(jīng)常使用 Kθ(X,X)+σ2I，使得噪聲方差 σ2被添加到對(duì)角線K(X,X), 顯著改善其條件。如果噪聲方差很小，或者我們正在進(jìn)行無噪聲回歸，通常的做法是在對(duì)角線上添加少量“抖動(dòng)”，數(shù)量級(jí)為 10?6, 改善調(diào)理。

18.3.4。從零開始的工作示例

讓我們創(chuàng)建一些回歸數(shù)據(jù)，然后用 GP 擬合數(shù)據(jù)，從頭開始實(shí)施每一步。我們將從中抽樣數(shù)據(jù)

(18.3.8)y(x)=sin?(x)+12sin?(4x)+?,

和?～N(0,σ2). 我們希望找到的無噪聲函數(shù)是 f(x)=sin?(x)+12sin?(4x). 我們將從使用噪聲標(biāo)準(zhǔn)偏差開始σ=0.25.

def data_maker1(x, sig):
  return np.sin(x) + 0.5 * np.sin(4 * x) + np.random.randn(x.shape[0]) * sig

sig = 0.25
train_x, test_x = np.linspace(0, 5, 50), np.linspace(0, 5, 500)
train_y, test_y = data_maker1(train_x, sig=sig), data_maker1(test_x, sig=0.)

d2l.plt.scatter(train_x, train_y)
d2l.plt.plot(test_x, test_y)
d2l.plt.xlabel("x", fontsize=20)
d2l.plt.ylabel("Observations y", fontsize=20)
d2l.plt.show()

在這里，我們將嘈雜的觀察結(jié)果視為圓圈，將我們希望找到的藍(lán)色無噪聲函數(shù)視為圓圈。

現(xiàn)在，讓我們?cè)跐撛跓o噪聲函數(shù)之前指定一個(gè) GP， f(x)～GP(m,k). 我們將使用均值函數(shù) m(x)=0和 RBF 協(xié)方差函數(shù)（核）

(18.3.9)k(xi,xj)=a2exp?(?12?2||x?x′||2).

mean = np.zeros(test_x.shape[0])
cov = d2l.rbfkernel(test_x, test_x, ls=0.2)

我們從 0.2 的長(zhǎng)度尺度開始。在我們擬合數(shù)據(jù)之前，重要的是要考慮我們是否指定了合理的先驗(yàn)。讓我們想象一下這個(gè)先驗(yàn)的一些樣本函數(shù)，以及 95% 的可信集（我們相信真正的函數(shù)有 95% 的可能性在這個(gè)區(qū)域內(nèi)）。

prior_samples = np.random.multivariate_normal(mean=mean, cov=cov, size=5)
d2l.plt.plot(test_x, prior_samples.T, color='black', alpha=0.5)
d2l.plt.plot(test_x, mean, linewidth=2.)
d2l.plt.fill_between(test_x, mean - 2 * np.diag(cov), mean + 2 * np.diag(cov),
         alpha=0.25)
d2l.plt.show()

這些樣本看起來合理嗎？函數(shù)的高級(jí)屬性是否與我們嘗試建模的數(shù)據(jù)類型一致？

現(xiàn)在讓我們?cè)谌我鉁y(cè)試點(diǎn)形成后驗(yàn)預(yù)測(cè)分布的均值和方差x?.

(18.3.10)fˉ?=K(x,x?)T(K(x,x)+σ2I)?1y

(18.3.11)V(f?)=K(x?,x?)?K(x,x?)T(K(x,x)+σ2I)?1K(x,x?)

在我們做出預(yù)測(cè)之前，我們應(yīng)該學(xué)習(xí)我們的內(nèi)核超參數(shù) θ和噪聲方差σ2. 讓我們將長(zhǎng)度尺度初始化為 0.75，因?yàn)榕c我們正在擬合的數(shù)據(jù)相比，我們之前的函數(shù)看起來變化太快了。我們還將猜測(cè)噪聲標(biāo)準(zhǔn)偏差σ0.75。

為了學(xué)習(xí)這些參數(shù)，我們將最大化這些參數(shù)的邊際似然。

(18.3.12)log?p(y|X)=log?∫p(y|f,X)p(f|X)df

(18.3.13)log?p(y|X)=?12yT(K(x,x)+σ2I)?1y?12log?|K(x,x)+σ2I|?n2log?2π

也許我們之前的功能變化太快了。我們假設(shè)長(zhǎng)度尺度為 0.4。我們還將猜測(cè)噪聲標(biāo)準(zhǔn)偏差為 0.75。這些只是超參數(shù)初始化——我們將從邊際似然中學(xué)習(xí)這些參數(shù)。

ell_est = 0.4
post_sig_est = 0.5

def neg_MLL(pars):
  K = d2l.rbfkernel(train_x, train_x, ls=pars[0])
  kernel_term = -0.5 * train_y @ 
    np.linalg.inv(K + pars[1] ** 2 * np.eye(train_x.shape[0])) @ train_y
  logdet = -0.5 * np.log(np.linalg.det(K + pars[1] ** 2 * 
                     np.eye(train_x.shape[0])))
  const = -train_x.shape[0] / 2. * np.log(2 * np.pi)

  return -(kernel_term + logdet + const)


learned_hypers = optimize.minimize(neg_MLL, x0=np.array([ell_est,post_sig_est]),
                  bounds=((0.01, 10.), (0.01, 10.)))
ell = learned_hypers.x[0]
post_sig_est = learned_hypers.x[1]

在本例中，我們學(xué)習(xí)了 0.299 的長(zhǎng)度尺度和 0.24 的噪聲標(biāo)準(zhǔn)差。請(qǐng)注意，學(xué)習(xí)到的噪聲非常接近真實(shí)噪聲，這有助于表明我們的 GP 非常明確地解決了這個(gè)問題。

一般來說，仔細(xì)考慮選擇內(nèi)核和初始化超參數(shù)是至關(guān)重要的。雖然邊際似然優(yōu)化對(duì)初始化來說相對(duì)穩(wěn)健，但它也不能避免糟糕的初始化。嘗試使用各種初始化運(yùn)行上面的腳本，看看你發(fā)現(xiàn)了什么結(jié)果。

現(xiàn)在，讓我們用這些學(xué)習(xí)到的超級(jí)計(jì)算機(jī)進(jìn)行預(yù)測(cè)。

K_x_xstar = d2l.rbfkernel(train_x, test_x, ls=ell)
K_x_x = d2l.rbfkernel(train_x, train_x, ls=ell)
K_xstar_xstar = d2l.rbfkernel(test_x, test_x, ls=ell)

post_mean = K_x_xstar.T @ np.linalg.inv((K_x_x + 
        post_sig_est ** 2 * np.eye(train_x.shape[0]))) @ train_y
post_cov = K_xstar_xstar - K_x_xstar.T @ np.linalg.inv((K_x_x + 
        post_sig_est ** 2 * np.eye(train_x.shape[0]))) @ K_x_xstar

lw_bd = post_mean - 2 * np.sqrt(np.diag(post_cov))
up_bd = post_mean + 2 * np.sqrt(np.diag(post_cov))

d2l.plt.scatter(train_x, train_y)
d2l.plt.plot(test_x, test_y, linewidth=2.)
d2l.plt.plot(test_x, post_mean, linewidth=2.)
d2l.plt.fill_between(test_x, lw_bd, up_bd, alpha=0.25)
d2l.plt.legend(['Observed Data', 'True Function', 'Predictive Mean', '95% Set on True Func'])
d2l.plt.show()

我們看到橙色的后驗(yàn)均值幾乎完全匹配真正的無噪聲函數(shù)！請(qǐng)注意，我們顯示的 95% 可信集是針對(duì)無潛在噪聲（真實(shí)）函數(shù)的，而不是數(shù)據(jù)點(diǎn)。我們看到這個(gè)可信集完全包含了真實(shí)的功能，并且看起來并不過分寬泛或狹窄。我們不希望也不期望它包含數(shù)據(jù)點(diǎn)。如果我們希望有一個(gè)可靠的觀察集，我們應(yīng)該計(jì)算

lw_bd_observed = post_mean - 2 * np.sqrt(np.diag(post_cov) + post_sig_est ** 2)
up_bd_observed = post_mean + 2 * np.sqrt(np.diag(post_cov) + post_sig_est ** 2)

不確定性有兩個(gè)來源，認(rèn)知不確定性，代表可減少的不確定性，以及任意或不可減少的 不確定性。這里的認(rèn)知不確定性表示無噪聲函數(shù)真實(shí)值的不確定性。隨著我們遠(yuǎn)離數(shù)據(jù)點(diǎn)，這種不確定性應(yīng)該會(huì)增加，因?yàn)檫h(yuǎn)離數(shù)據(jù)有更多種類的函數(shù)值與我們的數(shù)據(jù)一致。隨著我們觀察到越來越多的數(shù)據(jù)，我們對(duì)真實(shí)函數(shù)的信念變得更加自信，認(rèn)知不確定性消失。這種情況下的任意 不確定性是觀察噪聲，因?yàn)閿?shù)據(jù)是通過這種噪聲提供給我們的，并且無法減少。

數(shù)據(jù)中的認(rèn)知不確定性由潛在無噪聲函數(shù) np.diag(post_cov) 的方差捕獲。任意不 確定性由噪聲方差 post_sig_est**2 捕獲。

不幸的是，人們常常不關(guān)心他們?nèi)绾伪硎静淮_定性，許多論文顯示的誤差條完全未定義，不清楚我們是在可視化認(rèn)知不確定性還是任意不確定性或兩者兼而有之，并且將噪聲方差與噪聲標(biāo)準(zhǔn)差混淆，標(biāo)準(zhǔn)差與標(biāo)準(zhǔn)誤差、可信集的置信區(qū)間等。如果不明確不確定性代表什么，它本質(zhì)上是沒有意義的。

本著密切關(guān)注我們的不確定性代表什么的精神，重要的是要注意我們對(duì)無噪聲函數(shù)的方差估計(jì)取平方根的 兩倍。由于我們的預(yù)測(cè)分布是高斯分布，這個(gè)數(shù)量使我們能夠形成一個(gè) 95% 的可信集，代表我們對(duì) 95% 可能包含基本事實(shí)函數(shù)的區(qū)間的信念。噪聲 方差存在于完全不同的尺度上，并且更難解釋。

最后，讓我們來看看20個(gè)后驗(yàn)樣本。這些樣本告訴我們我們認(rèn)為哪些類型的函數(shù)可能適合我們的數(shù)據(jù)，這是后驗(yàn)的。

post_samples = np.random.multivariate_normal(post_mean, post_cov, size=20)
d2l.plt.scatter(train_x, train_y)
d2l.plt.plot(test_x, test_y, linewidth=2.)
d2l.plt.plot(test_x, post_mean, linewidth=2.)
d2l.plt.plot(test_x, post_samples.T, color='gray', alpha=0.25)
d2l.plt.fill_between(test_x, lw_bd, up_bd, alpha=0.25)
plt.legend(['Observed Data', 'True Function', 'Predictive Mean', 'Posterior Samples'])
d2l.plt.show()

在基本回歸應(yīng)用中，最常見的做法是分別使用后驗(yàn)預(yù)測(cè)均值和標(biāo)準(zhǔn)差作為不確定性的點(diǎn)預(yù)測(cè)器和度量。在更高級(jí)的應(yīng)用程序中，例如使用蒙特卡洛采集函數(shù)的貝葉斯優(yōu)化或基于模型的 RL 的高斯過程，通常需要采用后驗(yàn)樣本。然而，即使在基本應(yīng)用程序中沒有嚴(yán)格要求，這些樣本也能讓我們更直觀地了解我們對(duì)數(shù)據(jù)的擬合程度，并且通常有助于將其包含在可視化中。

18.3.5。使用 GPyTorch 讓生活更輕松

正如我們所見，完全從頭開始實(shí)現(xiàn)基本的高斯過程回歸實(shí)際上非常容易。然而，一旦我們想要探索各種內(nèi)核選擇，考慮近似推理（這甚至是分類所必需的），將 GP 與神經(jīng)網(wǎng)絡(luò)相結(jié)合，甚至擁有大于約 10,000 個(gè)點(diǎn)的數(shù)據(jù)集，那么從頭開始的實(shí)現(xiàn)就變成了笨重和麻煩。一些最有效的可擴(kuò)展 GP 推理方法，例如 SKI（也稱為 KISS-GP），可能需要數(shù)百行代碼來實(shí)現(xiàn)高級(jí)數(shù)值線性代數(shù)例程。

在這些情況下，GPyTorch庫(kù)將使我們的生活變得更加輕松。我們將在未來關(guān)于高斯過程數(shù)值和高級(jí)方法的筆記本中更多地討論 GPyTorch。GPyTorch 庫(kù)包含 許多示例。為了感受這個(gè)包，我們將通過簡(jiǎn)單的回歸示例，展示如何使用 GPyTorch 對(duì)其進(jìn)行調(diào)整以重現(xiàn)我們的上述結(jié)果。這看起來像是簡(jiǎn)單地重現(xiàn)上面的基本回歸的大量代碼，從某種意義上說，它確實(shí)如此。但是我們可以立即使用各種內(nèi)核、可擴(kuò)展推理技術(shù)和近似推理，只需從下面更改幾行代碼，而不是編寫可能有數(shù)千行的新代碼。

# First let's convert our data into tensors for use with PyTorch
train_x = torch.tensor(train_x)
train_y = torch.tensor(train_y)
test_y = torch.tensor(test_y)

# We are using exact GP inference with a zero mean and RBF kernel
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.ZeroMean()
    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)

此代碼塊將數(shù)據(jù)置于 GPyTorch 的正確格式中，并指定我們正在使用精確推理，以及我們要使用的均值函數(shù)（零）和核函數(shù)（RBF）。我們可以很容易地使用任何其他內(nèi)核，例如，通過調(diào)用 gpytorch.kernels.matern_kernel() 或 gpyotrch.kernels.spectral_mixture_kernel()。到目前為止，我們只討論了精確推理，其中可以在不進(jìn)行任何近似的情況下推斷出預(yù)測(cè)分布。對(duì)于高斯過程，我們只有在具有高斯似然時(shí)才能進(jìn)行精確推理；更具體地說，當(dāng)我們假設(shè)我們的觀察結(jié)果是由高斯過程和高斯噪聲表示的無噪聲函數(shù)生成的。在未來的筆記本中，我們將考慮其他設(shè)置，例如分類，我們無法做出這些假設(shè)。

# Initialize Gaussian likelihood
likelihood = gpytorch.likelihoods.GaussianLikelihood()
model = ExactGPModel(train_x, train_y, likelihood)
training_iter = 50
# Find optimal model hyperparameters
model.train()
likelihood.train()
# Use the adam optimizer, includes GaussianLikelihood parameters
optimizer = torch.optim.Adam(model.parameters(), lr=0.1)
# Set our loss as the negative log GP marginal likelihood
mll = gpytorch.mlls.ExactMarginalLogLikelihood(likelihood, model)

在這里，我們明確指定我們想要使用的可能性（高斯），我們將用于訓(xùn)練內(nèi)核超參數(shù)的目標(biāo)（這里是邊際可能性），以及我們想要用于優(yōu)化該目標(biāo)的過程（在本例中，Adam ). 我們注意到，雖然我們使用的是“隨機(jī)”優(yōu)化器 Adam，但在本例中，它是全批次 Adam。因?yàn)檫呺H似然不會(huì)分解數(shù)據(jù)實(shí)例，所以我們不能對(duì)“小批量”數(shù)據(jù)使用優(yōu)化器并保證收斂。GPyTorch 也支持其他優(yōu)化器，例如 L-BFGS。與標(biāo)準(zhǔn)深度學(xué)習(xí)不同，做好邊際似然優(yōu)化工作與良好的泛化密切相關(guān)，這通常使我們傾向于使用像 L-BFGS 這樣強(qiáng)大的優(yōu)化器，前提是它們的成本并不高得令人望而卻步。

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()
  if i % 10 == 0:
    print(f'Iter {i+1:d}/{training_iter:d} - Loss: {loss.item():.3f} '
       f'squared lengthscale: '
       f'{model.covar_module.base_kernel.lengthscale.item():.3f} '
       f'noise variance: {model.likelihood.noise.item():.3f}')
  optimizer.step()

Iter 1/50 - Loss: 0.973 squared lengthscale: 0.693 noise variance: 0.693
Iter 11/50 - Loss: 0.684 squared lengthscale: 0.511 noise variance: 0.313
Iter 21/50 - Loss: 0.422 squared lengthscale: 0.533 noise variance: 0.128
Iter 31/50 - Loss: 0.304 squared lengthscale: 0.535 noise variance: 0.056
Iter 41/50 - Loss: 0.320 squared lengthscale: 0.522 noise variance: 0.040

在這里，我們實(shí)際運(yùn)行優(yōu)化程序，每 10 次迭代輸出一次損失值。

# Get into evaluation (predictive posterior) mode
test_x = torch.tensor(test_x)
model.eval()
likelihood.eval()
observed_pred = likelihood(model(test_x))

上面的代碼塊使我們能夠?qū)ξ覀兊臏y(cè)試輸入進(jìn)行預(yù)測(cè)。

with torch.no_grad():
  # Initialize plot
  f, ax = d2l.plt.subplots(1, 1, figsize=(4, 3))
  # Get upper and lower bounds for 95% credible set (in this case, in
  # observation space)
  lower, upper = observed_pred.confidence_region()
  ax.scatter(train_x.numpy(), train_y.numpy())
  ax.plot(test_x.numpy(), test_y.numpy(), linewidth=2.)
  ax.plot(test_x.numpy(), observed_pred.mean.numpy(), linewidth=2.)
  ax.fill_between(test_x.numpy(), lower.numpy(), upper.numpy(), alpha=0.25)
  ax.set_ylim([-1.5, 1.5])
  ax.legend(['True Function', 'Predictive Mean', 'Observed Data',
        '95% Credible Set'])

最后，我們繪制擬合圖。

我們看到擬合幾乎相同。需要注意的幾件事：GPyTorch 正在處理平方長(zhǎng)度尺度和觀察噪聲。例如，我們?cè)?for scratch 代碼中學(xué)習(xí)到的噪聲標(biāo)準(zhǔn)差約為 0.283。GPyTorch 找到的噪聲方差是 0.81≈0.2832. 在 GPyTorch 圖中，我們還展示了觀察空間中的可信集而不是潛在函數(shù)空間，以證明它們確實(shí)覆蓋了觀察到的數(shù)據(jù)點(diǎn)。

18.3.6。概括

我們可以將高斯過程先驗(yàn)與數(shù)據(jù)結(jié)合起來形成后驗(yàn)，我們用它來進(jìn)行預(yù)測(cè)。我們還可以形成邊際似然，這對(duì)于核超參數(shù)的自動(dòng)學(xué)習(xí)很有用，核超參數(shù)控制高斯過程的變化率等屬性。為回歸形成后驗(yàn)和學(xué)習(xí)內(nèi)核超參數(shù)的機(jī)制很簡(jiǎn)單，涉及大約十幾行代碼。對(duì)于任何想要快速“啟動(dòng)并運(yùn)行”高斯過程的讀者來說，這個(gè)筆記本都是一個(gè)很好的參考。我們還介紹了 GPyTorch 庫(kù)。盡管用于基本回歸的 GPyTorch 代碼相對(duì)較長(zhǎng)，但可以針對(duì)其他內(nèi)核函數(shù)或我們將在以后的筆記本中討論的更高級(jí)功能（例如可擴(kuò)展推理或用于分類的非高斯似然）進(jìn)行簡(jiǎn)單修改。

18.3.7。練習(xí)

我們強(qiáng)調(diào)了學(xué)習(xí)內(nèi)核超參數(shù)的重要性，以及超參數(shù)和內(nèi)核對(duì)高斯過程泛化特性的影響。嘗試跳過我們學(xué)習(xí) hypers 的步驟，而是猜測(cè)各種長(zhǎng)度尺度和噪聲方差，并檢查它們對(duì)預(yù)測(cè)的影響。當(dāng)你使用大長(zhǎng)度尺度時(shí)會(huì)發(fā)生什么？一個(gè)小的長(zhǎng)度尺度？噪聲方差大？噪聲方差?。?/p>

我們已經(jīng)說過邊際似然不是凸目標(biāo)，但是可以在 GP 回歸中可靠地估計(jì)長(zhǎng)度尺度和噪聲方差等超參數(shù)。這通常是正確的——事實(shí)上，邊際似然在學(xué)習(xí)長(zhǎng)度尺度超參數(shù)方面比空間統(tǒng)計(jì)中的傳統(tǒng)方法要好得多，后者涉及擬合經(jīng)驗(yàn)自相關(guān)函數(shù)（“協(xié)方差圖”）?？梢哉f，至少在最近關(guān)于可擴(kuò)展推理的工作之前，機(jī)器學(xué)習(xí)對(duì)高斯過程研究的最大貢獻(xiàn)是為超參數(shù)學(xué)習(xí)引入了邊際似然。

然而，即使這些參數(shù)的不同配對(duì)也為許多數(shù)據(jù)集提供了可解釋的不同合理解釋，從而導(dǎo)致我們目標(biāo)中的局部最優(yōu)。如果我們使用大長(zhǎng)度尺度，那么我們假設(shè)真正的基礎(chǔ)函數(shù)正在緩慢變化。如果觀察到的數(shù)據(jù)變化很大，那么我們唯一可以合理地?fù)碛写箝L(zhǎng)度尺度的就是大的噪聲方差。另一方面，如果我們使用較小的長(zhǎng)度尺度，我們的擬合將對(duì)數(shù)據(jù)的變化非常敏感，幾乎沒有空間來解釋噪聲變化（任意不確定性）。

嘗試看看您是否可以找到這些局部最優(yōu)值：使用具有大噪聲的非常大的長(zhǎng)度尺度和具有小噪聲的小長(zhǎng)度尺度進(jìn)行初始化。你會(huì)收斂到不同的解決方案嗎？

我們已經(jīng)說過，貝葉斯方法的一個(gè)基本優(yōu)勢(shì)是自然地表示認(rèn)知不確定性。在上面的例子中，我們無法完全看到認(rèn)知不確定性的影響。嘗試使用 進(jìn)行預(yù)測(cè)。當(dāng)您的預(yù)測(cè)超出數(shù)據(jù)范圍時(shí)，95% 的可信集會(huì)發(fā)生什么變化？它是否涵蓋了該區(qū)間內(nèi)的真實(shí)函數(shù)？如果你只想象那個(gè)區(qū)域的任意不確定性，會(huì)發(fā)生什么？test_x = np.linspace(0, 10, 1000)

嘗試運(yùn)行上述示例，但使用 10,000、20,000 和 40,000 個(gè)訓(xùn)練點(diǎn)，并測(cè)量運(yùn)行時(shí)間。培訓(xùn)時(shí)間如何安排？或者，運(yùn)行時(shí)間如何隨測(cè)試點(diǎn)的數(shù)量變化？預(yù)測(cè)均值和預(yù)測(cè)方差是否不同？通過從理論上計(jì)算出訓(xùn)練和測(cè)試的時(shí)間復(fù)雜度，以及使用不同數(shù)量的點(diǎn)運(yùn)行上面的代碼來回答這個(gè)問題。

嘗試使用不同的協(xié)方差函數(shù)運(yùn)行 GPyTorch 示例，例如 Matern 內(nèi)核。結(jié)果如何變化？在 GPyTorch 庫(kù)中找到的光譜混合內(nèi)核怎么樣？有些比其他更容易訓(xùn)練邊際可能性嗎？一些對(duì)于長(zhǎng)期預(yù)測(cè)和短期預(yù)測(cè)更有價(jià)值嗎？

在我們的 GPyTorch 示例中，我們繪制了包括觀察噪聲在內(nèi)的預(yù)測(cè)分布，而在我們的“從頭開始”示例中，我們僅包括認(rèn)知不確定性。重做 GPyTorch 示例，但這次只繪制認(rèn)知不確定性，并與從頭開始的結(jié)果進(jìn)行比較。預(yù)測(cè)分布現(xiàn)在看起來一樣嗎？（他們應(yīng)該。）

Discussions