【摘 要】高斯过程作为一种用于预测的非参数模型,可以用于回归任务,也可以用于分类任务,本文主要介绍其在回归任务中的主要原理和方法。《机器学习中的高斯过程》一书是高斯过程研究领域的扛鼎之作,本文主要节选自该书的第二章。

【原 文】 Rasmussen, C.E. and Williams, C.K. (2006), Chapter 2 of Gaussian processes for machine learning. Cambridge, Mass: MIT press Cambridge, MA (3).

【提 醒】 本文所有内容均是在假设协方差(核)函数已知的情况下进行的讨论。因此,无论是从权重视角还是从函数视角,关注的主要是(权重或函数的)先验以及(权重或函数的)后验推断。

第2章 高斯过程回归

监督学习可以分为回归和分类问题。分类的输出是离散的类标签,而回归与连续量的预测有关。例如,在金融应用程序中,人们可能会尝试根据利率、货币汇率、可用性和需求来预测商品价格。在本章中,我们描述了回归问题的高斯过程方法;分类问题在第 3 章讨论

有多种方法可以解释高斯过程 (GP) 回归模型。可以将高斯过程视为定义函数的分布,并直接在函数空间(函数空间视角)中进行推断。尽管这种视角很吸引人,但最初可能难以理解,因此:

  • 2.1 节 中,从许多人熟悉和易于理解的权重空间视角开始阐述
  • 2.2 节 中继续介绍函数空间视角
  • 高斯过程通常具有可以通过设置某些参数来改变的特性,在 2.3 节 中将讨论这些参数的变化对高斯过程特性的影响。
  • 高斯过程模型的预测会给出完整的预测分布形式,但实际工作中,有时需要点估计,在 2.4 节 中,我们将讨论如何使用决策理论将损失函数与预测分布结合,进而给出最佳点预测。
  • 2.5 节 中介绍了一个涉及学习机器人手臂的逆向动力学的实际比较示例。
  • 2.6 节中给出了关于高斯过程回归平滑性质的一些理论分析。
  • 2.7 节 中讨论了如何将显式均值函数组合到模型中,并进行相应的推断和预测。
  • 本章中大部分材料是标准的,我们将在 第 2.8 节 中介绍历史概述。

2.1 权重空间视角

在已被广泛研究和使用的简单线性回归模型中,输出是输入的线性组合。其主要优点是实施简单和可解释性;缺点是灵活性有限;如果输入和输出之间的关系不能用线性函数合理地近似,则该模型会给出较差的预测。

在本节中,我们先讨论线性模型的贝叶斯处理。然后,将输入投影到高维特征空间,并在该空间中应用线性模型来实现对简单线性模型的增强。事实表明:在某些特殊的空间中,可以应用 “核技巧” 在原始输入空间中实现高维空间的隐式计算;在特征空间维度较大时,这一点能够有效节省计算量。

在介绍具体方法之前,我们先对场景进行形式化的描述:

我们有一个包含 nn 个观测值的训练集 D\mathcal{D}D={(xi,yi)i=1,,n}\mathcal{D} = \{(\mathbf{x}_i, y_i) | i = 1, \ldots, n\},其中 x\mathbf{x} 表示维度为 DD输入向量yy 表示标量输出或目标。所有 nn 个案例的输入列向量可以聚合在一个 D×nD × n设计矩阵 XX 中,而目标则可以收集在向量 y\mathbf{y} 中,因此我们可以写成向量形式 D=(X,y)D = (X, \mathbf{y})。在回归任务中,目标是实数值。我们的兴趣是推断输入和目标之间的关系,即在给定输入的情况下目标的条件分布(注:在回归任务中,我们通常对输入的分布并不感兴趣)。

2.1.1 简单线性模型的贝叶斯分析

(1)模型定义

首先回顾含高斯噪声的标准线性回归模型:

f(x)=xw,y=f(x)+ε(2.1)f(\mathbf{x}) = \mathbf{x}^{\top} \mathbf{w}, \qquad y = f(\mathbf{x}) + \varepsilon \tag{2.1}

其中 x\mathbf{x} 是输入向量,w\mathbf{w} 是线性模型的权重(参数)向量,ff 是函数值,yy 是观测到的目标值。通常模型中还包含偏置(或偏移量),但这可以通过使用值始终为 11 的附加元素增广输入向量 x\mathbf{x} 来实现,因此没有明确表示。我们假设观测值 yy 与函数值 f(x)f(\mathbf{x}) 的差异在于加性噪声,并进一步假设该噪声服从均值为零、方差为 σn2\sigma^2_n 的独立同分布高斯分布。

εN(0,σn2)(2.2)\varepsilon \sim \mathcal{N}(0, \sigma^2_n) \tag{2.2}

(2)似然函数

这种噪声假设与模型一起产生了似然,即给定参数时观测值的概率密度,并且可以被分解为训练集中单个案例的似然之积(或对数似然之和)(由于独立性假设):

p(yX,w)=i=1np(yixi,w)=i=1n12πσnexp((yixiw)22σn2)=1(2πσn2)n/2exp(12σn2yXw2)=N(Xw,σn2I)(2.3)\begin{align*} p(\mathbf{y}|X, \mathbf{w}) &= \prod^{n}_{i=1} p(y_i|\mathbf{x}_i, \mathbf{w}) = \prod^{n}_{i=1} \frac{1}{\sqrt{2π} \sigma_n} \exp (− \frac{(y_i − \mathbf{x}^{\top}_i \mathbf{w})^2}{2 \sigma^2_n} ) \\ &= \frac{1}{(2π \sigma^2_n )^{n/2}} \exp ( − \frac{1}{2 \sigma^2_n} |\mathbf{y} − X^{\top} \mathbf{w}|^2) = \mathcal{N}(X^{\top} \mathbf{w}, \sigma^2_n I ) \end{align*} \tag{2.3}

其中符号 |\cdot| 表示向量的欧几里得长度。

(3)先验分布

在贝叶斯形式中,需要在参数上指定一个先验,用于表达我们在看到观测之前对参数的信念。我们在模型的权重参数上放置一个协方差矩阵为 Σp\Sigma_p 的零均值高斯先验。

wN(0,Σp)(2.4)\mathbf{w} \sim \mathcal{N}(\mathbf{0}, \Sigma_p) \tag{2.4}

此先验的作用和性质将在 第 2.2 节中讨论;现在我们按照指定先验继续推导。

(4)后验分布

贝叶斯线性模型的推断基于权重参数的后验分布,根据贝叶斯规则计算:

posterior=likelihood×priormarginal likelihood,p(wy,X)=p(yX,w)p(w)p(yX)(2.5)\text{posterior} = \frac{\text{likelihood} \times \text{prior}}{ \text{marginal likelihood}}, \qquad p(\mathbf{\mathbf{w}}|\mathbf{y},X) = \frac{p(\mathbf{y}|X, \mathbf{w}) p(\mathbf{w})}{p(\mathbf{y}|X)} \tag{2.5}

(5)边缘似然

上式中 p(yX)p(\mathbf{y}|X) 为归一化常数 ,也称为边缘似然(或证据),是一个与权重参数无关的量(权重被边缘化掉了),由下式给出:

p(yX)=p(yX,w)p(w)dw(2.6)p(\mathbf{y}|X) = \int p(\mathbf{y}|X, \mathbf{w}) p(\mathbf{w}) d \mathbf{w} \tag{2.6}

(6)后验的解析形式

式 (2.5) 中的后验结合了似然和先验,并捕获了关于参数的一切。如果只考虑贝叶斯公式右侧的分子项目(即似然项和先验项),结合高斯分布,我们可以发现权重参数的后验也是一个高斯分布:

p(wX,y)p(yX,w)p(w)exp(12σn2(yXw)(yXw))exp(12wΣp1w)exp(12(wwˉ)(1σn2XX+Σp1)(wwˉ)) where wˉ=σn2(σn2XX+Σp1)1Xy.(2.7)\begin{align*} p(\mathbf{\mathbf{w}}|X, \mathbf{y}) &\propto p(\mathbf{y}| X, \mathbf{w})p(\mathbf{w}) \\ &\propto \exp (− \frac{1}{ 2\sigma^2_n} (\mathbf{y} − X^{\top} \mathbf{\mathbf{w}})^{\top} (\mathbf{y} − X^{\top} \mathbf{\mathbf{w}})) \exp (− \frac{1}{2} \mathbf{\mathbf{w}}^{\top} \boldsymbol{\Sigma}^{-1}_p \mathbf{\mathbf{w}}) \\ &\propto \exp \left(− \frac{1}{2} ( \mathbf{\mathbf{w}} −\bar{\mathbf{\mathbf{w}}})^{\top} ( \frac{1}{\sigma^2_n}X X^{\top} + \boldsymbol{\Sigma}^{-1}_p )(\mathbf{\mathbf{w}} −\bar{\mathbf{\mathbf{w}}}) \right)\\ \\ &\text{ where } \bar{\mathbf{\mathbf{w}}} = \sigma^{-2}_n (\sigma^{-2}_n XX^{\top} + \boldsymbol{\Sigma}^{-1}_p )^{−1} X \mathbf{y}. \end{align*} \tag{2.7}

通过进一步推导(略),可以得到:权重参数的后验是一个均值为 wˉ\bar{\mathbf{\mathbf{w}}} 、协方差为 A1A^{-1} 的高斯分布( AA 也被称为 精度矩阵):

