Да, ядро 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 (каппа), поэтому, задав 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, не имеет значения, поскольку мы устанавливаем их все l будет таким же, как с нашим выбором W и kappa).
Если ваш ввод был более чем одномерным, вы можете установить active_dims=list(range(X.shape[1] - 1)) для первого ядра (ядер) и active_dims=[X.shape[1]-1] для ядра Coregion.