Secure your code as it's written. Use Snyk Code to scan source code in minutes - no build needed - and fix issues immediately.
X: Training data matrix with inputs of size NxNx, with Nx number of
inputs to the GP.
Y: Training data matrix with outpyts of size (N x Ny).
hyper: Array with hyperparameters [ell_1 .. ell_Nx sf sn].
inputmean: Mean from the last GP iteration of size (1 x Nx)
inputcov: Covariance matrix from the last GP iteration of size (Nx x Nx)
# Returns
mean: Array of the output mean of size (Ny x 1).
covariance: Covariance matrix of size (Ny x Ny).
"""
hyper = ca.log(hyper)
Ny = len(invK)
N, Nx = ca.MX.size(X)
mean = ca.MX.zeros(Ny, 1)
beta = ca.MX.zeros(N, Ny)
log_k = ca.MX.zeros(N, Ny)
v = X - ca.repmat(inputmean, N, 1)
covariance = ca.MX.zeros(Ny, Ny)
#TODO: Fix that LinsolQr don't work with the extended graph?
A = ca.SX.sym('A', inputcov.shape)
[Q, R2] = ca.qr(A)
determinant = ca.Function('determinant', [A], [ca.exp(ca.trace(ca.log(R2)))])
for a in range(Ny):
beta[:, a] = ca.mtimes(invK[a], Y[:, a])
iLambda = ca.diag(ca.exp(-2 * hyper[a, :Nx]))
R = inputcov + ca.diag(ca.exp(2 * hyper[a, :Nx]))
iR = ca.mtimes(iLambda, (ca.MX.eye(Nx) - ca.solve((ca.MX.eye(Nx)
if sym.name() in child.symbol_dict:
name = child.symbol_dict[sym.name()][1]
v = x_i[label][name]
ind = self.q_i[child][name]
sym2 = MX.zeros(sym.size())
sym2[ind] = v
sym2 = reshape(sym2, sym.shape)
c = substitute(c, sym, sym2)
break
for nghb in self.q_ij.keys():
for label, child in nghb.group.items():
if sym.name() in child.symbol_dict:
name = child.symbol_dict[sym.name()][1]
v = z_ij[nghb.label, label, name]
ind = self.q_ij[nghb][child][name]
sym2 = MX.zeros(sym.size())
sym2[ind] = v
sym2 = reshape(sym2, sym.shape)
c = substitute(c, sym, sym2)
break
for name, s in self.par_global.items():
if s.name() == sym.name():
c = substitute(c, sym, par[name])
lb, ub = con[1], con[2]
self.define_constraint(c, lb, ub)
# construct problem
prob, buildtime = self.father_updx.construct_problem(
self.options, str(self._index), problem)
self.problem_upd_xz = prob
self.father_updx.init_transformations(self.problem.init_primal_transform,
self.problem.init_dual_transform)
self.init_var_dd()
def gp_casadi(invK, hyp, X, Y, z):
E = len(invK)
n = ca.MX.size(X[:, 1])[0]
D = ca.MX.size(X[1, :])[1]
mean = ca.MX.zeros(E, 1)
var = ca.MX.zeros(E, 1)
for a in range(E):
ell = ca.MX(hyp[a, 0:D])
sf2 = ca.MX(hyp[a, D]**2)
m = 0 #hyp[a, D + 2]
kss = covSEard2(z, z, ell, sf2)
ks = ca.MX.zeros(n, 1)
for i in range(n):
ks[i] = covSEard2(X[i, :] - m, z - m, ell, sf2)
#ks = repmat()
ksK = ca.mtimes(ks.T, invK[a])
mean[a] = ca.mtimes(ksK, Y[:, a] - m) + m
var[a] = kss - ca.mtimes(ksK, ks)
return mean, var
def calc_NLL_casadi(hyper, X, Y):
""" Objective function """
# Calculate NLL
n, D = ca.MX.size(X)
ell = hyper[:D]
sf2 = hyper[D]**2
lik = hyper[D + 1]**2
m = hyper[D + 2]
K = ca.MX.zeros(n, n)
for i in range(n):
for j in range(n):
K[i, j] = covSEard2(X[i, :] - m, X[j, :] - m, ell, sf2)
K = K + lik * ca.MX.eye(n)
K = (K + K.T) * 0.5 # Make sure matrix is symmentric
A = ca.SX.sym('A', ca.MX.size(K))
#L = ca.chol(A) # Should check for PD !!!
cholesky = ca.Function('Cholesky', [A], [ca.chol(A)])
L = cholesky(K).T
logK = 2 * ca.sum1(ca.MX.log(ca.fabs(ca.diag(L))))
invLy = ca.solve(L, Y - m)
def gp_taylor_approx(invK, X, F, hyper, D, inputmean, inputcov):
#hyper = ca.MX.log(hyper)
E = len(invK)
n = ca.MX.size(F[:, 1])[0]
mean = ca.MX.zeros(E, 1)
var = ca.MX.zeros(D, 1)
v = X - ca.repmat(inputmean, n, 1)
covar = ca.MX.zeros(E, E)
covariance = ca.MX.zeros(E, E)
d_mean = ca.MX.zeros(E, 1)
dd_var = ca.MX.zeros(E, E)
for a in range(E):
ell = hyper[a, :D]
w = 1 / ell**2
sf2 = ca.MX(hyper[a, D]**2)
iK = ca.MX(invK[a])
alpha = ca.mtimes(iK, F[:, a])
kss = sf2
ks = ca.MX.zeros(n, 1)
def GP_noisy_input(invK, X, F, hyper, D, inputmean, inputcov):
hyper = ca.MX.log(hyper)
E = len(invK)
n = ca.MX.size(F[:, 1])[0]
mean = ca.MX.zeros(E, 1)
beta = ca.MX.zeros(n, E)
log_k = ca.MX.zeros(n, E)
v = X - ca.repmat(inputmean, n, 1)
#invK = MX(invK)
covariance = ca.MX.zeros(E, E)
A = ca.SX.sym('A', inputcov.shape)
[Q, R2] = ca.qr(A)
determinant = ca.Function('determinant', [A], [ca.exp(ca.trace(ca.log(R2)))])
for a in range(E):
beta[:, a] = ca.mtimes(invK[a], F[:, a])
iLambda = ca.diag(ca.exp(-2 * hyper[a, :D]))
R = inputcov + ca.diag(ca.exp(2 * hyper[a, :D]))
iR = ca.mtimes(iLambda, (ca.MX.eye(D) - ca.solve((ca.MX.eye(D) + ca.mtimes(inputcov, iLambda)), (ca.mtimes(inputcov, iLambda)))))
def gp_taylor_approx(invK, X, F, hyper, D, inputmean, inputcov):
#hyper = ca.MX.log(hyper)
E = len(invK)
n = ca.MX.size(F[:, 1])[0]
mean = ca.MX.zeros(E, 1)
var = ca.MX.zeros(D, 1)
v = X - ca.repmat(inputmean, n, 1)
covar = ca.MX.zeros(E, E)
covariance = ca.MX.zeros(E, E)
d_mean = ca.MX.zeros(E, 1)
dd_var = ca.MX.zeros(E, E)
for a in range(E):
ell = hyper[a, :D]
w = 1 / ell**2
sf2 = ca.MX(hyper[a, D]**2)
iK = ca.MX(invK[a])
alpha = ca.mtimes(iK, F[:, a])
kss = sf2
ks = ca.MX.zeros(n, 1)
for i in range(n):
ks[i] = covSEard2(X[i, :], inputmean, ell, sf2)
for child in q_j.keys():
for name in q_j[child].keys():
z = z_ij[str(nghb), child.label, name]
l = l_ij[str(nghb), child.label, name]
obj -= mtimes(l.T, z)
self.define_objective(obj)
# construct constraints
for con in self.constraints:
c = con[0]
for sym in symvar(c):
for label, child in self.group.items():
if sym.name() in child.symbol_dict:
name = child.symbol_dict[sym.name()][1]
v = x_i[label][name]
ind = self.q_i[child][name]
sym2 = MX.zeros(sym.size())
sym2[ind] = v
sym2 = reshape(sym2, sym.shape)
c = substitute(c, sym, sym2)
break
for nghb in self.q_ij.keys():
for label, child in nghb.group.items():
if sym.name() in child.symbol_dict:
name = child.symbol_dict[sym.name()][1]
v = z_ij[nghb.label, label, name]
ind = self.q_ij[nghb][child][name]
sym2 = MX.zeros(sym.size())
sym2[ind] = v
sym2 = reshape(sym2, sym.shape)
c = substitute(c, sym, sym2)
break
for name, s in self.par_i.items():
def GP_noisy_input(invK, X, F, hyper, D, inputmean, inputcov):
hyper = ca.MX.log(hyper)
E = len(invK)
n = ca.MX.size(F[:, 1])[0]
mean = ca.MX.zeros(E, 1)
beta = ca.MX.zeros(n, E)
log_k = ca.MX.zeros(n, E)
v = X - ca.repmat(inputmean, n, 1)
#invK = MX(invK)
covariance = ca.MX.zeros(E, E)
A = ca.SX.sym('A', inputcov.shape)
[Q, R2] = ca.qr(A)
determinant = ca.Function('determinant', [A], [ca.exp(ca.trace(ca.log(R2)))])
for a in range(E):
beta[:, a] = ca.mtimes(invK[a], F[:, a])
iLambda = ca.diag(ca.exp(-2 * hyper[a, :D]))
R = inputcov + ca.diag(ca.exp(2 * hyper[a, :D]))
iR = ca.mtimes(iLambda, (ca.MX.eye(D) - ca.solve((ca.MX.eye(D) + ca.mtimes(inputcov, iLambda)), (ca.mtimes(inputcov, iLambda)))))
T = ca.mtimes(v, iR)
c = ca.exp(2 * hyper[a, D]) / ca.sqrt(determinant(R)) * ca.exp(ca.sum2(hyper[a, :D]))