docs: def _build_predict(self, Xnew, full_cov=False):
          """
          The posterior variance of F is given by
              q(f) = N(f | K alpha + mean, [K^-1 + diag(lambda**2)]^-1)
          Here we project this to F*, the values of the GP at Xnew which is given
          by
             q(F*) = N ( F* | K_{*F} alpha + mean, K_{**} - K_{*f}[K_{ff} +
                                             diag(lambda**-2)]^-1 K_{f*} )
          """
  
          # compute kernel things
          Kx = self.kern.K(self.X, Xnew)
          K = self.kern.K(self.X)
  
          # predictive mean
          f_mean = tf.matmul(Kx, self.q_alpha, transpose_a=True) + self.mean_function(Xnew)
  
          # predictive var
          A = K + tf.matrix_diag(tf.transpose(1. / tf.square(self.q_lambda)))
          L = tf.cholesky(A)
          Kx_tiled = tf.tile(tf.expand_dims(Kx, 0), [self.num_latent, 1, 1])
          LiKx = tf.matrix_triangular_solve(L, Kx_tiled)
          if full_cov:
              f_var = self.kern.K(Xnew) - tf.matmul(LiKx, LiKx, transpose_a=True)
          else:
              f_var = self.kern.Kdiag(Xnew) - tf.reduce_sum(tf.square(LiKx), 1)
          return f_mean, tf.transpose(f_var)