WillFeng's Blog
0%

（转载）简单线性回归-基于高斯核函数

发表于 分类于

高斯过程

高斯过程 Gaussian Processes 是概率论和数理统计中随机过程的一种，是多元高斯分布的扩展，被应用于机器学习、信号处理等领域。本文对高斯过程进行公式推导、原理阐述、可视化以及代码实现，介绍了以高斯过程为基础的高斯过程回归 Gaussian Process Regression 基本原理、超参优化、高维输入等问题。

  • 一元高斯分布公式
    其概率密度公式如下：

  • 核函数(协方差函数)
    核函数是一个高斯过程的核心，核函数决定了一个高斯过程的性质。核函数在高斯过程中起的作用是生成一个协方差矩阵（相关系数矩阵），衡量任意两个点之间的“距离”。最常用的一个核函数为高斯核函数，也成为径向基函数 RBF。其基本形式如下。其中 σ 和 h 是高斯核的超参数。

代码实现

1.定义核函数

1
2
3
4
5
6
7
8
# 定义核函数
   def kernel(self, x1, x2):
       m,n = x1.shape[0], x2.shape[0]
       dist_matrix = np.zeros((m,n), dtype=float)
       for i in range(m):
           for j in range(n):
               dist_matrix[i][j] = np.sum((x1[i]-x2[j])**2)
       return np.exp(-0.5/self.h**2*dist_matrix)

2.主要代码

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
# 导入相关库
import numpy as np
import matplotlib.pyplot as plt
# 定义高斯过程类
class GPR:
   def __init__(self, h):
       self.is_fit = False
       self.train_x, self.train_y = None, None
       self.h = h
       
   def fit(self, x, y):
       self.train_x = np.asarray(x)
       self.train_y = np.asarray(y)
       self.is_fit = True
       
   def predict(self, x):
       if not self.is_fit:
           print("Sorry! GPR Model can't fit!")
           return
       
       x = np.asarray(x)
       kff = self.kernel(x,x)
       kyy = self.kernel(self.train_x, self.train_x)
       kfy = self.kernel(x, self.train_x)
       kyy_inv =  np.linalg.inv(kyy + 1e-8*np.eye(len(self.train_x)))
       
       mu = kfy.dot(kyy_inv).dot(self.train_y)
       return mu
       
       
   # 定义核函数
   def kernel(self, x1, x2):
       m,n = x1.shape[0], x2.shape[0]
       dist_matrix = np.zeros((m,n), dtype=float)
       for i in range(m):
           for j in range(n):
               dist_matrix[i][j] = np.sum((x1[i]-x2[j])**2)
       return np.exp(-0.5/self.h**2*dist_matrix)

# 创造训练集
train_x = np.arange(0,10).reshape(-1,1)
train_y = np.cos(train_x) + train_x
# 制造槽点
train_y = train_y + np.random.normal(0, 0.01, size=train_x.shape)
# 显示训练集的分布
plt.figure()
plt.scatter(train_x, train_y, label="train", c="red", marker="x")
plt.legend()
plt.show()

输出的图：

1

3.学习代码

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
# 创建训练集
test_x = np.arange(0, 10, 0.1).reshape(-1,1)

# 针对不同h值得到拟合图像
h=0.1
for i in range(10):
    gpr = GPR(h)
    gpr.fit(train_x, train_y)
    mu = gpr.predict(test_x)
    test_y = mu.ravel()
    plt.figure()
    plt.title("h=%.2f"%(h))
    plt.plot(test_x, test_y, label="predict")
    plt.scatter(train_x, train_y, label="train", c="red", marker="x")
    plt.legend()
    plt.show()
    h += 0.1

输出图像：

2

2

2

2

写在最后

本人是刚接触机器学习的小白，这篇文章让对老师布置的作业更好的理解。文章转载自CSDN（Vince_Cheng），文章地址转载仅作为本人学习使用的笔记，同时也是我的第一篇博客，相信日后我也能写出原创的文章。感谢这位大佬的分享，侵删！