高斯后验过程 (Python)

Question

我使用下面的代码创建并采样了平均值 = 0 的联合高斯先验：

\n\n

import numpy as np\nimport matplotlib.pyplot as plt \nfrom math import pi \nfrom scipy.spatial.distance import cdist\nimport scipy.stats as sts\n\nx_prior = np.linspace(-10,10,101)\nx_prior = x_prior.reshape(-1,1)\nmu = np.zeros(x_prior.shape)\n\n#defining the Kernel for the covariance function\n\ndef sec(a,b, length_scale , sigma) : \n    K = sigma * np.exp(-1/(2*length_scale) * cdist(a,b)**2)\n    return K \n\n#defining the Gaussian Process prior\n\ndef GP(a , b, mu , kernel , length_scale, sigma , samples ) :\n    f = np.random.multivariate_normal(mu.flatten(), kernel(a ,b , length_scale , sigma ) , samples)\n    return f\n\nprior = GP(x_prior ,x_prior, mu , sec , 100, 1 , 5)\n\nplt.figure()\nplt.grid()\nplt.title(\'samples from the Gaussian prior\')\nplt.plot(x_prior , prior.T)\nplt.show()\n

\n\n

\n\n

然后，当添加一些“观察到的”数据时，我希望计算这些点的后验，但这就是我陷入困境的地方。

\n\n

这是我引入新数据的代码：

\n\n

x_train = np.array([-10,-8,5,-1,2])\nx_train = x_train.reshape(-1,1)\ndef straight_line(m , x , c):\n    y = 5*x + c\n    return y\nytrain = straight_line(5 , x_train , 0)\n

\n\n

据我了解，您可以根据与观察到的数据关联的先前 x 值和新 x 值来计算新数据的条件分布。

\n\n

那么您是否希望通过对平均值进行某种更改以包含新的 y 值来更新多元变量，使其成为后验变量？

\n\n

我已使用以下资源来尝试此操作：

\n\n

http://katbailey.github.io/post/gaussian-processes-for-dummies/ \n https://www.robots.ox.ac.uk/~mebden/reports/GPtutorial.pdf

\n\n

但我真的很想了解每个阶段会发生什么，以及为什么，这样当我得到后验（我做不到）时，我确切地知道我是如何到达那里的。

\n\n

这是我一直在尝试实施的一些解决方案，但到目前为止还没有成功：

\n\n

K_train = sec(x_train , x_train , 1,1)\nK_prior = sec(x_prior , x_prior , 1,1)\nK_pt =  sec(x_prior , x_train , 1,1)\nK_tp = sec(x_train , x_prior ,  1,1)  ##\xc2\xa0= k_tp transpose\nprior = sts.multivariate_normal(mu.flatten(), K_prior) \n#mean_test = np.dot(K_p , np.linalg.inv(K_prior))\nmean_function = np.dot(np.dot(K_tp ,np.linalg.inv(K_prior).T) , prior )\ncovariance_function = K_train - np.dot(np.dot(K_tp ,np.linalg.inv(K_prior).T) , K_pt) \n

\n

Answer 1

只是为所有看过此内容的人提供更新。我通过阅读这篇论文找到了解决方案：

https://arxiv.org/pdf/1711.10834.pdf

和以下代码：

mean_function = np.dot(np.dot(K_pt ,np.linalg.inv(K_train)), ytrain) 

covariance_function = K_prior - np.dot(np.dot(K_pt ,np.linalg.inv(K_train)) , K_tp) 

f = np.random.multivariate_normal(mean_function[:,0],covariance_function , 100)

其中 f 是从中采样的后关节高斯