8 #ifndef IMPISD_MULTIVARIATE_FNORMAL_SUFFICIENT_H
9 #define IMPISD_MULTIVARIATE_FNORMAL_SUFFICIENT_H
11 #include <IMP/isd/isd_config.h>
12 #include "internal/timer.h"
18 #include <Eigen/Dense>
19 #include <Eigen/Cholesky>
20 #include <IMP/isd/internal/cg_eigen.h>
22 IMPISD_BEGIN_NAMESPACE
23 using Eigen::MatrixXd;
24 using Eigen::VectorXd;
27 #define IMP_MVN_TIMER_NFUNCS 11
93 VectorXd FM_, Fbar_, epsilon_,Peps_;
94 double JF_,lJF_,norm_,lnorm_;
95 MatrixXd P_,W_,Sigma_,FX_,PW_;
99 Eigen::LDLT<MatrixXd, Eigen::Upper> ldlt_;
101 bool flag_FM_, flag_FX_, flag_Fbar_,
102 flag_W_, flag_Sigma_, flag_epsilon_,
103 flag_PW_, flag_P_, flag_ldlt_, flag_norms_,
107 bool use_cg_, first_PW_, first_PWP_;
110 IMP::Pointer<internal::ConjugateGradientEigen> cg_;
112 internal::CallTimer<IMP_MVN_TIMER_NFUNCS> timer_;
124 const VectorXd& FM,
const MatrixXd& Sigma,
double factor=1);
137 const VectorXd& FM,
int Nobs,
const MatrixXd& W,
138 const MatrixXd& Sigma,
double factor=1);
141 double density()
const;
144 double evaluate()
const;
147 VectorXd evaluate_derivative_FM()
const;
150 MatrixXd evaluate_derivative_Sigma()
const;
153 double evaluate_derivative_factor()
const;
156 MatrixXd evaluate_second_derivative_FM_FM()
const;
161 MatrixXd evaluate_second_derivative_FM_Sigma(
unsigned l)
const;
164 MatrixXd evaluate_second_derivative_Sigma_Sigma(
unsigned k,
unsigned l)
const;
168 void set_FX(
const MatrixXd& f);
169 MatrixXd get_FX()
const;
171 void set_Fbar(
const VectorXd& f);
172 VectorXd get_Fbar()
const;
174 void set_jacobian(
double f);
175 double get_jacobian()
const;
176 void set_minus_log_jacobian(
double f);
177 double get_minus_log_jacobian()
const;
179 void set_FM(
const VectorXd& f);
180 VectorXd get_FM()
const;
182 void set_W(
const MatrixXd& f);
183 MatrixXd get_W()
const;
185 void set_Sigma(
const MatrixXd& f);
186 MatrixXd get_Sigma()
const;
188 void set_factor(
double f);
189 double get_factor()
const;
195 void set_use_cg(
bool use,
double tol);
201 VectorXd get_Sigma_eigenvalues()
const;
204 double get_Sigma_condition_number()
const;
207 MatrixXd solve(MatrixXd B)
const;
210 double get_mean_square_residuals()
const;
216 double get_minus_exponent()
const;
222 double get_minus_log_normalization()
const;
226 out <<
"MultivariateFNormalSufficient: "
227 << N_ <<
" observations of "
228 << M_ <<
" variables " <<std::endl,
237 MatrixXd get_P()
const;
238 void set_P(
const MatrixXd& P);
241 MatrixXd get_PW()
const;
242 MatrixXd compute_PW_direct()
const;
243 MatrixXd compute_PW_cg()
const;
244 void set_PW(
const MatrixXd& PW);
247 VectorXd get_Peps()
const;
248 void set_Peps(
const VectorXd& Peps);
251 VectorXd get_epsilon()
const;
252 void set_epsilon(
const VectorXd& eps);
256 Eigen::LDLT<MatrixXd, Eigen::Upper> get_ldlt()
const;
257 void set_ldlt(
const Eigen::LDLT<MatrixXd, Eigen::Upper>& ldlt);
260 void set_norms(
double norm,
double lnorm);
261 std::vector<double> get_norms()
const;
264 double trace_WP()
const;
267 MatrixXd compute_PTP()
const;
270 MatrixXd compute_PWP()
const;
273 void compute_epsilon();