浅析高斯过程回归（Gaussian process regression）[通俗易懂]

大家好，欢迎来到IT知识分享网。

前言

高斯过程回归的和其他回归算法的区别是：一般回归算法给定输入X，希望得到的是对应的Y值，拟合函数可以有多种多样，线性拟合、多项式拟合等等，而高斯回归是要得到函数f(x)的分布，那么是如何实现的呢？

对于数据集 D:(X,Y) ，令 $f(x^{_{i}})=y_{i}$ ，从而得到向量 f=[f(x_1),f(x_2),...,f(x_n)] , 将所需要预测的 $x_{i}$ 的集合定义为，对应的预测值为, 根据贝叶斯公式有：

$p(f*|f)=\frac{p(f|f*)p(f*)}{p(f)}=\frac{p(f,f*)}{p(f)}$

高斯回归首先要计算数据集中样本之间的联合概率分布, $f\sim N(\mu ,K)$ ， $\mu$ 为 $f(x_{1}),f(x_{2}),...,f(x_{n})$ 的均值所组成的向量，K为其协方差矩阵，再根据需要预测的的先验概率分布 $f*\sim N(\mu* ,K*)$ 与 $f\sim N(\mu ,K)$ ，来计算出的后验概率分布。

其中共有两个核心问题：（1）如何计算和方差矩阵（2）如何具体如何计算的概率分布。

{
x_1,x_2,...,x_n }的协方差矩阵，一定要注意是自变量x的协方差矩阵，因为我们在进行对对应的 f(x*) 的预测时，是依据和{
}之间的协方差，来判断其y值之间的差异性，从而基于概率公式得出预测值。

协方差矩阵的计算

定义函数 m(x)=E(f(x)) ， $k(x,x^{^{T}})=K$ ，则由概率论基本公式可得：

$f(x)\sim N(m(x),k(x,x^{^{T}}))\\$

m(x)=E(f(x))

$k(x,x^{^{T}})=E(f(x)-m(x)(f(x)-m(x))^{T})$

利用这样的原始公式来计算协方差矩阵K是十分不方便的，我们先来理解一下高斯过程是如何利用Gaussian distribution 来描述样本的，先来看图1和图2:

浅析高斯过程回归（Gaussian process regression）[通俗易懂] — 图1

图1中很明显可以看出来y值的大小与自变量x值的取值相关性很小，对于这样的数据，我们可以给出f(x)的先验分布：

$f(x)\sim N(0,\begin{bmatrix} 1 & & \\ &... & \\& &1 \end{bmatrix})$

其协方差矩阵为对角阵，任意不同x值之间的协方差均设为0
对于图2，f(x)受x值影响较大，x值相近，y值也相近，呈现出较高的相关性，我们可以给出其先验分布， $f(x)\sim N(0,K)$ ，通过以上两个例子我们发现， f(x) 的协方差矩阵与 x_1,x_2,...,x_n 相关度很高，如果我们能够使得得到一个 $k(x,x^{^{T}})$ 函数，只需输入x，即可得到对应 f(x) 的协方差矩阵，将给我们的计算带来极大的便利。那么用什么样的函数来取代原始公式才是合理的呢？

我们需要知道以下两个定理：

（1）协方差矩阵必须是半正定阵

（2）kernel fountion都是半正定阵。这就意味着我们在学习SVM的时候所学过的核函数形式都可以用

当然应用最广的是RBF kernel，即如下式：

$k(x,{x}')=\alpha ^{2}exp(-\frac{1}{2l^{2}})(x-{x}')^2)$

其中 $\alpha$ 为超参数，是需要通过学习进行确定的参数，那么我们只需要通过监督学习的方式学习到合适的kernel，即可很方便的计算出的协方差矩阵。

