def sim_multitask_GP(times,length,noise_vars,K_f,trainfrac):
"""
draw from a multitask GP.
we continue to assume for now that the dim of the input space is 1, ie just time
M: number of tasks (labs/vitals/time series)
train_frac: proportion of full M x N data matrix Y to include
"""
M = np.shape(K_f)[0]
N = len(times)
n = N*M
K_t = OU_kernel_np(length,times) #just a correlation function
Sigma = np.diag(noise_vars)
K = np.kron(K_f,K_t) + np.kron(Sigma,np.eye(N)) + 1e-6*np.eye(n)
L_K = np.linalg.cholesky(K)
y = np.dot(L_K,np.random.normal(0,1,n)) #Draw normal
#get indices of which time series and which time point, for each element in y
ind_kf = np.tile(np.arange(M),(N,1)).flatten('F') #vec by column
ind_kx = np.tile(np.arange(N),(M,1)).flatten()
#randomly dropout some fraction of fully observed time series
perm = np.random.permutation(n)
n_train = int(trainfrac*n)
train_inds = perm[:n_train]
y_ = y[train_inds]
ind_kf_ = ind_kf[train_inds]
ind_kx_ = ind_kx[train_inds]
return y_,ind_kf_,ind_kx_
评论列表
文章目录