p(wX,y)N(wˉ,A1)wˉ=1σn2A1XyA=σn2XX+Σp1\begin{align*} p(\mathbf{\mathbf{w}}|X, \mathbf{y}) &\sim \mathcal{N}(\bar{\mathbf{w}}, A^{-1}) \tag{2.8} \\ \bar{\mathbf{w}} &= \frac{1}{\sigma^2_n} A^{-1}X \mathbf{y}\\ A &= \sigma^{-2}_n XX^{\top} + \boldsymbol{\Sigma}^{-1}_p \end{align*}

请注意:由于高斯分布的均值等于峰值,因此 wˉ\bar{\mathbf{w}} 也被称为 的最大后验 (MAP) 估计。

注: 非贝叶斯情况

在非贝叶斯设置中,先验(或其负对数形式)有时被认为是惩罚项,此时作出的 MAP 估计也被称为 w\mathbf{w} 的带惩罚最大似然估计。但请注意,在贝叶斯设置中, MAP 估计仅指后验概率为峰值时的选择,并没有其他特殊含义。本例中的带惩罚最大似然估计方法也被称为 “ 岭回归 ” [Hoerl 和 Kennard,1970],其惩罚来自于 式(2.7) 中的二次惩罚项 12wΣp1w\frac{1}{2} \mathbf{w}^{\top} \boldsymbol{\Sigma}^{-1}_p \mathbf{w}

(7)预测分布

为了对测试案例进行预测,贝叶斯方法会依据参数的后验概率对预测做出加权平均(注:这与非贝叶斯方法形成明显不同,在非贝叶斯方案中,通常由某些准则确定参数的单一点估计值)。也就是说,在贝叶斯方法中,测试点 x\mathbf{x}_* 处的预测结果 ff(x)f_* \triangleq f(\mathbf{x}_*) 不是一个点估计,而是一个关于后验分布进行积分(术语称之为边缘化)得到的、有关所有可能输出值的概率分布,我们称之为预测分布(Predictive Distribution):

p(fx,X,y)=p(fx,w)p(wX,y)dw=xwp(wX,y)dw=N(μ,Σ)μ=1σn2xA1XyΣ=xA1x\begin{align*} p(f_*|\mathbf{x}_*, X, \mathbf{y}) &= \int p(f_*|\mathbf{x}_*, \mathbf{w})p(\mathbf{w}|X, \mathbf{y}) d\mathbf{w} \\ &= \int \mathbf{x}_*^{\top}\mathbf{w} \, p(\mathbf{w}|X, \mathbf{y})d \mathbf{w} \\ &= \mathcal{N}(\mu_*,\Sigma_*) \tag{2.9}\\ \\ \mu_* &= \frac{1}{\sigma^2_n} \mathbf{x}_*^{\top} A^{-1} X \mathbf{y}\\ \Sigma_* &= \mathbf{x}_*^{\top} A^{-1} \mathbf{x}_* \end{align*}

显然,再次由于高斯分布,预测分布仍然具有解析形式,并且也是一个高斯分布,其中

  • 预测均值:由测试输入 x\mathbf{x}_* 乘以 式 (2.8) 中模型参数的后验均值给出,由于是线性模型,因此预测均值可以被视为训练目标 (y\mathbf{y}) 的线性组合,所以也常被称为 “最佳线性预测”。
  • 预测方差:是测试输入 x\mathbf{x}_* 关于后验协方差矩阵 A1A^{-1} 的二次型,这表明预测不确定性会随着测试输入的位置变化而变化,符合人们的预期。

(8)图示

图 2.1 给出了贝叶斯线性回归的一个例子。我们在此选择了一维输入空间,因此权重空间是二维的,主要是为了更容易地可视化。高斯先验的等值线显示在子图 (a) 中。数据在子图 (b) 中描绘为十字。这产生了子图 ( c ) 中显示的似然和子图 (d) 中的后验分布。预测分布及其误差条也在子图 (b) 中作出了标记。

Fig2.01

图 2.1:贝叶斯线性模型示例 f(x)=w1+w2xf(x) = w_1 + w_2 x,其中截距为 w1w_1,斜率为 w2w_2。 (a) 显示了先验分布 p(w)N(0,I)p(\mathbf{w}) \sim \mathcal{N}(\mathbf{0, I}) 的等值线,见式(2.4)。(b) 显示了三个用十字标记的训练数据点。图 ( c ) 显示了似然 p(yX,w)p(\mathbf{y} \mid X, \mathbf{w}) 的等高线,见式 (2.3),假设噪声等级为 σn=1σ_n = 1,因此,斜率会比截距更为 “确定”。 (d) 显示了后验 p(wX,y)p(\mathbf{w} |X, \mathbf{y}), 见式 (2.7);比较后验概率的最大值,可以看到截距更接近零了,而更 “确定” 的斜率几乎没有变化。所有等值线图都给出了 11 倍和 22 倍标准差的等概率等值线。叠加在 (b) 中数据上的是(无噪声)预测分布 p(fx,X,y)p(f_*|\mathbf{x}_*, X, \mathbf{y}) 的预测均值加/减两倍标准差(虚线),见式 (2.9)。

2.1.2 广义线性模型与核函数

(1)基本想法

在上一节中,我们回顾了贝叶斯简单线性模型,可以看出,该模型的表达能力有限。

为了提高表达能力,一个简单想法是用一组基函数将输入投影到某个高维空间中,然后在这个高维度空间中使用线性模型。 例如,标量型输入 x\mathbf{x} 可以被投影到 xx 的幂空间:ϕ(x)=(1,x,x2,x3,)\boldsymbol{\phi}(x) = (1, x, x^2,x^3,\ldots)^{\top} 中,从而实现了一个多项式回归。

只要投影运算独立于模型参数 w\mathbf{w},那么模型关于参数 w\mathbf{w} 仍然是线性的,因此在分析上也易于处理。

在本章,我们将暂时回避基函数的选择问题,假设基函数已经给定。不过,在 第 5 章 中我们将证明:高斯过程能够给出选择基函数的方法。

(2)模型

具体来说,我们引入了一组函数 ϕ(x)\boldsymbol{\phi}(\mathbf{x}),将 DD 维的输入向量 x\mathbf{x} 映射到 NN 维的特征空间。如果令矩阵 Φ(X)\Phi(X) 表示训练集中所有样例经 ϕ(x)\boldsymbol{\phi}(\mathbf{x}) 转换后的列向量构成的聚合矩阵。则新模型可以改写为基函数的线性形式:

f(x)=ϕ(x)w(2.10)f(\mathbf{x}) = \boldsymbol{\phi}(\mathbf{x})^{\top} \mathbf{\mathbf{w}} \tag{2.10}

式(2.1) 的不同之处在于:特征空间中的输入维度 NN (通常 N>DN>D )以及对应的参数向量 w\mathbf{w}。由于在特征空间中该模型仍然是一个简单线性模型,因此我们可以直接使用上一节的分析结果,在推断时,只需要用映射后的矩阵 Φ(X)\Phi(X) 代替原输入矩阵 XX 即可;而在预测时,只需要用映射后的 ϕ(x)\phi(\mathbf{x}_*) 代替了原始预测输入 x\mathbf{x}_*。根据 式(2.9),可以很容易得到如下预测分布:

fx,X,yN(μ,Σ)μ=1σn2ϕ(x)A1ΦyΣ=ϕ(x)A1ϕ(x)\begin{align*} f_*|\mathbf{x}_*, X, \mathbf{y} &\sim \mathcal{N}(\mu_*,\Sigma_*)\\ \\ \mu_* &= \frac{1}{\sigma^2_n} \boldsymbol{\phi} (\mathbf{x}_*)^{\top} A^{-1} \Phi \mathbf{y}\\ \Sigma_* &= \boldsymbol{\phi}(\mathbf{x}_*)^{\top} A^{-1} \boldsymbol{\phi}(\mathbf{x}_*) \tag{2.11} \end{align*}

式中的 Φ=Φ(X)\Phi = \Phi(X)A=σn2ΦΦ+Σp1A = \sigma^{-2}_n \Phi \Phi^{\top} + \boldsymbol{\Sigma}^{-1}_p

我们很快会发现,使用此方程进行预测,需要对大小为 N×NN × NAA 矩阵求逆,如果特征空间的维数 NN 很大就会很麻烦(但大多数时候,我们希望能将原始空间张成更高维度的特征空间,以便发现模式)。不过,我们可以用以下方式重写方程:

fx,X,yN(μ,Σ)μ=ϕΣpΦ(K+σn2I)1yΣ=ϕΣpϕϕΣpΦ(K+σn2I)1ΦΣpϕ\begin{align*} f_*|\mathbf{x}_*, X, \mathbf{y} &\sim \mathcal{N}(\mu_*,\Sigma_*)\\ \\ \mu_* &= \boldsymbol{\phi}^{\top}_* \Sigma_p \Phi(K + \sigma^2_n I )^{−1} \mathbf{y}\\ \Sigma_* &= \boldsymbol{\phi}^{\top}_* \Sigma_p \boldsymbol{\phi}_* − \boldsymbol{\phi}^{\top}_* \Sigma_p \Phi(K + \sigma^2_n I )^{−1} \Phi^{\top} \Sigma_p \boldsymbol{\phi}_* \tag{2.12} \end{align*}

其中使用了简写 ϕ=ϕ(x)\boldsymbol{\phi}_* = \boldsymbol{\phi}(\mathbf{x}_*) 并定义了一个新的符号 K=ΦΣpΦK = \Phi^{\top} \Sigma_p \Phi

