高斯过程(Gaussian Process)回归预测,例子,代码及可视化展示

news2024/11/15 12:43:38

高斯过程指的是一组随机变量的集合,这个集合里面的任意有限个随机变量都服从联合正态分布。(联合正态分布是指多个随机变量的联合分布满足正态分布。联合分布是指多个随机变量同时满足的概率分布,一个常见的例子是考虑两个随机变量:X 表示一个人的年龄,Y 表示他们的身高。假设 X 和 Y 的联合分布是正态分布。在这种情况下,联合分布可以描述不同年龄的人在不同身高下出现的概率。我们可以使用联合概率密度函数来计算某个年龄范围内的人在某个身高范围内出现的概率。)

高斯过程可以用于回归拟合。先来看一个简单的例子。一个地点的气温,我们有过去10年6月1号的气温数据,想要预测明年6月1日的气温,那么假定那一天的气温服从高斯分布,根据过去的数据,可以算出均值\mu和方差\sigma,则明年那一天的气温t\sim N(\mu, \sigma)满足这样一个正态分布,最有可能的气温就是均值\mu

那么如果是一个一元5次函数比如

f(x)=0.03x^{5}+0.2x^{4}-0.1x^{3}-2.4x^{2}-2.5x+6

画出图像和几个已知点

import numpy as np
import matplotlib.pyplot as plt


def f(x):
    coefs = [6, -2.5, -2.4, -0.1, 0.2, 0.03]
    total = 0
    for exp, coef in enumerate(coefs):
       total += coef * (x ** exp)
    return total


x_obs = np.array([-4, -3.6, -2.0, -1.5, 0])
y_obs = f(x_obs)
x_s = np.linspace(-5, 1, 80)

fig, ax = plt.subplots()
legend_list = []
ax.plot(x_s, f(x_s))
legend_list.append("true")
ax.scatter(x_obs, y_obs)
legend_list.append("sample points")

plt.legend(legend_list)
plt.show()

在已知几个点的情况下,如何通过高斯过程回归预测这个函数呢?

首先回顾一下高斯概率分布,一维高斯分布概率密度函数。

 x有95%的可能性出现在95%置信区间[-1.96, 1.96]。调整均值和方差,96%置信区间就会移动。

如果是二维高斯分布概率密度函数,下图。

代码如下

import numpy as np
import matplotlib.pyplot as plt

# 设置均值和协方差矩阵
mu = np.array([0, 0])   # 均值
cov = np.array([[1, 0], [0, 1]])   # 协方差矩阵

# 生成二维高斯分布的概率密度函数
x, y = np.meshgrid(np.linspace(-3, 3, 100), np.linspace(-3, 3, 100))
pos = np.dstack((x, y))
inv_cov = np.linalg.inv(cov)
pdf = np.exp(-0.5 * np.einsum('...k,kl,...l->...', pos - mu, inv_cov, pos - mu)) / (2 * np.pi * np.sqrt(np.linalg.det(cov)))

# 绘制二维高斯分布图
plt.contourf(x, y, pdf, levels=20)  # 使用contourf函数填充等高线图
plt.colorbar()  # 添加颜色条

# 设置图形标题和轴标签
plt.xlabel('x')
plt.ylabel('y')

# 显示图形
plt.show()

那么我们进行一次采样,得到的点大概率是出现在黄色区域,如果我们把这个点表示为p(y_{1},y_{2}),y1,y2分别表示横纵坐标。通过调整协方差矩阵,黄色区域就会移动,我们就能得到不同的y。

如果是n维高斯分布,取样一次得到1个点的n个维度值y_{1},...,y_{n}由此可以联想到,是否可以通过对n维高斯分布进行一次取样,来得到前文中那个5次函数的回归预测值呢?当然是可以的。

现在有一系列的x值

x_s = np.linspace(-5, 1, 80)

通过对高维高斯分布进行一次取样,可以得到一系列y_s,从x_s映射到y_s的关键就在于协方差矩阵,通过核函数,可以构造不同的协方差矩阵。这里我们选择最简单的核函数。

import numpy as np


def k(xs, ys, sigma=1, l=1):
    dx = np.expand_dims(xs, 1) - np.expand_dims(ys, 0)
    tmp = (sigma ** 2) * np.exp(-((dx / l) ** 2) / 2)
    return tmp


x_s = np.linspace(-5, 1, 80)
kernal = k(x_s, x_s)

预备知识告一段落,下面来看高斯过程回归是如何实现的:

先说明字母含义,x和y是我们已知的采样点,对应代码x_obs, y_obs。

x_obs = np.array([-4, -3.6, -2.0, -1.5, 0])
y_obs = f(x_obs)

x_{*}对应我们代码的x_s,y_{*}对应代码y_s。m(x)表示求均值,K是核函数运算得到的矩阵。

高斯过程回归的原理是:yy_{*}的联合概率分布(先验分布)符合高斯分布。

通过采样得到观测值y_obs以后,则y_{*}的后验分布满足均值为\mu_{*},方差为\Sigma _{*}的高斯分布