如何学习kernel的参数？很简单kernel k(x,{x}') 优劣的评价标准就是在要 $f(x)\sim N(m(x),k(x,x^{^{T}}))\\$ 的条件下，让 p(Y|X) 最大，为了方便求导我们将目标函数设为 $logp(Y|X)=log N(\mu ,K_{y})$ ，接下来利用梯度下降法来求最优值即可：

$\frac{\partial logp(Y|X) }{\theta }=\frac{1}{2}y^{T}K_{y}^{-1}\frac{\partial K_{y} }{\theta }K_{y}^{-1}y-\frac{1}{2}tr(K_{y}^{-1}\frac{\partial K_{y} }{\theta })$

后验概率分布的计算

图3

如图所示，红色的点代表已知的数据点，即训练集 $D={(x_{1},f(x_1)),(x_{2},f(x_2)),...,(x_{n},f(x_n)}$ ，给出的先验分布：

$f(x)\sim N(\mu,K)$ ，其中K为协方差矩阵，我们以RBF kernel来计算： $K_{ij}=\alpha ^{2}exp((-\frac{1}{2l^{2}})(x_i-x_j)^2)$

绿色的叉代表需要估计的点的X值，我们给定其先验分布 $f(x*)\sim N(\mu *,K(x*,x*))$

已知 $f(x)\sim N(\mu,K)$ ， $f(x*)\sim N(\mu *,K(x*,x*))$ ，可计算其其联合概率分布的先验：

$\bigl(\begin{smallmatrix} f\\ f* \end{smallmatrix}\bigr)\sim (\begin{pmatrix} \mu \\ \mu*\end{pmatrix},\begin{pmatrix} K& K* \\ K*^{T}&K**\end{pmatrix})$

其中 K** 为 f(x*) 的协方差矩阵， K**=k(X*,X*) ， K*=k(X,X*)

有了 p(f) 的先验分布，以及上面计算的 p(f,f*) ，根据贝叶斯公式可以计算 p(f*|f) 的后验概率:

$p(f*|f)=\frac{p(f|f*)p(f*)}{p(f)}=\frac{p(f,f*)}{p(f)}$

从而得出对于的估计， $f*\sim (\mu {}',K{}')$

$\mu{}'=K^TK^{-1}f$

$K{}'=K*^TK^{-1}K*+K**$

图3拟合的结果如图5所示，图5展示的拟合曲线是仅显示均值的曲线。图6、图7展现了高斯过程回归对于模型误差的估计能力，在图6中在后半部分数据波动较大时，其估计的方差表征了在这一部分的波动情况；在图7中也较好的展示了，在位置数据点，高斯过程回归对于该点函数值均值和方差的估计。

求解高斯过程的工具sklearn.gaussian_process

在sklearn中提供了十分方便的gaussian_process库，可以用来进行高斯过程求解，我们可以调用GaussianProcessRegressor来进行高斯过程回归的求解，函数具体介绍如下：

参数：

参数：	kernel：内核对象内核指定GP的协方差函数。如果通过None，则默认使用内核“1.0 * RBF（1.0）”。请注意，内核的超参数在拟合期间进行了优化。 alpha：float或array-like，可选（默认值：1e-10）在拟合期间将值添加到核矩阵的对角线。较大的值对应于观察中增加的噪声水平。通过确保计算值形成正定矩阵，这还可以防止拟合期间潜在的数值问题。如果传递数组，则它必须具有与用于拟合的数据相同的条目数，并用作与数据点相关的噪声级别。请注意，这相当于添加带有c = alpha的WhiteKernel。允许直接将噪声级别指定为参数主要是为了方便和与Ridge的一致性。优化器：字符串或可调用，可选（默认值：“fmin_l_bfgs_b”）可以是内部支持的优化器之一，用于优化由字符串指定的内核参数，也可以是作为可调用方传递的外部定义的优化器。n_restarts_optimizer：int，optional（默认值：0）优化程序重新启动的次数，用于查找最大化对数边际可能性的内核参数。优化程序的第一次运行是从内核的初始参数执行的，其余的那些（如果有的话）来自采样的log-uniform，从允许的the-value的空间中随机均匀。如果大于0，则所有边界必须是有限的。请注意，n_restarts_optimizer == 0表示执行了一次运行。 normalize_y：boolean，optional（默认值：False）目标值y是否被归一化，即，观察到的目标值的平均值变为零。如果预期目标值的平均值与零相差很大，则此参数应设置为True。启用后，归一化有效地根据数据修改GP的先验，这与可能性原则相矛盾; 因此，默认情况下禁用标准化。 copy_X_train：bool，optional（默认值：True）如果为True，则将训练数据的持久副本存储在对象中。否则，仅存储对训练数据的引用，如果在外部修改数据，则可能导致预测发生变化。 random_state：int，RandomState实例或None，可选（默认值：None）用于初始化中心的生成器。如果是int，则random_state是随机数生成器使用的种子; 如果是RandomState实例，则random_state是随机数生成器; 如果没有，随机数生成器所使用的RandomState实例np.random。
属性：	X_train_：类似数组，shape =（n_samples，n_features）训练数据中的特征值（也是预测所需） y_train_：array-like，shape =（n_samples，[n_output_dims]）训练数据中的目标值（也是预测所需） kernel_：内核对象内核用于预测。内核的结构与作为参数传递的内核相同，但具有优化的超参数 L_：类似数组，shape =（n_samples，n_samples）下三角形Cholesky分解内核 `X_train_` alpha_：array-like，shape =（n_samples，）核空间中训练数据点的双重系数 log_marginal_likelihood_value_：float 对数边际可能性 `self.kernel_.theta`

kernel：内核对象

内核指定GP的协方差函数。如果通过None，则默认使用内核“1.0 * RBF（1.0）”。请注意，内核的超参数在拟合期间进行了优化。

alpha：float或array-like，可选（默认值：1e-10）

在拟合期间将值添加到核矩阵的对角线。较大的值对应于观察中增加的噪声水平。通过确保计算值形成正定矩阵，这还可以防止拟合期间潜在的数值问题。如果传递数组，则它必须具有与用于拟合的数据相同的条目数，并用作与数据点相关的噪声级别。请注意，这相当于添加带有c = alpha的WhiteKernel。允许直接将噪声级别指定为参数主要是为了方便和与Ridge的一致性。

优化器：字符串或可调用，可选（默认值：“fmin_l_bfgs_b”）

可以是内部支持的优化器之一，用于优化由字符串指定的内核参数，也可以是作为可调用方传递的外部定义的优化器。n_restarts_optimizer：int，optional（默认值：0）

优化程序重新启动的次数，用于查找最大化对数边际可能性的内核参数。优化程序的第一次运行是从内核的初始参数执行的，其余的那些（如果有的话）来自采样的log-uniform，从允许的the-value的空间中随机均匀。如果大于0，则所有边界必须是有限的。请注意，n_restarts_optimizer == 0表示执行了一次运行。

normalize_y：boolean，optional（默认值：False）

目标值y是否被归一化，即，观察到的目标值的平均值变为零。如果预期目标值的平均值与零相差很大，则此参数应设置为True。启用后，归一化有效地根据数据修改GP的先验，这与可能性原则相矛盾; 因此，默认情况下禁用标准化。

copy_X_train：bool，optional（默认值：True）

如果为True，则将训练数据的持久副本存储在对象中。否则，仅存储对训练数据的引用，如果在外部修改数据，则可能导致预测发生变化。

random_state：int，RandomState实例或None，可选（默认值：None）

用于初始化中心的生成器。如果是int，则random_state是随机数生成器使用的种子; 如果是RandomState实例，则random_state是随机数生成器; 如果没有，随机数生成器所使用的RandomState实例np.random。

属性：

X_train_：类似数组，shape =（n_samples，n_features）

训练数据中的特征值（也是预测所需）

y_train_：array-like，shape =（n_samples，[n_output_dims]）

训练数据中的目标值（也是预测所需）

kernel_：内核对象

内核用于预测。内核的结构与作为参数传递的内核相同，但具有优化的超参数

L_：类似数组，shape =（n_samples，n_samples）

下三角形Cholesky分解内核 X_train_

alpha_：array-like，shape =（n_samples，）

核空间中训练数据点的双重系数

log_marginal_likelihood_value_：float

对数边际可能性 self.kernel_.theta

`fit`（X，y）	拟合高斯过程回归模型。
`get_params`（[深]）	获取此估算工具的参数。
`log_marginal_likelihood`（[theta，eval_gradient]）	返回训练数据的theta的log-marginal似然。
`predict`（X [，return_std，return_cov]）	使用高斯过程回归模型进行预测
`sample_y`（X [，n_samples，random_state]）	从高斯过程中抽取样本并在X处进行评估。
`score`（X，y [，sample_weight]）	返回预测的确定系数R ^ 2。
`set_params`（** PARAMS）	设置此估算器的参数。

提供了一个简单的小程序，大家可以用来做一些小实验：

import numpy as np
import matplotlib.pyplot as plt
import random
from sklearn.gaussian_process import GaussianProcessRegressor
a=np.random.random(50).reshape(50,1)
b=a*2+np.random.random(50).reshape(50,1)
plt.scatter(a,b,marker = 'o', color = 'r', label='3', s = 15)
plt.show()
gaussian=GaussianProcessRegressor()
fiting=gaussian.fit(a,b)
c=np.linspace(0.1,1,100)
d=gaussian.predict(c.reshape(100,1))
plt.scatter(a,b,marker = 'o', color = 'r', label='3', s = 15)
plt.plot(c,d)
plt.show()

IT知识分享网

免责声明：本站所有文章内容,图片，视频等均是来源于用户投稿和互联网及文摘转载整编而成，不代表本站观点，不承担相关法律责任。其著作权各归其原作者或其出版社所有。如发现本站有涉嫌抄袭侵权/违法违规的内容,侵犯到您的权益，请在线联系站长,一经查实,本站将立刻删除。本文来自网络,若有侵权，请联系删除，如若转载，请注明出处：https://yundeesoft.com/7487.html

浅析高斯过程回归（Gaussian process regression）[通俗易懂]

前言

协方差矩阵的计算

后验概率分布的计算

求解高斯过程的工具sklearn.gaussian_process

相关推荐

发表回复