Class interpm_krige_optim (o2scl)

O2scl : Class List

template<class vec_t = boost::numeric::ublas::vector<double>, class mat_x_t = o2scl::const_matrix_view_table<>, class mat_x_row_t = const const_matrix_row_gen<o2scl::const_matrix_view_table<>>, class mat_y_t = o2scl::matrix_view_table<>, class mat_y_col_t = const matrix_column_gen<o2scl::matrix_view_table<>>, class mat_inv_kxx_t = boost::numeric::ublas::matrix<double>, class mat_inv_t = o2scl_linalg::matrix_invert_det_cholesky<boost::numeric::ublas::matrix<double>>, class vec3_t = std::vector<std::vector<std::vector<double>>>>
class interpm_krige_optim : public o2scl::interpm_base<boost::numeric::ublas::vector<double>, o2scl::const_matrix_view_table<>, o2scl::matrix_view_table<>>

Multi-dimensional interpolation using Kriging (a.k.a. Gaussian process interpolation)

See also the Higher-dimensional Interpolation section of the User’s guide.

Note

The type mat_y_t should be const_matrix_view_table, but the code would need to be modified because when the data is rescaled the user copy of the data is modified. This should be fixed in the future, as there’s no need to change the user-specified data, it would be better to just modify local storage.

Settings

bool timing

If true, output timing results.

bool full_min

If true, use the full minimizer.

bool skip_optim

If true, skip optimization.

bool rescale

If true, then the data will be automatically rescaled (default true)

bool keep_matrix

If true, keep \( K^{-1} \) (default true)

Minimizers and settings

mmin_simp2<multi_funct, ubvector> def_mmin

Default minimizer.

diff_evo_adapt alt_mmin

Alternate minimizer.

bool use_alt_mmin

If true, use the alternate minimizer.

inline void set_mmin(mmin_base<multi_funct, multi_funct, ubvector> &mb)

Set the minimizer to use.

size_t mode

Function to minimize (default mode_loo_cv)

size_t loo_npts

Number of points to test for cross validation (default 100)

static const size_t mode_loo_cv_bf = 1

Leave-one-out cross validation (brute force version)

static const size_t mode_max_lml = 2

Minus Log-marginal-likelihood.

static const size_t mode_loo_cv = 3

Leave-one-out cross validation (default)

static const size_t mode_final = 10

No optimization (for internal use)

Covariance function

std::vector<std::shared_ptr<mcovar_base<vec_t, mat_x_row_t>>> cf

Pointer to the covariance function.

inline int set_covar(std::vector<std::shared_ptr<mcovar_base<vec_t, mat_x_row_t>>> covar, vec3_t &param_lists)

Set the covariance function and parameter lists.

Other functions

inline virtual int addl_const(size_t iout, double &ret)

Additional constraints to add to the fit.

inline virtual double qual_fun(size_t iout, int &success)

Function to optimize the covariance parameters.

inline double min_fun(size_t iout, size_t n, const ubvector &v, double max_val)

Minimization function for the covariance parameters.

Constructor and destructor

inline interpm_krige_optim()
inline virtual ~interpm_krige_optim()

Basic interpolation functions

inline virtual int set_data(size_t n_in, size_t n_out, size_t n_pts, mat_x_t &user_x, mat_y_t &user_y)

Initialize the data for the interpolation.

inline virtual int eval(const vec_t &xp, vec_t &yp) const

Evaluate the interpolation at point x, returning y.

inline virtual int eval_unc(const vec_t &xp, vec_t &yp, vec_t &y_unc) const

Evaluate the interpolation at point x, returning y and the uncertainties in y_unc.

template<class vec4_t>
inline int eval_tl(const vec_t &x0, vec4_t &y0) const

Given input vector x store the result of the interpolation in y.

Uncertainties and derivatives

template<class vec4_t>
inline void sigma(const vec_t &x0, vec4_t &dy0) const

Return the interpolation uncertainty from the Gaussian process.

template<class vec2_t, class vec4_t>
inline void deriv(const vec2_t &x0, vec4_t &y0, size_t ix)

Given input vector x store the result of the interpolation in y.

template<class vec2_t, class vec4_t>
inline void deriv2(const vec2_t &x0, vec4_t &y0, size_t ix, size_t iy)

Given input vector x store the result of the interpolation in y.

Public Types

typedef boost::numeric::ublas::matrix<double> ubmatrix

Internal matrix type for the Nelder-Mead simplex.

typedef interpm_krige_optim<vec_t, mat_x_t, mat_x_row_t, mat_y_t, mat_y_col_t, mat_inv_kxx_t, mat_inv_t, vec3_t> class_t

Typedef for this type.

Protected Attributes

std::vector<ubvector> Kinvf

Inverse covariance matrix times function vector.

std::vector<mat_inv_kxx_t> inv_KXX

The inverse of the covariance matrix for each output quantity.

mat_inv_t mi

The matrix inversion object.

vec3_t plists

The list of parameter values.

The first index runs over output quantities, the second over different sets of values, and the third over parameters.

std::vector<double> qual

The quality factor of the optimization for each output function.

mmin_base<multi_funct, multi_funct, ubvector> *mp

Pointer to the user-specified minimizer.

mat_x_t x

The input data.

mat_y_t y

The output data.

bool data_set

True if the data has been specified.

ubvector mean_y

The output means for rescaling.

ubvector std_y

The output standard deviations for rescaling.