m和K的求法:

K = k(x_obs, x_obs)
K_s = k(x_obs, x_s)
K_ss = k(x_s, x_s)

K_sTKinv = np.matmul(K_s.T, np.linalg.pinv(K))

mu_s = m(x_s) + np.matmul(K_sTKinv, y_obs - m(x_obs))
Sigma_s = K_ss - np.matmul(K_sTKinv, K_s)

高斯过程回归结果如下图

解释一下这张图,红色实线是函数真实值,蓝色实线是高斯过程回归的结果。三根虚线表示:对80维高斯分布(因为我们的x_s是80维数组)进行三次采样,得到的3个点,每个点80个维度坐标值组成的曲线。每个维度上的值,都是满足均值为蓝色实线,方差为灰色背景范围的高斯分布。在已知观测点上,方差为0,而距离观测点越远,方差范围越大,可信度也就越低。

完整代码如下:

import numpy as np
import matplotlib.pyplot as plt


def f(x):
    coefs = [6, -2.5, -2.4, -0.1, 0.2, 0.03]
    total = 0
    for exp, coef in enumerate(coefs):
       total += coef * (x ** exp)
    return total


def m(x):
    tmp = np.zeros_like(x)
    return tmp


def k(xs, ys, sigma=1, l=1):
    dx = np.expand_dims(xs, 1) - np.expand_dims(ys, 0)
    tmp = (sigma ** 2) * np.exp(-((dx / l) ** 2) / 2)
    return tmp


x_obs = np.array([-4, -3.6, -2.0, -1.5, 0])
y_obs = f(x_obs)
x_s = np.linspace(-5, 1, 80)
K = k(x_obs, x_obs)
K_s = k(x_obs, x_s)
K_ss = k(x_s, x_s)

K_sTKinv = np.matmul(K_s.T, np.linalg.pinv(K))

mu_s = m(x_s) + np.matmul(K_sTKinv, y_obs - m(x_obs))
Sigma_s = K_ss - np.matmul(K_sTKinv, K_s)

stds = np.sqrt(Sigma_s.diagonal())
tmp = np.flip(x_s, 0)
err_xs = np.concatenate((x_s, np.flip(x_s, 0)))
err_ys = np.concatenate((mu_s + 2 * stds, np.flip(mu_s - 2 * stds, 0)))

fig, ax = plt.subplots()
legend_list = []
ax.plot(x_s, f(x_s), color='r')
legend_list.append("real")
ax.scatter(x_obs, y_obs)
legend_list.append("sample points")
ax.plot(x_s, mu_s, color='b')
legend_list.append("mean")
ax.fill_between(x_s, mu_s-2*stds, mu_s+2*stds, alpha=0.3, color='gray')
legend_list.append("err")
for i in range(3):
    ax.plot(x_s, np.random.multivariate_normal(mu_s, Sigma_s), linestyle='--')
    legend_list.append("test" + str(i))

plt.legend(legend_list)
plt.show()

本文来自互联网用户投稿,该文观点仅代表作者本人,不代表本站立场。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如若转载,请注明出处:http://www.coloradmin.cn/o/658122.html

如若内容造成侵权/违法违规/事实不符,请联系多彩编程网进行投诉反馈,一经查实,立即删除!

相关文章

【2023年最新】提高分类模型指标的六大方案详解

文章目录 数据增强特征选择调整模型参数模型集成迁移学习模型解释完结 当今,机器学习模型得到了广泛的应用,其中分类模型是其中最常见和重要的一种。在实际应用中,如何提高分类模型的指标,使其在不同场景下表现更佳并且具有更好的…

【Web网页制作】影视主题网页制作web页面开发(附源码)

【写在前面】 其实之前我就写过一篇关于我的家乡的页面,也有不少网友和我私下反馈,让我多出一些关于页面制作的,于是乎我就今天晚上抽出点时间来整理一篇关于影视内容的web页面制作,希望能够得到大家的喜欢。 【涉及内容】 web网页…

什么是进程?程序又是什么?程序运行在操作系统来看是怎么做的?

相信很多人都看到过课本上写的进程的概念,那么真的理解了吗? 课本上是这样讲的,课本概念:程序的一个执行实例,正在执行的程序等。 那么进程到底是什么?我先把内核层面上的概念拿出来:内核观点&a…

【算法基础】基础算法

快速排序 模板题:785. 快速排序 - AcWing题库 思路: 定义一个x(一般喜欢用中间的),我们快速排序,让x左边的都比它小,同时让右边的都比它大。然后像二分一样不断细分,缩小范围进行同…

flink主要组件及高可用配置

背景 flink不论运行在哪种环境,例如Yarn,Mesos,Kebernute以及独立集群,每个应用都会包含重要的几个组件,本文就来讲述下flink的主要组件以及如何实现flink的高可用配置 flink主要组件 如图所示,flink主要…

TiDB 7.1 资源管控特性试用

