gpflow - GPFlow 同一 GP 的多个独立实现，不规则的采样时间/长度

Question

在 GPflow 中，我有多个时间序列，并且采样时间没有跨时间序列对齐，并且时间序列可能有不同的长度（纵向数据）。我假设它们是来自同一个 GP 的独立实现。用 svgp 处理这个问题的正确方法是什么，更普遍的是用 GPflow？我需要使用共同区域化吗？共同区域化笔记本假设相关轨迹，而我想要共享均值/内核但独立。

Answer 1

是的，Coregion在 GPflow 中实现的内核可以用来解决您的问题。

让我们从您描述的生成模型中设置一些数据，时间序列具有不同的长度：

import numpy as np
import gpflow
import matplotlib.pyplot as plt

Ns = [80, 90, 100]  # number of observations for three different realizations
Xs = [np.random.uniform(0, 10, size=N) for N in Ns]  # observation locations

# three different draws from the same GP:
k = gpflow.kernels.Matern52(variance=2.0, lengthscales=0.5)  # kernel

Ks = [k(X[:, None]) for X in Xs]
Ls = [np.linalg.cholesky(K) for K in Ks]
vs = [np.random.randn(N, 1) for N in Ns]
fs = [(L @ v).squeeze(axis=-1) for L, v in zip(Ls, vs)]

要实际设置 gpflow GP 模型的训练数据：

# output indicator for the observations: which timeseries is this?
os = [o * np.ones(N) for o, N in enumerate(Ns)]  # [0 ... 0, 1 ... 1, 2 ... 2]

# now assemble the three timeseries in single data set:
allX = np.concatenate(Xs)
allo = np.concatenate(os)
allf = np.concatenate(fs)
X = np.c_[allX, allo]
Y = allf[:, None]
assert X.shape == (sum(Ns), 2)
assert Y.shape == (sum(Ns), 1)

# now let's set up a copy of the original kernel:
k2 = gpflow.kernels.Matern52(active_dims=[0])  # the same as k above, but with different hyperparameters

# and a Coregionalization kernel that effectively says they are all independent:
kc = gpflow.kernels.Coregion(output_dim=len(Ns), rank=1, active_dims=[1])
kc.W.assign(np.zeros(kc.W.shape))
kc.kappa.assign(np.ones(kc.kappa.shape))
gpflow.set_trainable(kc, False)  # we want W and kappa fixed

内核定义了Coregion一个协方差矩阵 B = W Wᵀ + diag(kappa)，因此通过设置 W=0，我们规定零相关（独立实现）和 kappa=1（实际上是默认值）确保原始副本的方差超参数内核仍然是可解释的。

现在构建实际模型并优化超参数：

k2c = k2 * kc

m = gpflow.models.GPR((X, Y), k2c, noise_variance=1e-5)

opt = gpflow.optimizers.Scipy()
opt.minimize(m.training_loss, m.trainable_variables, compile=False)

它很好地恢复了初始方差和长度尺度超参数。

如果你想预测，你必须在Xnew参数中提供额外的“输出”列m.predict_f()，例如如下：

Xtest = np.linspace(0, 10, 100)
Xtest_augmented = np.c_[Xtest, np.zeros_like(Xtest)]
f_mean, f_var = m.predict_f(Xtest_augmented)

（无论您将输出列设置为 0、1 还是 2 都无关紧要，因为我们将它们都设置为与我们选择的W和相同kappa）。

如果您的输入不止一维，您可以设置 active_dims=list(range(X.shape[1] - 1))第一个内核和active_dims=[X.shape[1]-1]内核Coregion。