49 using namespace Eigen;
66 void CKLCovarianceInferenceMethod::init()
69 "V is L'*V=diag(sW)*K",
72 "A is A=I-K*diag(sW)*inv(L)'*inv(L)*diag(sW)",
78 "Square root of noise matrix W",
81 "the gradient of the variational expection wrt sigma2",
84 "the gradient of the variational expection wrt mu",
104 Map<VectorXd> eigen_result(result.
vector, len);
107 eigen_result=eigen_alpha;
119 Map<VectorXd> eigen_mean(mean.
vector, mean.
vlen);
135 eigen_W=(2.0*eigen_log_neg_lambda.array().exp()).matrix();
137 Map<VectorXd> eigen_sW(m_sW.
vector, m_sW.
vlen);
138 eigen_sW=eigen_W.array().sqrt().matrix();
145 eigen_V=eigen_L.triangularView<Upper>().adjoint().solve(eigen_sW.asDiagonal()*eigen_K*
CMath::sq(
m_scale));
149 eigen_s2=(eigen_K.diagonal().array()*
CMath::sq(
m_scale)-(eigen_V.array().pow(2).colwise().sum().transpose())).abs().matrix();
158 "The length of gradients (%d) should the same as the length of parameters (%d)\n",
162 Map<VectorXd> eigen_sW(m_sW.
vector, m_sW.
vlen);
177 Map<VectorXd> eigen_dv(m_dv.
vector, m_dv.
vlen);
181 Map<VectorXd> eigen_df(m_df.
vector, m_df.
vlen);
183 MatrixXd eigen_U=eigen_L.triangularView<Upper>().adjoint().solve(MatrixXd(eigen_sW.asDiagonal()));
186 eigen_A=MatrixXd::Identity(len, len)-eigen_V.transpose()*eigen_U;
191 Map<VectorXd> eigen_dnlz_alpha(gradient.
vector, len);
192 Map<VectorXd> eigen_dnlz_log_neg_lambda(gradient.
vector+len, len);
200 eigen_dnlz_log_neg_lambda=(eigen_Sigma.array().pow(2)*2.0).matrix()*eigen_dv+eigen_s2;
201 eigen_dnlz_log_neg_lambda=eigen_dnlz_log_neg_lambda-(eigen_Sigma.array()*eigen_A.array()).rowwise().sum().matrix();
202 eigen_dnlz_log_neg_lambda=(eigen_log_neg_lambda.array().exp()*eigen_dnlz_log_neg_lambda.array()).matrix();
214 Map<VectorXd> eigen_mean(mean.
vector, mean.
vlen);
222 MatrixXd eigen_t=eigen_L.triangularView<Upper>().adjoint().solve(MatrixXd::Identity(eigen_L.rows(),eigen_L.cols()));
224 for(
index_t idx=0; idx<eigen_t.rows(); idx++)
225 trace +=(eigen_t.col(idx).array().pow(2)).sum();
228 float64_t result=-a+eigen_L.diagonal().array().
log().sum();
229 result+=0.5*(-eigen_K.rows()+eigen_alpha.dot(eigen_mu-eigen_mean)+trace);
238 Map<VectorXd> eigen_sW(m_sW.
vector, m_sW.
vlen);
243 Map<VectorXd> eigen_dv(m_dv.
vector, m_dv.
vlen);
244 Map<VectorXd> eigen_df(m_df.
vector, m_df.
vlen);
247 MatrixXd AdK=eigen_A*eigen_dK;
250 VectorXd z=AdK.diagonal()+(eigen_A.array()*AdK.array()).rowwise().sum().matrix()
251 -(eigen_A.transpose().array()*AdK.array()).colwise().sum().transpose().matrix();
254 return eigen_alpha.dot(eigen_dK*(eigen_alpha/2.0-eigen_df))-z.dot(eigen_dv);
270 MatrixXd::Identity(eigen_K.rows(), eigen_K.cols()));
271 MatrixXd LL=llt.matrixU();
272 MatrixXd tt=LL.triangularView<Upper>().adjoint().solve(MatrixXd::Identity(LL.rows(),LL.cols()));
274 for(
index_t idx=0; idx<tt.rows(); idx++)
275 trace+=(tt.col(idx).array().pow(2)).sum();
277 MatrixXd eigen_V=LL.triangularView<Upper>().adjoint().solve(eigen_K*
CMath::sq(
m_scale));
279 Map<VectorXd> eigen_s2(s2_tmp.
vector, s2_tmp.
vlen);
280 eigen_s2=(eigen_K.diagonal().array()*
CMath::sq(
m_scale)-(eigen_V.array().pow(2).colwise().sum().transpose())).abs().matrix();
287 nlml_def=-a+LL.diagonal().array().log().sum();
288 nlml_def+=0.5*(-eigen_K.rows()+trace);
290 if (nlml_new<=nlml_def)
virtual void update_approx_cov()
SGVector< float64_t > m_alpha
virtual SGVector< float64_t > get_variational_first_derivative(const TParameter *param) const =0
static SGMatrix< float64_t > get_choleksy(SGVector< float64_t > W, SGVector< float64_t > sW, SGMatrix< float64_t > kernel, float64_t scale)
virtual void update_deriv()
virtual void get_gradient_of_nlml_wrt_parameters(SGVector< float64_t > gradient)
The class Labels models labels, i.e. class assignments of objects.
virtual int32_t get_num_labels() const =0
static T sum(T *vec, int32_t len)
return sum(vec)
The variational Gaussian Likelihood base class. The variational distribution is Gaussian.
static SGMatrix< float64_t > get_inverse(SGMatrix< float64_t > L, SGMatrix< float64_t > kernel, SGVector< float64_t > sW, SGMatrix< float64_t > V, float64_t scale)
void log()
natural logarithm of vector elements
TParameter * get_parameter(int32_t idx)
virtual SGVector< float64_t > get_mean_vector(const CFeatures *features) const =0
An abstract class of the mean function.
virtual void set_variational_distribution(SGVector< float64_t > mu, SGVector< float64_t > s2, const CLabels *lab)
virtual float64_t get_derivative_related_cov(Eigen::MatrixXd eigen_dK)
SGMatrix< float64_t > m_Sigma
virtual ~CKLCovarianceInferenceMethod()
SGMatrix< float64_t > m_L
virtual float64_t lbfgs_optimization()
virtual SGVector< float64_t > get_variational_expection()=0
virtual void lbfgs_precompute()
The KL approximation inference method class.
The class Features is the base class of all feature objects.
virtual void update_alpha()
SGVector< float64_t > m_mu
SGVector< float64_t > m_s2
static float64_t log(float64_t v)
virtual SGVector< float64_t > get_diagonal_vector()
virtual CVariationalGaussianLikelihood * get_variational_likelihood() const
virtual SGVector< float64_t > get_alpha()
CKLCovarianceInferenceMethod()
virtual float64_t get_negative_log_marginal_likelihood_helper()
virtual bool parameter_hash_changed()
The Likelihood model base class.
SGMatrix< float64_t > m_ktrtr
virtual void update_chol()