下面对公式做一简单推导:

  • 对于均值:根据 AAKK 的定义,有 σn2Φ(K+σn2I)=σn2Φ(ΦΣpΦ+σn2I)=AΣpΦ\sigma^{-2}_n \Phi (K + \sigma^2_n I) = \sigma^{-2}_n \Phi(\Phi^{\top} \Sigma_p \Phi + \sigma^2_n I) = A \Sigma_p \Phi。现在从左边乘以 A1A^{-1} 并从右边乘以 (K+σn2I)1(K + \sigma^2_nI)^{-1},即可以得到 σn2A1Φ=ΣpΦ(K+σn2I)1\sigma^{-2}_n A^{-1}\Phi = \Sigma_p \Phi(K + \sigma^2_nI)^{-1},表明式(2.11)式 (2.12)中均值表达式的等价性。
  • 对于方差:可以直接使用 式 (A.9) 的矩阵求逆引理得到,只需要设 Z1=Σp2Z^{-1} = \Sigma^2_pW1=σn2IW^{-1} = \sigma^2_nIV=U=ΦV = U = \Phi

需要注意的是:在式 (2.12) 中,我们需要对大小为 n×nn × n 的矩阵求逆,当 n<Nn < N 时这更为方便。在几何上,请注意 nn 个数据点最多可以张成特征空间中的 nn 个维度。

(3) 核与核技巧

