使用python+sklearn實現高斯過程回歸

本示例以兩種不同的方式來實現簡單的一維回歸：

在這兩種情況下，核函數的參數都是用極大似然原理估計的。

下圖以逐點95％置信區間的形式說明了高斯過程模型的插值性質(interpolating property)及其概率性質。

注意，參數alpha 用作訓練點之間假定協方差的Tikhonov正則化。

print(__doc__)

# 作者: Vincent Dubourg <[e-mail]>
#         Jake Vanderplas <[e-mail]>
#         Jan Hendrik Metzen <[e-mail]>s
# 許可證: BSD 3 clause

import numpy as np
from matplotlib import pyplot as plt

from sklearn.gaussian_process import GaussianProcessRegressor
from sklearn.gaussian_process.kernels import RBF, ConstantKernel as C

np.random.seed(1)


def f(x):
    """預測函數"""
    return x * np.sin(x)

# 
#  首先是無噪聲的情況
X = np.atleast_2d([1., 3., 5., 6., 7., 8.]).T

# 觀測結果
y = f(X).ravel()

# 對輸入空間進行網格劃分，以評估實函數，預測和MSE
x = np.atleast_2d(np.linspace(0, 10, 1000)).T

# 實例化高斯過程模型
kernel = C(1.0, (1e-3, 1e3)) * RBF(10, (1e-2, 1e2))
gp = GaussianProcessRegressor(kernel=kernel, n_restarts_optimizer=9)

# 使用參數的最大似然估計來擬合數據
gp.fit(X, y)

# 在網格x軸上進行預測（也要求MSE）
y_pred, sigma = gp.predict(x, return_std=True)

# 根據MSE繪製函數，預測和95％置信區間
plt.figure()
plt.plot(x, f(x), 'r:', label=r'$f(x) = x\,\sin(x)$')
plt.plot(X, y, 'r.', markersize=10, label='Observations')
plt.plot(x, y_pred, 'b-', label='Prediction')
plt.fill(np.concatenate([x, x[::-1]]),
         np.concatenate([y_pred - 1.9600 * sigma,
                        (y_pred + 1.9600 * sigma)[::-1]]),
         alpha=.5, fc='b', ec='None', label='95% confidence interval')
plt.xlabel('$x$')
plt.ylabel('$f(x)$')
plt.ylim(-10, 20)
plt.legend(loc='upper left')

# 
# 現在是有噪聲的情況
X = np.linspace(0.1, 9.9, 20)
X = np.atleast_2d(X).T

# 觀測結果和噪音
y = f(X).ravel()
dy = 0.5 + 1.0 * np.random.random(y.shape)
noise = np.random.normal(0, dy)
y += noise

# 實例化高斯過程模型
gp = GaussianProcessRegressor(kernel=kernel, alpha=dy ** 2,
                              n_restarts_optimizer=10)

# 使用參數的最大似然估計來擬合數據
gp.fit(X, y)

# 在網格x軸上進行預測（也要求MSE）
y_pred, sigma = gp.predict(x, return_std=True)

# 根據MSE繪製函數，預測和95％置信區間
plt.figure()
plt.plot(x, f(x), 'r:', label=r'$f(x) = x\,\sin(x)$')
plt.errorbar(X.ravel(), y, dy, fmt='r.', markersize=10, label='Observations')
plt.plot(x, y_pred, 'b-', label='Prediction')
plt.fill(np.concatenate([x, x[::-1]]),
         np.concatenate([y_pred - 1.9600 * sigma,
                        (y_pred + 1.9600 * sigma)[::-1]]),
         alpha=.5, fc='b', ec='None', label='95% confidence interval')
plt.xlabel('$x$')
plt.ylabel('$f(x)$')
plt.ylim(-10, 20)
plt.legend(loc='upper left')

plt.show()

腳本的總運行時間： ( 0 分 0.606 秒)

估計的內存使用量： 8 MB

下載python原始碼: plot_gpr_noisy_targets.py下載Jupyter notebook原始碼: plot_gpr_noisy_targets.ipynb☆☆☆為方便大家查閱，小編已將scikit-learn學習路線專欄文章統一整理到公眾號底部菜單欄，同步更新中，關注公眾號，點擊左下方「系列文章」，如圖：歡迎大家和我一起沿著scikit-learn文檔這條路線，一起鞏固機器學習算法基礎。（添加微信：mthler，備註：sklearn學習，一起進【sklearn機器學習進步群】開啟打怪升級的學習之旅。）