作者: 啦啦啦啦啦 原文来源: https://tidb.net/blog/3ddb423a 一.背景 印象里 2022 年初的时候就有小伙伴在 asktug 提出 TiDB 未来是否会有多租户功能的问题了,没想到这么快就已经 GA 了。资源管控特性(Resource Control&…

【论文阅读】Weighted Boxes Fusion(WBF)模型融合原理解读

论文地址:https://arxiv.org/pdf/1910.13302.pdf 代码地址:GitHub - ZFTurbo/Weighted-Boxes-Fusion: Set of methods to ensemble boxes from different object detection models, including implementation of "Weighted boxes fusion (WBF)"…

chatgpt赋能python:Python中如何转化大小写

Python中如何转化大小写 在Python编程中,转化字符串的大小写是一个常见的操作。Python提供了内置函数和字符串方法来实现此操作。本文将介绍如何使用Python中的这些函数和方法来转换字符串的大小写。 使用内置函数str.upper()和str.lower() str.upper()函数将字符…

共创开源生态 | 小米肖翔荣获“2023中国开源优秀人物”奖

6月15-16日,以“开源创新 数字化转型 智能化重构”为主题的“第十八届开源中国・开源世界高峰论坛”在北京成功召开。小米工程师肖翔凭借其在 Apache 基金会的开源贡献及在操作系统领域内的技术突破,荣获“2023中国开源优秀人物”奖。 Xiaomi …

一文读懂候选边界框Selective Search、AnchorBased、Anchor Free

目标检测是计算机视觉中的一项重要任务,主要目的是在图像或视频中识别并定位感兴趣的对象。为了实现这一目标,目标检测算法通常会生成一系列候选边界框,这些框包围了图像中可能存在的目标对象。候选边界框技术对于减少目标检测的计算复杂度和…

Vue中的数据可视化词云展示与词云生成

Vue中的数据可视化词云展示与词云生成 数据可视化是现代Web应用程序中的一个重要组成部分,它使得数据更加易于理解和分析。词云是一种非常流行的数据可视化形式,它可以用来展示文本数据中的主题和关键字。在本文中,我们将介绍如何在Vue中使用…

chatgpt赋能python:Python怎么转化数据类型?

Python怎么转化数据类型? Python是一种高级编程语言,它已经成为了许多程序员的首选语言。在Python中,数据类型是非常重要的一部分。但是,当我们需要将数据从一种类型转换为另一种类型时,该怎么做呢?在本文…

深度学习:探索人工智能的新前沿

第一章:引言 人工智能(Artificial Intelligence,AI)作为一项前沿技术,在近年来取得了巨大的进展。其中,深度学习(Deep Learning)作为人工智能领域的一个重要分支,更是引…

如何在VMware上安装CentOS7?

目录 一、器材准备 二、创建一个虚拟机 三、安装Centos7系统 一、器材准备 1. Centos7及以上版本的iso镜像 链接:centos7镜像 提取码:ao3n 2. VMware15及以上版本的软件工具包 链接:VMware16安装包以及激活码 提取码:40pe 二、创…

Framework - Zygote

一、概念 Zygote是 Android 中的第一个进程,负责孵化(fork)其它进程,而它自己由 Linux 内核启动的用户级进程 Init 创建。 二、作用 应用程序不能直接以本地进程的形态运行,必须在一个独立的虚拟机中运行,一…

Springboot实现数据传输加解密

前言 先给大家看下效果,原本我们的请求是这样子的 加密后的数据传输是这样子的 加解密步骤: 1.前端请求前进行加密,然后发送到后端 2.后端收到请求后解密 3.后端返回数据前进行加密 4.前端拿到加密串后,解密数据 加解密算法&…

搭建TiDB负载均衡环境-LVS+KeepAlived实践

作者: 我是咖啡哥 原文来源: https://tidb.net/blog/f614b200 昨天,发了一篇使用HAproxyKP搭建TiDB负载均衡环境的文章,今天我们再用LVSKP来做个实验。 环境信息 TiDB版本:V7.1.0 haproxy版本:2.6.2 …

【EXCEL】如何查找特殊字符 问号‘?’星号 ‘*’

目录 0.环境 1.适用场景 1)直接搜索问号的结果: 2)修改【查找内容】后,搜索结果变为精准定位: 2.具体做法 0.环境 windows wps(或excel,这里试了,此问题wps和excel表格是通用…

chatgpt赋能python:Python如何计算圆周率π

Python如何计算圆周率π 圆周率,又称π,是数学中一个重要的常数,它与圆的周长和直径的比值始终保持不变。在计算机编程中,计算圆周率π也是一个颇具挑战的问题。本文介绍了使用Python编程语言来计算圆周率π的方法,希…

C语言进阶---指针的进阶

前言 指针的主题,我们在初级阶段的《指针》章节已经接触过了。我们直到指针的概念。 ​ 1、指针就是个变量,用来存放地址,地址唯一标识一块内存空间。 ​ 2、指针的大小是固定的4/8个字节(32为平台/64位平台) ​ 3、指…