仔细观察 式 (2.12),我们会发现,原始输入总是以 ΦΣpΦ\Phi^{\top} \Sigma_p \Phi, ϕΣpΦ\boldsymbol{\phi}^{\top}_* \Sigma_p \Phi, 或 ϕΣpϕ\boldsymbol{\phi}^{\top}_* \Sigma_p \boldsymbol{\phi}_* 的方式进入特征空间,也就是说,在特征空间中的运算总是具有 ϕ(x)Σpϕ(x)\boldsymbol{\phi}(\mathbf{x})^{\top} \Sigma_p \boldsymbol{\phi}(\mathbf{x'}) 的形式,其中 x\mathbf{x}x\mathbf{x'} 要么在训练集中,要么在测试集中。因此,我们会思考:是否能够直接用原始空间中的计算来替换特征空间中的计算呢?如果能够这样的话,理论上在特征空间的模式分析任务,都可以转换成原始空间中的直接计算,我们就可以不用再考虑 NN 的大小问题了。

事实表明,这样做是可行的,但正如上面的分析一样,需要具备一个条件,即: 特征空间中的模式分析任务必须是基于内积的计算

此时,我们可以定义一个函数 k(x,x)=ϕ(x)Σpϕ(x)k(\mathbf{x}, \mathbf{x'}) = \boldsymbol{\phi}(\mathbf{x})^{\top} \Sigma_p \boldsymbol{\phi}(\mathbf{x'}) 来表示这种形式,此函数就是我们熟知的 协方差函数

请注意,ϕ(x)Σpϕ(x)\boldsymbol{\phi}(\mathbf{x})^{\top} \Sigma_p \boldsymbol{\phi}(\mathbf{x}') 本质上是一个相对于 Σp\Sigma_p内积。下面是简单证明:

假设 Σp\Sigma_p 是正定矩阵,则我们可以定义 Σp1/2\Sigma^{1/2}_p 使得 (Σp1/2)2=Σp(\Sigma^{1/2}_p)^2 = \Sigma_p (注: 当 Σp\Sigma_p 的奇异值分解(SVD)为 Σp=UDU\Sigma_p = UDU^{\top},其中 DD 是对角矩阵,则 Σp1/2\Sigma^{1/2}_p 可以写成 UD1/2UU D^{1/2}U^{\top} 的形式)。进一步地,我们定义 ψ(x)=Σp1/2ϕ(x)\boldsymbol{\psi}(\mathbf{x}) = \Sigma^{1/2}_p \boldsymbol{\phi}(\mathbf{x}),则可以得到一个简单的点积表示 k(x,x)=ψ(x)ψ(x)k(\mathbf{x}, \mathbf{x}') = \boldsymbol{\psi}(\mathbf{x}) \cdot \boldsymbol{\psi}(\mathbf{x}')。也就是说,只要确保 Σp\Sigma_p 是正定矩阵,就可以仅根据输入空间中的某种计算来定义和计算特征空间中的内积,这种情况也被称为 核技巧。核技巧使 “核” 成为大家关注的主要对象,特征空间反而因此被放在了次要位置。

2.2 函数空间视角

2.2.1 函数的分布

通过直接在函数空间中考虑推断,可以获得与上一节相同结果的等效方法。为此,我们引入 高斯过程 (GP) 来描述函数的分布:

【定义 2.1】 高斯过程是随机变量的集合,该集合中任意有限数量的随机变量服从联合高斯分布。

与权重视角下用 权重的均值向量权重的协方差矩阵 来定义模型的方式不同,高斯过程由两个函数完全指定:一个是 均值函数 ,另外一个是 协方差函数 ,也称 。我们将真实过程 f(x)f(\mathbf{x}) 的均值函数 m(x)m(\mathbf{x}) 和协方差函数 k(x,x)k(\mathbf{x}, \mathbf{x}') 定义为:

m(x)=E[f(x)]k(x,x)=E[(f(x)m(x))(f(x)m(x))](2.13)\begin{align*} m(\mathbf{x}) &= \mathbb{E}[f (\mathbf{x})]\\ k(\mathbf{x}, \mathbf{x}') &= \mathbb{E}[(f(\mathbf{x}) − m(\mathbf{x}))(f (\mathbf{x}') − m(\mathbf{x}'))] \end{align*} \tag{2.13}

并将高斯过程写成

f(x)GP(m(x),k(x,x))(2.14)f(\mathbf{x}) \sim \mathcal{GP}(m(\mathbf{x}), k(\mathbf{x}, \mathbf{x}')) \tag{2.14}

通常,为了表示简单起见,我们将均值函数取为零(尽管不需要这样做,参见第 2.7 节)。

在上面的例子中,随机变量代表了位置 x\mathbf{x}f(x)f(\mathbf{x}) 的函数值。高斯过程的常见应用是随时间定义的,即随机变量的索引集是时间。不过,在更高维度空间中使用高斯过程也很常见,因此,我们在这里假设索引集 X\mathcal{X} 是所有可能的输入构成的集合,例如 XRd\mathcal{X} \subset \mathbb{R}^d

为了符号方便,我们使用训练集中案例的(任意)枚举来识别随机变量,即 fif(xi)f_i \triangleq f(\mathbf{x}_i) 是对应于案例 (xi,yi)(\mathbf{x}_i, y_i) 的随机变量。

根据定义,高斯过程是随机变量的集合。该定义自动隐含了 一致性要求(或称 边缘化性质): 如果指定了一个高斯过程,如 (y1,y2)N(μ,Σ)(\mathbf{y}_1, \mathbf{y}_2) \sim \mathcal{N}(\boldsymbol{μ}, \Sigma) ,那么它同时也指定了边缘分布 y1N(μ1,Σ11)\mathbf{y }_1 \sim \mathcal{N}(\boldsymbol{μ}_1, \Sigma_{11}),其中 Σ11\Sigma_{11}Σ\Sigma 的子矩阵,参见式(A.6)。请注意: 如果协方差矩阵的元素由某个协方差函数指定,则一致性要求会自动得到满足。上述定义并不排除具有有限索引集的高斯过程,只是这种情况仅仅是一个高斯分布而已。

2.2.2 均值函数 – 再看线性回归模型

对于上一节中的贝叶斯线性回归模型:f(x)=ϕ(x)wf(\mathbf{x}) = \boldsymbol{\phi} (\mathbf{x})^{\top}\mathbf{w},从权重视角来看,其权重参数的先验为 wN(0,Σp)\mathbf{w} \sim \mathcal{N}(\mathbf{0}, \Sigma_p),模型输出的均值和协方差为:

E[f(x)]=ϕ(x)E[w]=0E[f(x)f(x)]=ϕ(x)E[ww]ϕ(x)=ϕ(x)Σpϕ(x)(2.15)\begin{align*} \mathbb{E}[f (\mathbf{x})] &= \boldsymbol{\phi}(\mathbf{x})^{\top} \mathbb{E}[\mathbf{w}] = 0\\ \mathbb{E}[f (\mathbf{x})f(\mathbf{x'})] &= \boldsymbol{\phi} (\mathbf{x})^{\top}\mathbb{E}[\mathbf{w}\mathbf{w}^{\top}]\boldsymbol{\phi} (\mathbf{x'}) = \boldsymbol{\phi} (\mathbf{x})^{\top}\Sigma_p\boldsymbol{\phi} (\mathbf{x'}) \end{align*} \tag{2.15}

也就是说,f(x)f(\mathbf{x})f(x)f (\mathbf{x'}) 是联合的高斯分布,其均值为零,协方差由 ϕ(x)Σpϕ(x)\boldsymbol{\phi} (\mathbf{x})^{\top}\Sigma_p\boldsymbol{\phi} (\mathbf{x'}) 给出。实际上,对于高斯过程来说,任意数量 nn 输入点对应的函数值 f(x1),,f(xn)f(\mathbf{x}_1),\ldots , f (\mathbf{x}_n),都服从联合高斯分布。尽管当 N<nN < n 时,此联合高斯分布是奇异的(注:NN 是特征空间的维度,nn 为案例的数量;此处的奇异性,是因为在特征空间中的联合协方差矩阵的秩为 NN )。

2.2.3 协方差函数 – 再看协方差矩阵

协方差函数指定了两个随机变量之间的协方差,在本章中,协方差函数将主要以 平方指数 (SE) 协方差函数 式(2.16))为例;其他协方差函数将在第 4 章 中讨论。

cov(f(xp),f(xq))=k(xp,xq)=exp(12xpxq2)(2.16)\operatorname{cov}(f(\mathbf{x}_p), f(\mathbf{x}_q)) = k(\mathbf{x}_p, \mathbf{x}_q) = \exp (-\frac{1}{2} |\mathbf{x}_p − \mathbf{x}_q|^2) \tag{2.16}

这里有几点需要注意:

(1)输出之间的协方差被定义为输入的函数。对于这个特定的协方差函数,我们可以看到:输入之间非常接近的两个变量之间的协方差几乎是 11,随着两者在输入空间中距离的增加,协方差不断减小。

(2)平方指数协方差函数对应于具有无限多个基函数的贝叶斯线性回归模型 (详细论证见 第 4.3.1 节)。事实上,所有正定的协方差函数 k(,)k(·,·),都存在相应的基函数展开形式(可能是无限的,参见 第 4.3 节 中的 Mercer 定理)。

(3)指定了协方差函数就意味着确定了函数的分布。为了看到这一点,我们可以在任意数量的输入处,从函数的分布中抽取样本;具体来说,我们选择一些输入点 XX_*,然后用 式(2.16) 填充相应协方差矩阵中的元素,并用这个协方差矩阵生成一个随机高斯向量:

fN(0,K(X,X))(2.17)\mathbf{f}_* \sim \mathcal{N}(\mathbf{0}, K(X_*, X_*)) \tag{2.17}

进一步将生成的高斯向量绘制成输入的函数形式。图 2.2(a) 显示了三个这样的示例。

Fig2.2

图 2.2: (a) 显示了从高斯过程先验中随机抽取的三个函数;点表示实际生成的 yy 值(其中一个函数);通过将大量评估点连接起来,就形成了函数的线状形式(其他两个函数)。 (b) 显示了从后验中抽取的三个随机函数,即五个无噪声观测为条件的先验。在这两个图中,阴影区域分别代表先验和后验的逐点平均值加上和减去两倍标准差(对应于 95%95\% 置信区域)。

图 2.2 的示例中,输入值是等距的。请注意,函数看起来很平滑。事实上,平方指数协方差函数是无限可微的,导致该高斯过程是无限均方可微的(见 4.1 节)。我们还看到函数似乎具有特征长度尺度(即变程超参数),请参阅第 4.2.1 节。对于式 (2.16),特征长度尺度在一个单位左右。如果用 xpxq/|\mathbf{x}_p −\mathbf{x}_q|/\ell 替换 式(2.16) 中的 xpxq|\mathbf{x}_p −\mathbf{x}_q|\ell 为某个正常数),我们就可以改变该过程的特征长度尺度。此外,随机函数的总体方差可以通过 式(2.16) 中指数函数之前的正值前置因子来控制 。我们将在 2.3 节 中更多地讨论这些因素如何影响预测,并在 第 5 章 中更多地讨论如何设置这些尺度参数。

2.2.4 函数空间视角下的用无噪声观测做预测

我们通常对从先验中抽取的具体随机函数并不感兴趣,而是希望能够结合训练数据,得到关于函数的知识(如:某个输入处最可能的输出值,不确定性等)。我们先考虑无噪声观测的简单情形,即知道训练点处的真实过程值 {(xi,fi)i=1,n}\{(\mathbf{x}_i, f_i)| i = 1,\ldots ,n\}(当然,此情形是现实世界中不可能发生的)。此时,训练输出 f\mathbf{f} 和测试输出 f\mathbf{f}_* 的联合先验分布是:

[ff]N(0,[K(X,X)K(X,X)K(X,X)K(X,X)])(2.18)\left[\begin{array}{l} \mathbf{f} \\ \mathbf{f}_{*} \end{array}\right] \sim \mathcal{N}\left(\mathbf{0},\left[\begin{array}{ll} K(X, X) & K\left(X, X_{*}\right) \\ K\left(X_{*}, X\right) & K\left(X_{*}, X_{*}\right) \end{array}\right]\right) \tag{2.18}

如果有 nn 个训练点和 nn_* 个测试点,则 K(X,X)K(X, X_*) 表示在所有训练点和测试点对处计算的协方差的 n×nn × n_* 矩阵,其他项 K(X,X)K(X, X ), K(X,X)K(X_*, X_*)K(X,X)K(X_*, X) 也类似。

为了获得函数的后验分布,需要限制该联合先验只能包含那些与已观测数据点一致的函数。参考 图 2.2 ,最简单的方法可能是从先验抽取很多函数样本,然后拒绝那些与观测结果不一致的函数样本(这种方法看起来很蠢,但很能说明问题)。 幸运的是,从概率角度来看并不需要使用如此笨的方法,因为,训练点的预测直接对应于以观测为条件的先验(更多详细信息请参见第 A.2 节):

fX,X,fN(K(X,X)K(X,X)1f,K(X,X)K(X,X)K(X,X)1K(X,X))(2.19)\mathbf{f}_*|X_*, X, \mathbf{f} \sim \mathcal{N}(K(X_*,X) K(X, X)^{-1} \mathbf{f}, \quad K(X_*, X_*) − K(X_*, X)K(X, X)^{−1} K(X, X_*)) \tag{2.19}

或等价的:

fˉ=K(X,X)K(X,X)1fcov(f)=K(X,X)K(X,X)K(X,X)1K(X,X)\begin{align*} \bar{\mathbf{f}}_* &= K(X_*,X) K(X, X)^{-1} \mathbf{f}\\ \operatorname{cov}(\mathbf{f}_*) &= K(X_*, X_*) − K(X_*, X) K(X, X)^{−1} K(X, X_*) \end{align*}

注: 这里其实使用了三个机制:

  • 一是以观测数据为条件的高斯过程先验就是高斯过程后验,即如果高斯过程的先验为 p(f,f)p(\mathbf{f,f_*}),则其后验为 p(f,f{(xi,fi)})p(\mathbf{f,f_*}|\{(\mathbf{x}_i, f_i)\})
  • 二是后验也是联合高斯分布,根据边缘性质,其条件分布依然是高斯分布 p(ff)N(μ,Σ)p(\mathbf{f}_*|\mathbf{f}) \sim \mathcal{N}(\mu,\Sigma)
  • 三是上述所有高斯分布,在高斯过程设置中都受协方差函数 k(,)k(\cdot,\cdot) 控制,而协方差函数与输入 X,XX,X_* 有关,因此上述公式中均隐含了对 X,XX,X_* 的依赖。

函数值 f\mathbf{f}_*(对应于测试输入 XX_*)可以通过采样方式获得,其分布对应于 式 (2.19) 中均值和协方差定义的联合高斯分布(注:采样方法可参考 第 A.2 节 )。

图 2.2(b) 显示了这种计算的结果,图中给出了标有 ‘+’ 符号的五个训练数据点。请注意,将这种计算扩展到多维输入非常简单,只需要根据 式(2.16) 改变协方差函数的计算(不过生成的函数可能更难于用图形方式显示)。

2.2.5 函数空间视角下的用含噪声观测做预测

现实场景中,我们无法访问真实函数值本身,而只能访问其含噪声版本 y=f(x)+εy = f (\mathbf{x}) + \varepsilon。假设加性独立同分布高斯噪声 ε\varepsilon 具有方差 σn2\sigma^2_n,含噪声观测的先验变成

cov(yp,yq)=k(xp,xq)+σn2δpq or cov(y)=K(X,X)+σn2I(2.20)\operatorname{cov}(y_p, y_q) = k(\mathbf{x}_p, \mathbf{x}_q) + \sigma^2_n \delta_{pq} \quad \text{ or } \quad \operatorname{cov}(\mathbf{y}) = K(X, X) + \sigma^2_n I \tag{2.20}

其中 δpq\delta_{pq} 是 Kronecker delta,当且仅当 p=qp = q11,否则为 00,代表了噪声的独立性假设。与式 (2.16)的无噪声情况相比,上式添加了对角矩阵。进一步的,在式(2.18) 基础上引入噪声项,可以将已观测目标值和测试位置处函数值的联合先验分布写成:

[yf]N(0,[K(X,X)+σn2IK(X,X)K(X,X)K(X,X)])(2.21)\left[\begin{array}{l} \mathbf{y} \\ \mathbf{f}_{*} \end{array}\right] \sim \mathcal{N}\left(\mathbf{0},\left[\begin{array}{ll} K(X, X) + \sigma^2_n I & K\left(X, X_{*}\right) \\ K(X_{*}, X) & K(X_{*}, X_{*}) \end{array}\right]\right) \tag{2.21}

对应于 式 (2.19) 的条件分布。我们得出高斯过程回归的关键预测方程:

fX,y,XN(fˉ,cov(f)) where fˉE[fX,y,X]=K(X,X)[K(X,X)+σn2I]1ycov(f)=K(X,X)K(X,X)[K(X,X)+σn2I]1K(X,X)\begin{align*} \mathbf{f}_* | X, \mathbf{y}, X_* &\sim \mathcal{N}(\bar{\mathbf{f}}_*, \operatorname{cov}(\mathbf{f}_*)) \text{ where } \tag{2.22}\\ \bar{\mathbf{f}}_* &\triangleq \mathbb{E}[\mathbf{f}_*|X, \mathbf{y}, X_*] = K(X_*, X)[K(X, X) + \sigma^2_n I]^{−1} \mathbf{y} \tag{2.23}\\ \operatorname{cov}(\mathbf{f}_*) &= K(X_*, X_*) − K(X_*, X)[K(X, X) + \sigma^2_n I]^{−1} K(X, X_*) \tag{2.24} \end{align*}

2.2.6 函数-权重空间视角的等价性

请注意: 当明确 K(C,D)=Φ(C)ΣpΦ(D)K(C, D) = \Phi(C)^{\top} \Sigma_p \Phi(D) 时( 其中 C,DC, D 代表 XXXX_*),式(2.22)式(2.24) 与权重视角下的 式 (2.12) 完全对应。这意味着:给定任意一组基函数,我们都可以将其对应的协方差函数写为 k(xp,xq)=ϕ(xp)Σpϕ(xq)k(\mathbf{x}_p, \mathbf{x}_q) = \boldsymbol{\phi} (\mathbf{x}_p)^{\top}\Sigma_p\boldsymbol{\phi} (\mathbf{x}_q);反之,如果给定(正定的)协方差函数 kk,我们也可以导出某种基函数展开形式(可能是无限的、不唯一的),请参见第 4.3 节

K(X,X)K(X, X), K(X,X)K(X, X_*)K(X,X)K(X_*, X_*) 等表达形式看起来相当恼火,因此从现在起,我们引入一种紧凑的符号表示: K=K(X,X)K = K(X, X)K=K(X,X)K_* = K(X, X_*)。在只有一个测试点 x\mathbf{x}_* 的情况下,我们用 k(x)=k\mathbf{k}(\mathbf{x}_*) = \mathbf{k}_* 来表示测试点和 nn 个训练点之间的协方差向量。使用这种紧凑符号,对于单个测试点 x\mathbf{x}_*式(2.23)式(2.24) 可以简写为:

fˉ=k(K+σn2I)1yV[f]=k(x,x)k(K+σn2I)1k\begin{align*} \bar{f}_* &= \mathbf{k}_*^{\top} (K + \sigma^2_n I)^{-1} \mathbf{y} \tag{2.25}\\ \mathbb{V}[f_*] &= k(\mathbf{x}_*, \mathbf{x}_*) − \mathbf{k}_*^{\top} (K + \sigma^2_n I)^{-1} \mathbf{k}_* \tag{2.26} \end{align*}

仔细观察 式(2.25)式(2.26) 会发现: 式 (2.25) 中的预测均值是已观测目标 y\mathbf{y} 的线性组合;这也是它被称为线性预测器的原因。也可以采用另一个角度来看 式 (2.25),将预测均值看成是 nn 个训练案例处核函数值的线性组合:

fˉ(x)=i=1nαik(xi,x)(2.27)\bar{f}(\mathbf{x}_*) = \sum^{n}_{i=1} α_i k(\mathbf{x}_i, \mathbf{x}_*) \tag{2.27}

其中,权重为 α=(K+σn2I)1y\boldsymbol{α} = (K + \sigma^2_n I)^{-1} \mathbf{y}

式 (2.27) 中关于预测均值的这种形式是 表征定理(representer theorem) 的一种表现;有关这一点的更多信息,请参见 第 6.2 节

我们可以直观地理解这个结果。首先,尽管高斯过程被定义为所有 yy 变量的联合高斯分布(索引集 X\mathcal{X} 中的所有位置都对应一个,因此是无穷的),但为了在 x\mathbf{x}_* 处进行预测,我们只用关心由 nn 个训练点和 11 个测试点定义的 (n+1)(n+1) 维分布即可。其次,多元高斯分布的特性使得其边缘分布(或条件分布),都可以被解析地写成另外一个高斯形式,而且只需要处理联合协方差矩阵中相关的块即可(参阅第 A.2 节)。因此,以观测为条件的 (n+1)(n + 1) 维分布能够为我们提供所需的结果。图 2.3 给出了含噪声观测高斯过程的图模型表示。

Fig2.3

图 2.3:用于回归的高斯过程的图模型(链图)。正方形代表观测到的变量,圆圈代表未知量。粗水平条表示一组完全连接的节点。请注意,给定隐变量 fif_i 时,观测 yiy_i 条件独立于所有其他节点。由于高斯过程的边缘化特性,添加更多输入 xx、隐变量 ff 和未观测的目标 yy_* 并不会改变任何其他变量的分布。

另请特别注意:

预测方差不依赖于已观测目标,而只依赖于测试输入、训练输入和超参数(见 式(2.24) );这是高斯过程的一个性质。方差代表了两项之间的差异:第一项 K(X,X)K(X_*, X_*) 是先验协方差;从中减去一个代表了数据给我们的信念的(正)项。通过向隐变量值添加方差的形式,我们可以很容易地得到含噪声测试目标 y\mathbf{y}_* 的预测分布,即将 σn2I\sigma^2_n \mathbf{I} 添加到 cov(f)\operatorname{cov}(f_*) 表达式的方差中。

事实上,高斯过程模型的预测分布不仅仅给出了 式 (2.26) 的逐点误差,当 XX_* 表示多个测试输入时,式(2.24) 依然成立;此时能够同时计算测试目标的协方差(其对角矩阵元素仍是逐点方差)。事实上,式 (2.23) 是均值函数,式(2.24) 是后验高斯过程的协方差函数。后验协方差如 图 2.4(b) 所示。

Fig2.4

图 2.4:(a) 与图 2.2(b) 相同,显示了从后验抽取的三个随机函数。 (b) 显示了 f(x)f(\mathbf{x})f(x)f(\mathbf{x'}) 之间的后验协方差,用于相同数据的三个不同 x\mathbf{x'} 值。请注意,封闭点处的协方差很高,在训练点下降到零(因为是无噪声过程,所以此处方差为零),然后变为负值…。发生这种情况是:如果平滑函数在数据点的一侧小于均值,它会另一侧倾向于超过均值,导致数据点处协方差符号的反转。请注意对比,先验协方差只是高斯形状并且从不为负。

2.2.7 边缘似然

现在引入重要的概念:边缘似然(或证据)p(yX)p(\mathbf{y}|X),因为它在模型选择和超参数优化方面意义重大,经常被用作确定超参数的目标函数(见 第 5 章 )。

边缘似然是似然乘以先验的积分:

p(yX)=p(yf,X)p(fX)df(2.28)p(\mathbf{y}|X) = \int p(y|f , X)p(f |X) df \tag{2.28}

可以看出,在函数空间视角下,边缘似然主要指函数值 f\mathbf{f} 的边缘化。

(1)先验

首先看一下先验。在高斯过程模型下,先验是高斯的,fXN(0,K)\mathbf{f} |X \sim \mathcal{N}(\mathbf{0}, K),其具体的对数形式为:

logp(fX)=12fK1f12logKn2log2π(2.29)\log p(\mathbf{f} |X) = -\frac{1}{2} \mathbf{f}^{\top} K^{−1} \mathbf{f} -\frac{1}{2} \log |K| − \frac{n}{2} \log 2π \tag{2.29}

(2)似然

高斯过程的似然也是高斯的 yfN(f,σn2I)\mathbf{y}|\mathbf{f} \sim \mathcal{N}(\mathbf{f}, \sigma^2_n I),而且可以被分解为各训练案例似然的组合。

(3)边缘似然

利用附录的 式 A.7式 A.8 可以完成对先验的积分,并产生如下对数边缘似然形式:

logp(yX)=12y(K+σn2I)1y12logK+σn2In2log2π(2.30)\log p(\mathbf{y}|X) = -\frac{1}{2} y^{\top}(K + \sigma^2_n I)^{-1} \mathbf{y} -\frac{1}{2} \log |K + \sigma^2_n I| − \frac{n}{2} \log 2π \tag{2.30}

该边缘似然可以直接通过训练观测 y\mathbf{y} 获得(注:根据前面的假设,yN(0,K+σn2I)\mathbf{y} \sim \mathcal{N}(\mathbf{0}, K + \sigma^2_n I ) )。

2.2.8 预测算法

当给定协方差函数的形式和超参数时,算法 2.1 给出了高斯过程回归 (GPR) 预测的实际实现,该算法的输出中不仅包含预测均值和预测方差,还包含边缘似然值。该算法使用 Cholesky 分解实现矩阵的求逆运算,因为其具有快速和数值稳定的特点,参见附录 第 A.4 节。该算法返回无噪声测试案例的预测均值和预测方差,要计算含噪声测试数据 yy_* 的预测分布,只需将噪声方差 σn2\sigma^2_n 添加到 ff_* 的预测方差中即可。

Alg2.1

算法 2.1:高斯过程回归的预测及其对数边缘似然。该算法使用 Cholesky 分解实现了 式(2.25)式 (2.26) 所需的矩阵求逆,参见 A.4 节。对于多个测试案例,只需重复第 4-6 行即可。式 (2.30) 中所需的对数行列式是根据 Cholesky 分解计算的(对于大 nn,可能无法表示行列式本身,但可以给出值)。第 2 行中 Cholesky 分解的计算复杂度为 n3/6n^3/6,第 3 行和第 5 行中(对于每个测试案例)三角矩阵求解的计算复杂度为 n2/2n^2/2

2.3 直观地认识超参数

注:

获得最佳的协方差函数形式和超参数值,是高斯过程推断的主要任务,也是高斯过程应用 pipeline 中最核心环节。这里面可能需要理清几个术语:

高斯过程训练: 获得最佳协方差函数形式和超参数值的过程;

高斯过程推断: 根据观测数据过得高斯过程后验的过程;

高斯过程预测: 根据推断的后验结果,获得测试点处预测均值和预测方差的过程。

通常,协方差函数会有一些自由参数。例如,一维的平方指数协方差函数具有以下形式

ky(xp,xq)=σf2exp(122(xpxq)2)+σn2δpq(2.31)k_y(x_p, x_q) = \sigma^2_f \exp (-\frac{1}{2\ell^2} (x_p − x_q)^2) + \sigma^2_n \delta_{pq} \tag{2.31}

协方差表示为 kyk_y,因为它是针对含噪声目标 yy 而不是底层真实函数 ff 的。式中的 长度尺度 \ell信号方差 σf2\sigma^2_f噪声方差 σn2\sigma^2_n 都是可以变化的自由参数,通常我们将其称为超参数。

第 5 章中,我们将考虑从训练数据中确定超参数的方法。但在本节中,我们的目标比较简单,只是直观地探索一下改变超参数会对高斯过程预测产生哪些影响。

考虑图 2.5(a)+ 号所示的数据。这是从具有平方指数核的高斯过程生成的,其中 (,σf,σn)=(1,1,0.1)(\ell, σ_f , \sigma_n) = (1, 1, 0.1)。图中还显示了使用这些超参数值获得的 22 倍标准差预测条。请注意,对于远离任何训练点的输入值,误差会变大。事实上,如果 xx 轴被延长,人们会发现,随着远离数据,标准差逐步偏离了过程方差( σfσ_f )。

Fig2.5

图 2.5:(a) 数据是从具有超参数 (,σf,σn)=(1,1,0.1)(\ell, σ_f , σ_n) = (1, 1, 0.1) 的高斯过程生成的,如 + 符号所示。使用具有上述超参数值的高斯过程预测,可以获得真实函数 ff95%95\% 置信区域(以灰色显示)。 (b) 和 ( c ) 则分别针对了超参数值为 (0.31.080.00005)(0.3、1.08、0.00005)(3.01.160.89)(3.0、1.16、0.89) 的情况。

如果我们将长度尺度设置得更短,使 =0.3\ell = 0.3 并保持其他参数不变,那么我们将从该过程中生成类似于 图 2.5(a) 的图,不过 xx 轴应当缩放 0.30.3 倍;等效地,如果保持与 图 2.5(a) 中相同的 xx 轴,那么函数样本看起来会更加摇摆不定。如果我们对从 =1\ell = 1 过程生成的数据使用 =0.3\ell = 0.3 的过程进行预测,将获得 图 2.5(b) 中的结果。其余两个参数将像第 5 章所述那样,通过优化边缘似然的方法设置。在这种情况下,噪声参数减小到 σn=0.00005\sigma_n = 0.00005,因为 “信号” 的更大灵活性意味着可以降低噪声水平。这可以在图中的 x=2.5x = 2.5 附近的两个数据点处观测到。在图 2.5(a) (=1\ell = 1) 中,这些基本上被解释为具有不同噪声的相似函数值。然而,在 图 2.5(b) (=0.3\ell = 0.3) 中,噪声水平非常低,因此这两点必须用真实函数值 ff 的剧烈变化来解释。另请注意,较短的长度尺度意味着 图 2.5(b) 中的误差条在远离数据点时会快速增长。

相反,如果将长度尺度设置得更长一些,比如设置为 =3\ell=3,如 图2.5(c) 所示。同样,剩余的两个参数通过优化边缘似然来设置。此时噪声水平已增加到 σn=0.89\sigma_n = 0.89,我们看到数据现在由具有大量噪声的缓慢变化的函数来解释。

当然,我们可以将快速变化的低噪声信号或缓慢变化的高噪声信号的位置取到极端;前者会产生信号的白噪声过程模型,而后者会产生带有附加白噪声的恒定信号。在这两种模型下,产生的数据点应该看起来像白噪声。但是,通过研究 图 2.5(a),我们发现白噪声并不是一个令人信服的数据模型,因为 yy 的序列不会足够快地交替,底层真实函数的可变性会为其带来相关性。

上述分析在此处的一维案例中相对容易看到,但高维情况会比较复杂。 第 5 章 中即将讨论的边缘似然等方法,可以推广到更高维度,并允许我们对各种模型进行量化评价。例如,根据边缘似然,上述案例会明显倾向于采用 (,σf,σn)=(1,1,0.1)(\ell, σ_f , \sigma_n) = (1, 1, 0.1) 而不是其他备选方案。

2.4 高斯过程回归的决策理论

(1)决策问题的提出

高斯过程在测试点的预测输出表现为一个高斯形式的预测分布,而人们做决策时需要的有时是确定的单点估计(可以简单理解为给出 “是” 或 “否” 的确切答案,即便存在不确定性),此时就存在一个决策问题。 例如在前面部分中,测试输入 x\mathbf{x}_* 的输出 yy_∗ 就对应于一个预测分布。该分布也是高斯的,均值和方差由 式(2.25)式 (2.26) 给出。但在实际应用中,我们经常被迫需要作出一个决定,即一个在某种意义上最优的点估计。

为了提供决策支持,通常会选择一个准则或损失函数 L(ytrue,yguess)\mathcal{L}(y_{true}, y_{guess}),用于定义在真实值为 ytruey_{true} 时猜测值 yguessy_{guess} 带来的损失(或惩罚)。例如,损失函数可以等于猜测与事实之间的绝对偏差。

(2)贝叶斯与非贝叶斯决策的比较

需要注意的是: 高斯过程是在没有参考损失函数的情况下计算得出的预测分布,这与非贝叶斯方法存在显著不同。 在非贝叶斯方法中,通常会通过经验风险(或损失)最小化来训练模型(可以理解为训练和决策一体);而在贝叶斯方法中,似然函数(主要用于训练)和损失函数(主要用于决策)之间存在明显的区分

  • 似然函数用于描述噪声测量值偏离(假设的)无噪声真实函数值的程度;
  • 损失函数则用于捕获作出特定选择后的后果(在给定实际真实状态的情况下);
  • 似然函数和损失函数之间理论上不需要有任何共同点,尽管有时两者之间会存在联系。

(3)高斯过程的决策理论

决策的目标是让点预测 yguessy_{guess} 损失最小,但是当我们不知道 ytruey_{true} 真实值时该怎么决定损失大小呢?一个自然的想法是采用期望损失(或风险)最小化,即关于模型的所有可能选项进行平均,并将其结果视为真值:

R~L(yguessx)=L(y,yguess)p(yx,D)dy(2.32)\tilde{R}_{\mathcal{L}}(y_{guess} | \mathbf{x}_*) = \int \mathcal{L}(y_*, y_{guess}) p(y_*|\mathbf{x}_*, \mathcal{D}) d y_* \tag{2.32}

因此,从期望损失最小化的意义上说,最佳猜测应当是:

yoptimalx=argminyguessR~L(yguessx)(2.33)y_{optimal}|\mathbf{x}_* = \arg \min_{y_{guess}} \tilde{R}_{\mathcal{L}}(y_{guess}|\mathbf{x}_*) \tag{2.33}

一般来说,使期望风险 yguessy|y_{guess}− y_*| 最小化的 yguessy_{guess} 值会是预测分布 p(yx,D)p(y_*|\mathbf{x}_*, D) 的中值,但对于平方损失 (yguessy)2(y_{guess} − y_*)^2 而言,最优值会是该分布的均值。特别是:当预测分布是高斯分布时,均值和中值重合,此时的最优点估计就是均值。实际上,对于任何呈对称状的损失函数和预测分布,都会得出预测分布的均值是最佳 yguessy_{guess} 的结果。不过,在许多实际问题中,损失函数可能不对称,例如在严格安全的应用中,点预测可以直接从 式(2.32)式 (2.33) 计算 。在 Berger [1985] 中可以找到对决策理论的综述。

决策理论可以说是贝叶斯方法特有的一种方法,因为在非贝叶斯方法中只有点估计的概念,没有预测分布的概念,因此不存在决策问题。本节给出了如下知识点:

  • 强调了贝叶斯方法和非贝叶斯方法在训练和决策上的不同,即训练和决策是否分离;
  • 给出了利用期望损失最小化来解决决策问题的一个基本框架( 式(2.32)式 (2.33) );
  • 给出了一个常用的决策选择,即对于大多数场景,可以采用所以可能值的均值作为最优点估计。

2.5 示例应用

暂略。

2.6 高斯过程的平滑性质分析

高斯过程回归旨在通过去除污染噪声 ε\varepsilon 来重建底层真实信号 ff。为此,它计算含噪声观测 yy 的加权平均值 fˉ(x)=k(K+σn2I)1y\bar{f} (\mathbf{x}_*) = \mathbf{k}_*^{\top}(K +\sigma^2_nI)^{-1} \mathbf{y}。由于 fˉ(x)\bar{f} (\mathbf{x}_*)y\mathbf{y} 的线性组合,所以高斯过程回归是一种线性平滑器(有关详细信息,请参见 Hastie 和 Tibshirani [1990,第 2.8 节])。

在本节中,我们首先通过训练点处预测结果的矩阵分析来研究平滑性质,然后再通过等效核研究。

(1)训练点处的平滑性质

训练点的预测均值 fˉ\bar{\mathbf{f}} 由下式给出

fˉ=K(K+σn2I)1y(2.35)\bar{\mathbf{f}}= K(K + \sigma^2_n I)^{-1} \mathbf{y} \tag{2.35}

KK 具有特征分解 K=i=1nλiuiuiK = \sum^{n}_{i=1} \lambda_i \mathbf{u}_i \mathbf{u}^{\top}_i ,其中 λi\lambda_i 是第 ii 个特征值,ui\mathbf{u}_i 是对应的特征向量。由于 KK 是实数且对称半正定,其特征值是实数且非负,并且其特征向量相互正交。对于系数 γi=uiy\gamma_i = \mathbf{u}^{\top}_i \mathbf{y},有 y=i=1nγiui\mathbf{y} = \sum^{n}_{i=1} \gamma_i \mathbf{u}_i 。则预测均值可以被分解为:

fˉ=i=1nγiλiλi+σn2ui(2.36)\bar{\mathbf{f}} = \sum^{n}_{i=1} \frac{\gamma_i \lambda_i }{\lambda_i + \sigma^2_n} \mathbf{u}_i \tag{2.36}

请注意,如果 λi/(λi+σn2)1\lambda_i /(\lambda_i + \sigma^2_n) \ll 1,则 y\mathbf{y} 中沿 ui\mathbf{u}_i 的分量将被显著消去。

对于实践中使用的大多数协方差函数,变化较慢的特征向量(例如,较少的过零点)会有更大的特征值,这意味着 y\mathbf{y} 中的高频分量会被平滑掉。 有效参数的数量(或平滑器的自由度)被定义为 tr(K(K+σn2I)1)=i=1nλi/(λi+σn2)\operatorname{tr}(K(K + \sigma^2_n I)^{-1}) = \sum^{n}_{i=1} \lambda_i /(\lambda_i + \sigma^2_n),参见 Hastie 和 Tibshirani [1990,sec. 3.5]。请注意,这代表了未消去特征向量的数量。

注: 这一段没有太理解想解释什么。

(2)权重函数与等效核

我们可以定义一个函数的向量 h(x)=(K+σn2I)1k(x)\mathbf{h}(\mathbf{x}_*) = (K + \sigma^2_n I)^{-1} \mathbf{k}(\mathbf{x}_*)。因此有 fˉ(x)=h(x)y\bar{f}(\mathbf{x}_*) = \mathbf{h}(\mathbf{x}_*)^{\top} \mathbf{y},这清楚地表明: x\mathbf{x}_* 处的均值预测是目标值 y\mathbf{y} 的线性组合。对于固定的测试点 x\mathbf{x}_*h(x)\mathbf{h}(\mathbf{x}_*) 给出了应用于目标 y\mathbf{y} 的权重向量。 h(x)\mathbf{h}(\mathbf{x}_*) 被称为权重函数 [Silverman, 1984]。由于高斯过程回归是线性平滑器,因此权重函数不依赖于 y\mathbf{y}。请注意线性模型(预测是输入的线性组合)和线性平滑器(预测是训练集目标的线性组合)之间的区别。

K+σn2IK+\sigma^2_n I 矩阵的求逆以及 KK 取决于 nn 个数据点的具体位置,这使我们对权重函数形式的理解更为复杂。理想情况下,我们可以认为观测结果在 x\mathbf{x} 空间中以某种观测密度被 “抹掉” 了(注:此处翻译不是特别满意)。这样的话,我们可以使用如 第 7.1 节 所示的分析方法来解决此问题(注: 第 7.1 节 重点介绍了等效核这一理论工具)。Silverman [1984] 将权重函数与 核平滑 做类比,将理想化的权重函数称为 等效核;另见 Girosi 等 [1995 年,section 2.1]。

核平滑器将一个核函数 κκ 放置在 x\mathbf{x}_* 处,然后为每个数据点 (xi,yi)(\mathbf{x}_i, y_i) 计算 κi=κ(xix/)κ_i = κ(|\mathbf{x}_i − \mathbf{x}_*|/\ell),其中 \ell 为长度尺度。高斯是常用的核函数之一。 f(x)f(\mathbf{x}_*) 的预测被计算为 f^(x)=i=1nwiyi\hat{f}(\mathbf{x}_*) = \sum^{n}_{i=1} w_i y_i 其中 wi=κi/j=1nκjw_i = κ_i/ \sum^{n}_{j=1} κ_j。这也称为 Nadaraya-Watson 估计,参见例如 Scott [1992,section 8.1]。

对于一维输入变量 xx图 2.6 说明了高斯过程的权重函数和等效核。我们使用了平方指数协方差函数,并设置了长度尺度 =0.0632\ell = 0.0632(因此 2=0.004\ell^2 = 0.004)。沿 xx 轴随机分布了 n=50n = 50 个训练点。图 2.6(a)图 2.6(b) 分别显示了超参数 σn2=0.1\sigma^2_n = 0.1x=0.5x_* = 0.5x=0.05x_* = 0.05 的权重函数和等效核。图 2.6(c) 也用于 x=0.5x_* = 0.5,但超参数 σn2=10\sigma^2_n = 10。在所有情况下,点标记都对应于权重函数 h(x)\mathbf{h}(x_*),实线标记都对应于等效核。虚线标记则显示以测试点为中心的平方指数核,已经被缩放为与等效核的最大值相同的高度。图 2.6(d) 显示了等效核作为 nn 的函数时的变化情况,nn 是单位区间中的数据点数。

Fig2.6

图 2.6: (a)-( c ) 显示了对应于 n=50n = 50 个训练点的权重函数 h(x)\mathbf{h}(x_*)(点)、等效核(实线)和原始平方指数核(虚线)。 (d) 显示了两种不同数据密度的等效核。详情参见文本。测试点处的小十字用于在所有四个图中进行缩放。

从这些图中可以得出许多有趣的观察结果:

  • 首先观察到: 等效核具有与原始平方指数核完全不同的形状。在 图 2.6(a) 中,等效核显然是振荡的(具有负瓣)并且具有比原始核更高的空间频率。图 2.6(b) 显示了类似行为,不过相对于 图 2.6(a) 来说,由于存在边效应,等效核被截断了。在 图 2.6(c) 中,我们看到在较高的噪声水平下,负瓣减少了,等效核的宽度与原始核相似。
  • 此外,与 图 2.6(a)图 2.6(b) 相比,图 2.6(c) 中等效核的高度有所降低(它在更宽区域上取平均值)。

我们可以用上面的特征分解来分析和理解较低噪声水平时更加振荡的等效核:在更高的噪声水平下,仅保留了 y\mathbf{y} 的大 λλ(缓慢变化)分量;而对于更小的噪声水平,会保留更多的振荡分量。

图 2.6(d) 中,我们绘制了 [0,1][0, 1]n=10n = 10n=250n = 250 个数据点的等效核;注意等效核的宽度如何随着 nn 的增加而减小。我们将在 第 7.1 节 中进一步讨论这种行为。

图 2.6 中的等效核图是通过在 [0,1][0, 1] 上使用 ngridn_{grid} 点的密集网格然后计算更平滑的矩阵 K(K+σgrid2I)1K(K + σ^2_{grid} I)^{-1} 绘制的。该矩阵的每一行都是适当位置的等效核,因为这是对单位向量 y\mathbf{y} 的响应,该向量在除一点以外的所有点都为零。但是,为了获得正确的缩放比例,必须设置 σgrid2=σn2ngrid/nσ^2_{grid} = \sigma^2_n n_{grid}/n;对于 ngrid>nn_{grid}>n 这意味着每个 ngridn_grid 点的有效方差更大,但由于相应的点数更多,这种影响被抵消了。这可以通过想象如果在单个 xx 位置有方差为 σ2σ^2 网格的 ngrid/nn_{grid}/n 独立高斯观测的情况来理解;这相当于方差为 σn2\sigma^2_n 的一个高斯观测值。实际上,nn 个观测值已沿区间均匀平滑。如果我们达到连续极限并寻求平滑噪声函数,则可以解析地获得等效核的形式。第 7.1 节 给出了相关理论和一些示例等效核。

2.7 与显式均值函数的结合

考虑具有零均值函数的高斯过程很常见,但绝不是必要的,因为后验高斯过程的均值并不限于为零。人们可能希望对均值函数进行显式建模有几个原因,如:模型的可解释性、表达先验信息的便利性、为利于分析而做的必要约束等。使用显式基础函数是一种在函数上指定非零均值的方法,但还可以使用它们来实现其他有趣的效果。

2.7.1 固定的均值函数

使用固定的均值函数 m(x)m(\mathbf{x}) 是非常简单的:只需将零均值高斯过程应用于观测值和固定均值函数之间的差即可。即

f(x)GP(m(x),k(x,x))(2.37)f (\mathbf{x}) \sim \mathcal{GP}(m(\mathbf{x}), k(\mathbf{x, x'})) \tag{2.37}

此时的预测均值变为:

fˉ=m(X)+K(X,X)Ky1(ym(X))(2.38)\bar{\mathbf{f}}_* = \mathbf{m}(X_*) + K(X_*, X)K^{−1}_y (\mathbf{y} − \mathbf{m}(X)) \tag{2.38}

其中 Ky=K+σn2IK_y = K + \sigma^2_nI

预测方差由于和目标值无关,因此与式 (2.24) 保持一致。

2.7.2 不固定的均值函数

在实际工作中,通常很难指定一个固定的均值函数。在许多情况下,指定某一类由系数 β\boldsymbol{\beta} 确定的均值函数可能更方便。此时,需要将该系数从数据中推断出来。例如,考虑一个均值函数为线性的模型:

g(x)=f(x)+h(x)β, 其中 f(x)GP(0,k(x,x))(2.39)g(\mathbf{x}) = f(\mathbf{x}) + \mathbf{h}(\mathbf{x})^{\top} \boldsymbol{\beta},\text{ 其中 } f(\mathbf{x}) \sim \mathcal{GP}(0, k(\mathbf{x, x'})) \tag{2.39}

这里 f(x)f(\mathbf{x}) 是零均值高斯过程,h(x)\mathbf{h}(\mathbf{x}) 是一组固定的基函数,β\boldsymbol{\beta} 则是新增加的参数。该公式表示:数据接近于全局线性模型,其残差由高斯过程建模。早在 1975 年,Blight 和 Ott [1975] 就明确探索了这个想法,他们使用高斯过程对多项式回归的残差进行建模,即 h(x)=(1,x,x2,)\mathbf{h}(x) = (1, x, x^2, \ldots)

(1)推断

在拟合模型时,可以结合协方差函数的超参数对参数 β\boldsymbol{\beta} 进行优化。或者,也将 β\boldsymbol{\beta} 视为贝叶斯的,例如将 β\boldsymbol{\beta} 的先验取为高斯分布,βN(b,B)\boldsymbol{\beta} \sim \mathcal{N}(\mathbf{b}, B),并且采用贝叶斯方法的边缘化来完成后续任务(如利用对参数后验的边缘化,计算目标的预测分布)。

根据 O’Hagan [1978],我们事实上可以获得另一个高斯过程:

g(x)GP(h(x)b,k(x,x)+h(x)Bh(x)(2.40)g(\mathbf{x}) \sim \mathcal{GP}(\mathbf{h}(\mathbf{x})^{\top} \mathbf{b}, k(\mathbf{x, x'}) + \mathbf{h}(\mathbf{x})^{\top} B \mathbf{h}(\mathbf{x'}) \tag{2.40}

该高斯过程中的协方差函数中,增加了由参数的不确定性带来的贡献。

(2)预测

g(x)g(\mathbf{x}) 的均值函数和协方差函数代入 式(2.39)式 (2.24) 进行预测。重新排列后,可以得到:

gˉ(X)=Hβˉ+KKy1(yHβˉ)=fˉ(X)+Rβˉcov(g)=cov(f)+R(B1+HKy1H)1R(2.41)\begin{align*} \bar{\mathbf{g}}(X_*) &= H^{\top}_* \bar{\boldsymbol{\beta}} + K_* K^{-1}_y (\mathbf{y} − H \bar{\boldsymbol{\beta}}) =\bar{\mathbf{f}}(X_*) + R^{\top} \bar{\boldsymbol{\beta} }\\ \operatorname{cov}(\mathbf{g}_*) &= \operatorname{cov}(\mathbf{f}_*) + R^{\top}(B^{-1} + H K^{-1}_y H^{\top})^{-1} R \end{align*} \tag{2.41}

其中 HH 矩阵收集所有训练案例的 h(x)\mathbf{h}(\mathbf{x}) 向量(类似的 HH_* 表示所有测试),βˉ=(B1+HKy1H)1(HKy1y+B1b)\bar{\boldsymbol{\beta}} = (B^{-1} + HK^{-1}_y H^{\top})^{-1}(HK^{-1}_y \mathbf{y} + B^{-1}\mathbf{b}) , 和 R=HHKy1KR = H_* − HK^{-1}_y K_*

注意 式 (2.41) 第一行中均值表达式的漂亮解释,:βˉ\bar{\boldsymbol{\beta}} 是全局线性模型参数的均值,是数据项和先验项之间的折衷,而预测均值只是平均线性输出再加上高斯过程模型从残差中预测的值。协方差是通常的协方差项和新的非负贡献之和。

探索上述表达式的极限,随着 β\boldsymbol{\beta} 参数的先验变得模糊,B1OB^{-1} \rightarrow O(其中 OO 是零矩阵),我们将得到一个独立于 b\mathbf{b} 的预测分布

g(X)=fˉ(X)+Rβˉcov(g)=cov(f)+R(HKy1H)1R(2.42)\begin{align*} \mathbf{g}(X_*) &= \bar{\mathbf{f}}(X_*) + R^{\top} \bar{\boldsymbol{\beta}}\\ \operatorname{cov}(\mathbf{g}_*) &= \operatorname{cov}(\mathbf{f}_*) + R^{\top} (HK^{-1}_y H^{\top})^{-1} R \end{align*} \tag{2.42}

其中极限 βˉ=(HKy1H)1HKy1y\bar{\boldsymbol{\beta}} = (H K^{-1}_y H^{\top})^{-1} H K^{-1}_y \mathbf{y}

请注意,在极限 B1OB^{-1} \rightarrow O 下的预测,不应该通过将 式(2.40) 中修改后的协方差函数插入标准预测方程来实现。因为协方差函数的元素趋于无穷大,不适合数值实现。相反,式 (2.42) 必须被使用。即使对非极限情况感兴趣,式 (2.41) 在数值上也优于 式(2.40) 的直接实现 ,由于全局线性部分往往会在协方差矩阵中加入一些非常大的特征值,因此会影响其条件数。

(3)边缘似然

在这个简短的部分中,我们简要讨论 式(2.40) 中模型(均值函数的参数具有先验 βN(b,B)\boldsymbol{\beta} \sim \mathcal{N}(b, B) )的边缘似然,因为这在后面的章节会有用,特别是 第 6.3.1 节。我们可以利用 式(2.30) 来表示边缘似然:

logp(yX,b,B)=12(Hby)(Ky+HBH)1(Hby)12logKy+HBHn2log2π(2.43)\log p(\mathbf{y}|X, \mathbf{b}, B) = − \frac{1}{2} (H^{\top} \mathbf{b} − \mathbf{y})^{\top} (K_y + H^{\top} B H)^{-1}(H^{\top} \mathbf{b} − \mathbf{y}) − \frac{1}{2} \log |K_y + H^{\top} B H| − \frac{n}{2} \log 2π \tag{2.43}

其中已经包含了显式的均值。我们感兴趣的是当 B1OB^{-1} \rightarrow O 的极限情况下(即当先验模糊时)的表现。在此极限下,先验的均值是无关的(如 式(2.42)中的情况),因此在不失一般性的情况下(对于极限情况)我们假设均值为零( b=0\mathbf{b} = 0 ),从而给出:

logp(yX,b=0,B)=12yKy1y+12yCy12logKy12logB12logAn2log2π(2.44)\log p(\mathbf{y}|X, \mathbf{b = 0}, B) = − \frac{1}{2} \mathbf{y}^{\top} K^{-1}_y \mathbf{y} + \frac{1}{2} \mathbf{y}^{\top} C \mathbf{y} − \frac{1}{2} \log |K_y| − \frac{1}{2} \log |B| − \frac{1}{2} \log |A| − \frac{n}{2} \log 2π \tag{2.44}

其中 A=B1+HKy1HA = B^{-1} + H K^{-1}_y H^{\top}C=Ky1HA1HKy1C = K^{-1}_y H^{\top} A^{-1} H K^{-1}_y, 并且使用了 式(A.9)式 (A.10) 的矩阵求逆引理。

我们现在探索 β\boldsymbol{\beta} 上模糊先验极限下的对数边缘似然的行为。在此极限下,HH^{\top} 的列跨越的方向上的高斯方差将变得无穷大,显然这需要特殊处理。对数边缘似然由三项组成:y\mathbf{y} 的二次型、对数行列式项和涉及 log2π\log 2π 的项。执行协方差矩阵的特征分解,我们看到无限方差方向对二次型项的贡献为零。但对数行列式项倾向于负无穷大。在这种情况下,标准解 [Wahba,1985 年;Ansley 和 Kohn,1985 年] 是将 y\mathbf{y} 投影到与 HH^{\top} 张成的正交方向上,并计算该子空间中的边缘似然。如果 HH^{\top} 的秩为 mm,如 Ansley 和 Kohn [1985] 所示,则意味着我们必须丢弃 式 (2.44) 中的 12logBm2log2π− \frac{1}{2} \log |B| − \frac{m}{2} \log 2π 项,得到:

logp(yX)=12yKy1y+12yCy12logKy12logAnm2log2π(2.45)\log p(\mathbf{y}|X) = -\frac{1}{2} \mathbf{y}^{\top} K^{-1}_y \mathbf{y} + \frac{1}{2} \mathbf{y}^{\top} C \mathbf{y} -\frac{1}{2} \log |K_y| -\frac{1}{2} \log |A| − \frac{n−m}{2} \log 2π \tag{2.45}

其中 A=HKy1HA = H K^{-1}_y H^{\top}C=Ky1HA1HKy1C = K^{-1}_y H^{\top} A^{-1} H K^{-1}_y

2.8 历史及相关工作

使用高斯过程进行预测当然不是最近的话题,尤其是对于时间序列分析;基本理论至少可以追溯到 1940 年代 Wiener [1949] 和 Kolmogorov [1941] 的工作。事实上,Lauritzen [1981] 讨论了丹麦天文学家 T. N. Thiele 从 1880 年开始的相关工作。

高斯过程预测在地统计学领域(参见 Matheron,1973 年;Journel 和 Huijbregts,1978 年)也广为人知,在那里它被称为克里金法,在气象学领域 [Thompson,1956 年,Daley,1991 年],尽管该文献自然主要关注二维和三维输入空间。Whittle [1963,section 5.4] 也建议使用这种方法进行空间预测。 Ripley [1981] 和 Cressie [1993] 对空间统计中的高斯过程预测提供了有用的概述。

逐渐意识到高斯过程预测可以用于一般的回归环境。例如,O’Hagan [1978] 提出了方程 2.23 和 2.24 中给出的一般理论,并将其应用于许多一维回归问题。Sacks 等 [1989] 在计算机实验的背景下描述 GPR(其中观测值 yy 是无噪声的)并讨论了一些有趣的方向,例如协方差函数中参数的优化(参见我们的 第 5 章)和实验设计(即选择提供有关 ff 的最多信息的 x\mathbf{x} 点的数量)。作者描述了许多建模的计算机模拟,包括一个示例,其中响应变量是电路中的时钟异步,输入是六个晶体管宽度。Santner 等 [2003] 撰写了一本关于使用高斯过程进行计算机实验设计和分析的新书。

Williams 和 Rasmussen [1996] 描述了机器学习中的高斯过程回归,并描述了协方差函数中参数的优化,另请参见 Rasmussen [1996]。如第 4.2.3 节和 Neal [1996] 中所述,他们受到与无限神经网络连接的启发而使用高斯过程。上述线性岭回归的 “核化” 也称为核岭回归,参见例如 Saunders 等 [1998]。