interp_krige.h
Go to the documentation of this file.
1 /*
2  -------------------------------------------------------------------
3 
4  Copyright (C) 2017-2020, Andrew W. Steiner
5 
6  This file is part of O2scl.
7 
8  O2scl is free software; you can redistribute it and/or modify
9  it under the terms of the GNU General Public License as published by
10  the Free Software Foundation; either version 3 of the License, or
11  (at your option) any later version.
12 
13  O2scl is distributed in the hope that it will be useful,
14  but WITHOUT ANY WARRANTY; without even the implied warranty of
15  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16  GNU General Public License for more details.
17 
18  You should have received a copy of the GNU General Public License
19  along with O2scl. If not, see <http://www.gnu.org/licenses/>.
20 
21  -------------------------------------------------------------------
22 */
23 #ifndef O2SCL_INTERP_KRIGE_H
24 #define O2SCL_INTERP_KRIGE_H
25 
26 /** \file interp_krige.h
27  \brief One-dimensional interpolation by Kriging
28 */
29 
30 #include <iostream>
31 #include <string>
32 
33 #include <gsl/gsl_sf_erf.h>
34 
35 #include <boost/numeric/ublas/vector.hpp>
36 #include <boost/numeric/ublas/matrix.hpp>
37 #include <boost/numeric/ublas/operation.hpp>
38 #include <boost/numeric/ublas/vector_proxy.hpp>
39 
40 #include <o2scl/interp.h>
41 #include <o2scl/cholesky.h>
42 #include <o2scl/lu.h>
43 #include <o2scl/columnify.h>
44 #include <o2scl/vector.h>
45 #include <o2scl/vec_stats.h>
46 #include <o2scl/min_brent_gsl.h>
47 #include <o2scl/constants.h>
48 
49 #ifndef DOXYGEN_NO_O2NS
50 namespace o2scl {
51 #endif
52 
53  /** \brief Interpolation by Kriging with a user-specified
54  covariance function
55 
56  See also the \ref intp_section section of the \o2 User's guide.
57 
58  \note The function \ref set() stores a pointer to the covariance
59  function and its derivatives and integrals so they cannot go out
60  of scope before any of the interpolation functions are called.
61 
62  \note This class is experimental.
63  */
64  template<class vec_t, class vec2_t=vec_t,
65  class covar_func_t=std::function<double(double,double)>,
66  class covar_integ_t=std::function<double(double,double,double)> >
67  class interp_krige : public interp_base<vec_t,vec2_t> {
68 
69 #ifdef O2SCL_NEVER_DEFINED
70  }{
71 #endif
72 
73  public:
74 
77  typedef boost::numeric::ublas::matrix_column<ubmatrix> ubmatrix_column;
78 
79  protected:
80 
81  /** \brief Inverse covariance matrix times function vector
82  */
84 
85  /** \brief Pointer to user-specified covariance function
86  */
87  covar_func_t *f;
88 
89  /** \brief Pointer to user-specified derivative
90  */
91  covar_func_t *df;
92 
93  /** \brief Pointer to user-specified second derivative
94  */
95  covar_func_t *df2;
96 
97  /** \brief Pointer to user-specified integral
98  */
99  covar_integ_t *intp;
100 
101  public:
102 
103  interp_krige() {
104  this->min_size=2;
105  matrix_mode=matrix_cholesky;
106  }
107 
108  virtual ~interp_krige() {}
109 
110  /// \name Select matrix inversion method
111  //@{
112  /** \brief Method for matrix inversion
113  */
114  size_t matrix_mode;
115  /// Use Cholesky decomposition
116  static const size_t matrix_cholesky=0;
117  /// Use LU decomposition
118  static const size_t matrix_LU=1;
119  //@}
120 
121  /// Initialize interpolation routine
122  virtual void set(size_t size, const vec_t &x, const vec2_t &y) {
123  O2SCL_ERR2("Function set(size_t,vec_t,vec_t) unimplemented ",
124  "in interp_krige.",o2scl::exc_eunimpl);
125  return;
126  }
127 
128  /** \brief Initialize interpolation routine, specifying derivatives
129  and integrals [not yet implemented]
130  */
131  virtual int set_covar_di_noise(size_t n_dim, const vec_t &x,
132  const vec_t &y, covar_func_t &fcovar,
133  covar_func_t &fderiv,
134  covar_func_t &fderiv2,
135  covar_func_t &finteg,
136  double noise_var, bool err_on_fail=true) {
137  O2SCL_ERR("Function set_covar_di_noise not yet implemented.",
139  return 0;
140  }
141 
142  /// Initialize interpolation routine
143  virtual int set_covar_noise(size_t n_dim, const vec_t &x, const vec_t &y,
144  covar_func_t &fcovar, double noise_var,
145  bool err_on_fail=true) {
146 
147  if (n_dim<this->min_size) {
148  O2SCL_ERR((((std::string)"Vector size, ")+szttos(n_dim)+", is less"+
149  " than "+szttos(this->min_size)+" in interp_krige::"+
150  "set().").c_str(),exc_einval);
151  }
152 
153  // Store pointer to covariance function
154  f=&fcovar;
155 
156  // Construct the KXX matrix
157  ubmatrix KXX(n_dim,n_dim);
158  for(size_t irow=0;irow<n_dim;irow++) {
159  for(size_t icol=0;icol<n_dim;icol++) {
160  if (irow>icol) {
161  KXX(irow,icol)=KXX(icol,irow);
162  } else if (irow==icol) {
163  KXX(irow,icol)=fcovar(x[irow],x[icol])+noise_var;
164  } else {
165  KXX(irow,icol)=fcovar(x[irow],x[icol]);
166  }
167  }
168  }
169 
170  if (matrix_mode==matrix_LU) {
171 
172  // Construct the inverse of KXX
173  ubmatrix inv_KXX(n_dim,n_dim);
174  o2scl::permutation p(n_dim);
175  int signum;
176  o2scl_linalg::LU_decomp(n_dim,KXX,p,signum);
177  if (o2scl_linalg::diagonal_has_zero(n_dim,KXX)) {
178  if (err_on_fail) {
179  O2SCL_ERR2("Matrix singular (LU method) ",
180  "in interp_krige::set_covar_noise().",
182  }
183  return 1;
184  }
185  o2scl_linalg::LU_invert<ubmatrix,ubmatrix,ubmatrix_column>
186  (n_dim,KXX,p,inv_KXX);
187 
188  // Inverse covariance matrix times function vector
189  this->Kinvf.resize(n_dim);
190  boost::numeric::ublas::axpy_prod(inv_KXX,y,this->Kinvf,true);
191 
192  } else {
193 
194  // Construct the inverse of KXX
195  int cret=o2scl_linalg::cholesky_decomp(n_dim,KXX,false);
196  if (cret!=0) {
197  if (err_on_fail) {
198  O2SCL_ERR2("Matrix singular (Cholesky method) ",
199  "in interp_krige::set_covar_noise().",
201  }
202  return 2;
203  }
204  ubmatrix &inv_KXX=KXX;
205  o2scl_linalg::cholesky_invert<ubmatrix>(n_dim,inv_KXX);
206 
207  // Inverse covariance matrix times function vector
208  Kinvf.resize(n_dim);
209  boost::numeric::ublas::axpy_prod(inv_KXX,y,Kinvf,true);
210 
211  }
212 
213  // Set parent data members
214  this->px=&x;
215  this->py=&y;
216  this->sz=n_dim;
217 
218  return 0;
219  }
220 
221  /// Initialize interpolation routine
222  virtual int set_covar(size_t n_dim, const vec_t &x, const vec_t &y,
223  covar_func_t &fcovar, bool err_on_fail=true) {
224  return set_covar_noise(n_dim,x,y,fcovar,0.0,err_on_fail);
225  }
226 
227  /// Give the value of the function \f$ y(x=x_0) \f$ .
228  virtual double eval(double x0) const {
229 
230  // Initialize to zero to prevent uninit'ed var. warnings
231  double ret=0.0;
232  for(size_t i=0;i<this->sz;i++) {
233  ret+=(*f)(x0,(*this->px)[i])*Kinvf[i];
234  }
235 
236  return ret;
237  }
238 
239  /// Give the value of the derivative \f$ y^{\prime}(x=x_0) \f$ .
240  virtual double deriv(double x0) const {
241  return 0.0;
242  }
243 
244  /** \brief Give the value of the second derivative
245  \f$ y^{\prime \prime}(x=x_0) \f$
246  */
247  virtual double deriv2(double x0) const {
248  return 0.0;
249  }
250 
251  /// Give the value of the integral \f$ \int_a^{b}y(x)~dx \f$ .
252  virtual double integ(double a, double b) const {
253  return 0.0;
254  }
255 
256  /// Return the type, \c "interp_krige".
257  virtual const char *type() const { return "interp_krige"; }
258 
259 #ifndef DOXYGEN_INTERNAL
260 
261  private:
262 
267 
268 #endif
269 
270  };
271 
272 
273  /** \brief One-dimensional interpolation using an
274  optimized covariance function
275 
276  See also the \ref intp_section section of the \o2 User's guide.
277 
278  \note This class is experimental.
279  */
280  template<class vec_t, class vec2_t=vec_t>
281  class interp_krige_optim : public interp_krige<vec_t,vec2_t> {
282 
283  public:
284 
287  typedef boost::numeric::ublas::matrix_column<ubmatrix> ubmatrix_column;
288 
289  protected:
290 
291  /// Function object for the covariance
292  std::function<double(double,double)> ff;
293 
294  /// The covariance function length scale
295  double len;
296 
297  /// The quality factor of the optimization
298  double qual;
299 
300  /// The covariance function
301  double covar(double x1, double x2) {
302  return exp(-(x1-x2)*(x1-x2)/len/len/2.0);
303  }
304 
305  /// The derivative of the covariance function
306  double deriv(double x1, double x2) {
307  return -exp(-(x1-x2)*(x1-x2)/len/len/2.0)/len/len*(x1-x2);
308  }
309 
310  /// The second derivative of the covariance function
311  double deriv2(double x1, double x2) {
312  return ((x1-x2)*(x1-x2)-len*len)*
313  exp(-(x1-x2)*(x1-x2)/len/len/2.0)/len/len/len/len;
314  }
315 
316  /// The integral of the covariance function
317  double integ(double x, double x1, double x2) {
318  exit(-1);
319  // This is probably wrong
320  //return 0.5*len*sqrt(o2scl_const::pi)*
321  //(gsl_sf_erf((x2-x)/len)+gsl_sf_erf((x-x1)/len));
322  }
323 
324  /// Pointer to the user-specified minimizer
326 
327  /** \brief Function to optimize the covariance parameters
328  */
329  double qual_fun(double x, double noise_var, int &success) {
330 
331  len=x;
332  success=0;
333 
334  size_t size=this->sz;
335 
336  if (mode==mode_loo_cv) {
337 
338  qual=0.0;
339  for(size_t k=0;k<size;k++) {
340 
341  // Leave one observation out
342  ubvector x2(size-1);
343  o2scl::vector_copy_jackknife((*this->px),k,x2);
344  ubvector y2(size-1);
345  o2scl::vector_copy_jackknife((*this->py),k,y2);
346 
347  // Construct the KXX matrix
348  ubmatrix KXX(size-1,size-1);
349  for(size_t irow=0;irow<size-1;irow++) {
350  for(size_t icol=0;icol<size-1;icol++) {
351  if (irow>icol) {
352  KXX(irow,icol)=KXX(icol,irow);
353  } else {
354  KXX(irow,icol)=exp(-pow((x2[irow]-x2[icol])/len,2.0)/2.0);
355  if (irow==icol) KXX(irow,icol)+=noise_var;
356  }
357  }
358  }
359 
360  if (this->matrix_mode==this->matrix_LU) {
361 
362  // Construct the inverse of KXX
363  ubmatrix inv_KXX(size-1,size-1);
364  o2scl::permutation p(size-1);
365  int signum;
366  o2scl_linalg::LU_decomp(size-1,KXX,p,signum);
367  if (o2scl_linalg::diagonal_has_zero(size-1,KXX)) {
368  success=1;
369  return 1.0e99;
370  }
371  o2scl_linalg::LU_invert<ubmatrix,ubmatrix,ubmatrix_column>
372  (size-1,KXX,p,inv_KXX);
373 
374  // Inverse covariance matrix times function vector
375  this->Kinvf.resize(size-1);
376  boost::numeric::ublas::axpy_prod(inv_KXX,y2,this->Kinvf,true);
377 
378  } else {
379 
380  // Construct the inverse of KXX
381  int cret=o2scl_linalg::cholesky_decomp(size-1,KXX,false);
382  if (cret!=0) {
383  success=2;
384  return 1.0e99;
385  }
386  ubmatrix &inv_KXX=KXX;
387  o2scl_linalg::cholesky_invert<ubmatrix>(size-1,inv_KXX);
388 
389  // Inverse covariance matrix times function vector
390  this->Kinvf.resize(size-1);
391  boost::numeric::ublas::axpy_prod(inv_KXX,y2,this->Kinvf,true);
392 
393  }
394 
395  double ypred=0.0;
396  double yact=(*this->py)[k];
397  for(size_t i=0;i<size-1;i++) {
398  ypred+=exp(-pow(((*this->px)[k]-x2[i])/len,2.0)/2.0)*this->Kinvf[i];
399  }
400 
401  // Measure the quality with a chi-squared like function
402  qual+=pow(yact-ypred,2.0);
403 
404  }
405 
406  } else if (mode==mode_max_lml) {
407 
408  // Construct the KXX matrix
409  ubmatrix KXX(size,size);
410  for(size_t irow=0;irow<size;irow++) {
411  for(size_t icol=0;icol<size;icol++) {
412  if (irow>icol) {
413  KXX(irow,icol)=KXX(icol,irow);
414  } else {
415  KXX(irow,icol)=exp(-pow(((*this->px)[irow]-
416  (*this->px)[icol])/len,2.0)/2.0);
417  if (irow==icol) KXX(irow,icol)+=noise_var;
418  }
419  }
420  }
421 
422  // Note: We have to use LU here because O2scl doesn't yet
423  // have a lndet() function for Cholesky decomp
424 
425  // Construct the inverse of KXX
426  ubmatrix inv_KXX(size,size);
427  o2scl::permutation p(size);
428  int signum;
429  o2scl_linalg::LU_decomp(size,KXX,p,signum);
430  if (o2scl_linalg::diagonal_has_zero(size,KXX)) {
431  success=1;
432  return 1.0e99;
433  }
434  o2scl_linalg::LU_invert<ubmatrix,ubmatrix,ubmatrix_column>
435  (size,KXX,p,inv_KXX);
436 
437  double lndet=o2scl_linalg::LU_lndet<ubmatrix>(size,KXX);
438 
439  // Inverse covariance matrix times function vector
440  this->Kinvf.resize(size);
441  boost::numeric::ublas::axpy_prod(inv_KXX,*this->py,this->Kinvf,true);
442 
443  // Compute the log of the marginal likelihood, without
444  // the constant term
445  for(size_t i=0;i<size;i++) {
446  qual+=0.5*(*this->py)[i]*this->Kinvf[i];
447  }
448  qual+=0.5*lndet;
449  }
450 
451  return qual;
452  }
453 
454  public:
455 
457  nlen=20;
458  full_min=false;
459  mp=&def_min;
460  verbose=0;
462  }
463 
464  /// \name Function to minimize and various option
465  //@{
466  /// Leave-one-out cross validation
467  static const size_t mode_loo_cv=1;
468  /// Minus Log-marginal-likelihood
469  static const size_t mode_max_lml=2;
470  /// Function to minimize (default \ref mode_loo_cv)
471  size_t mode;
472  ///@}
473 
474  /// Verbosity parameter
475  int verbose;
476 
477  /** \brief Number of length scale points to try when full minimizer
478  is not used (default 20)
479  */
480  size_t nlen;
481 
482  /// Default minimizer
484 
485  /// If true, use the full minimizer
486  bool full_min;
487 
488  /// Initialize interpolation routine
489  virtual int set_noise(size_t size, const vec_t &x, const vec2_t &y,
490  double noise_var, bool err_on_fail=true) {
491 
492  // Set parent data members
493  this->px=&x;
494  this->py=&y;
495  this->sz=size;
496 
497  int success=0;
498 
499  if (full_min) {
500 
501  if (verbose>1) {
502  std::cout << "interp_krige_optim: full minimization"
503  << std::endl;
504  }
505 
506  // Choose first interval as initial guess
507  double len_opt=x[1]-x[0];
508 
509  funct mf=std::bind
510  (std::mem_fn<double(double,double,int &)>
512  std::placeholders::_1,noise_var,std::ref(success));
513 
514  mp->min(len_opt,qual,mf);
515  len=len_opt;
516 
517  if (success!=0) {
518  if (err_on_fail) {
519  O2SCL_ERR2("Minimization failed in ",
520  "interp_krige_optim::set_noise().",
522  }
523  }
524 
525  } else {
526 
527  if (verbose>1) {
528  std::cout << "interp_krige_optim: simple minimization"
529  << std::endl;
530  }
531 
532  // Compute a finite-difference array
533  std::vector<double> diff(size-1);
534  for(size_t i=0;i<size-1;i++) {
535  diff[i]=fabs(x[i+1]-x[i]);
536  }
537 
538  // Range of the length parameter
539  double len_min=o2scl::vector_min_value
540  <std::vector<double>,double>(size-1,diff)/3.0;
541  double len_max=fabs(x[size-1]-x[0])*3.0;
542  double len_ratio=len_max/len_min;
543 
544  if (verbose>1) {
545  std::cout << "len (min,max,ratio): " << len_min << " "
546  << len_max << " "
547  << pow(len_ratio,((double)1)/((double)nlen-1))
548  << std::endl;
549  }
550 
551  // Initialize to zero to prevent uninit'ed var. warnings
552  double min_qual=0.0, len_opt=0.0;
553 
554  if (verbose>1) {
555  std::cout << "ilen len qual fail min_qual" << std::endl;
556  }
557 
558  // Loop over the full range, finding the optimum
559  bool min_set=false;
560  for(size_t j=0;j<nlen;j++) {
561  len=len_min*pow(len_ratio,((double)j)/((double)nlen-1));
562 
563  int success=0;
564  qual=qual_fun(len,noise_var,success);
565 
566  if (success==0 && (min_set==false || qual<min_qual)) {
567  len_opt=len;
568  min_qual=qual;
569  min_set=true;
570  }
571 
572  if (verbose>1) {
573  std::cout << "interp_krige_optim: ";
574  std::cout.width(2);
575  std::cout << j << " " << len << " " << qual << " "
576  << success << " " << min_qual << " "
577  << len_opt << std::endl;
578  }
579 
580  }
581 
582  // Now that we've optimized the covariance function,
583  // just use the parent class to interpolate
584  len=len_opt;
585 
586  }
587 
588  ff=std::bind(std::mem_fn<double(double,double)>
590  std::placeholders::_1,std::placeholders::_2);
591 
592  this->set_covar_noise(size,x,y,ff,noise_var);
593 
594  return 0;
595  }
596 
597  /// Initialize interpolation routine
598  virtual void set(size_t size, const vec_t &x, const vec2_t &y) {
599 
600  // Use the mean absolute value to determine noise
601  double mean_abs=0.0;
602  for(size_t j=0;j<size;j++) {
603  mean_abs+=fabs(y[j]);
604  }
605  mean_abs/=size;
606 
607  set_noise(size,x,y,mean_abs/1.0e8,true);
608 
609  return;
610  }
611 
612 
613  };
614 
615 #ifndef DOXYGEN_NO_O2NS
616 }
617 #endif
618 
619 #endif
boost::numeric::ublas::matrix< double >
o2scl::exc_esing
@ exc_esing
apparent singularity detected
Definition: err_hnd.h:93
o2scl::interp_krige_optim::verbose
int verbose
Verbosity parameter.
Definition: interp_krige.h:475
o2scl::interp_krige< vec_t, vec_t >::matrix_LU
static const size_t matrix_LU
Use LU decomposition.
Definition: interp_krige.h:118
o2scl::interp_krige_optim::set_noise
virtual int set_noise(size_t size, const vec_t &x, const vec2_t &y, double noise_var, bool err_on_fail=true)
Initialize interpolation routine.
Definition: interp_krige.h:489
o2scl::interp_base< vec_t, vec_t >::py
const vec_t * py
Dependent vector.
Definition: interp.h:128
o2scl::permutation
A class for representing permutations.
Definition: permutation.h:70
o2scl::interp_krige::df2
covar_func_t * df2
Pointer to user-specified second derivative.
Definition: interp_krige.h:95
o2scl::interp_krige_optim::integ
double integ(double x, double x1, double x2)
The integral of the covariance function.
Definition: interp_krige.h:317
boost::numeric::ublas::vector< double >
o2scl::interp_krige::intp
covar_integ_t * intp
Pointer to user-specified integral.
Definition: interp_krige.h:99
o2scl::interp_krige::deriv2
virtual double deriv2(double x0) const
Give the value of the second derivative .
Definition: interp_krige.h:247
o2scl::exc_efailed
@ exc_efailed
generic failure
Definition: err_hnd.h:61
o2scl::interp_krige::eval
virtual double eval(double x0) const
Give the value of the function .
Definition: interp_krige.h:228
O2SCL_ERR2
#define O2SCL_ERR2(d, d2, n)
Set an error, two-string version.
Definition: err_hnd.h:281
o2scl::interp_krige::type
virtual const char * type() const
Return the type, "interp_krige".
Definition: interp_krige.h:257
o2scl::interp_krige_optim::mp
min_base * mp
Pointer to the user-specified minimizer.
Definition: interp_krige.h:325
o2scl
The main O<span style='position: relative; top: 0.3em; font-size: 0.8em'>2</span>scl O$_2$scl names...
Definition: anneal.h:42
o2scl::interp_krige_optim::len
double len
The covariance function length scale.
Definition: interp_krige.h:295
o2scl::vector_min_value
data_t vector_min_value(size_t n, const vec_t &data)
Compute the minimum of the first n elements of a vector.
Definition: vector.h:1196
o2scl_linalg::LU_decomp
int LU_decomp(const size_t N, mat_t &A, o2scl::permutation &p, int &signum)
Compute the LU decomposition of the matrix A.
Definition: lu_base.h:86
o2scl::exc_eunimpl
@ exc_eunimpl
requested feature not (yet) implemented
Definition: err_hnd.h:99
o2scl::interp_krige::matrix_mode
size_t matrix_mode
Method for matrix inversion.
Definition: interp_krige.h:114
o2scl::interp_krige_optim::ff
std::function< double(double, double)> ff
Function object for the covariance.
Definition: interp_krige.h:292
o2scl_inte_qng_coeffs::x2
static const double x2[5]
Definition: inte_qng_gsl.h:66
o2scl::interp_krige_optim::deriv
double deriv(double x1, double x2)
The derivative of the covariance function.
Definition: interp_krige.h:306
o2scl::interp_krige_optim::set
virtual void set(size_t size, const vec_t &x, const vec2_t &y)
Initialize interpolation routine.
Definition: interp_krige.h:598
o2scl::success
@ success
Success.
Definition: err_hnd.h:47
o2scl::interp_krige_optim::mode_loo_cv
static const size_t mode_loo_cv
Leave-one-out cross validation.
Definition: interp_krige.h:467
o2scl::interp_krige_optim
One-dimensional interpolation using an optimized covariance function.
Definition: interp_krige.h:281
o2scl::interp_krige::integ
virtual double integ(double a, double b) const
Give the value of the integral .
Definition: interp_krige.h:252
o2scl_inte_qng_coeffs::x1
static const double x1[5]
Definition: inte_qng_gsl.h:48
o2scl::exc_einval
@ exc_einval
invalid argument supplied by user
Definition: err_hnd.h:59
o2scl_linalg::cholesky_decomp
int cholesky_decomp(const size_t M, mat_t &A, bool err_on_fail=true)
Compute the in-place Cholesky decomposition of a symmetric positive-definite square matrix.
Definition: cholesky_base.h:62
o2scl::interp_krige::set_covar_di_noise
virtual int set_covar_di_noise(size_t n_dim, const vec_t &x, const vec_t &y, covar_func_t &fcovar, covar_func_t &fderiv, covar_func_t &fderiv2, covar_func_t &finteg, double noise_var, bool err_on_fail=true)
Initialize interpolation routine, specifying derivatives and integrals [not yet implemented].
Definition: interp_krige.h:131
o2scl::interp_krige_optim::deriv2
double deriv2(double x1, double x2)
The second derivative of the covariance function.
Definition: interp_krige.h:311
o2scl::interp_krige_optim::full_min
bool full_min
If true, use the full minimizer.
Definition: interp_krige.h:486
o2scl::interp_krige_optim::nlen
size_t nlen
Number of length scale points to try when full minimizer is not used (default 20)
Definition: interp_krige.h:480
o2scl::interp_krige_optim::mode_max_lml
static const size_t mode_max_lml
Minus Log-marginal-likelihood.
Definition: interp_krige.h:469
o2scl::interp_krige::set
virtual void set(size_t size, const vec_t &x, const vec2_t &y)
Initialize interpolation routine.
Definition: interp_krige.h:122
o2scl::interp_krige::deriv
virtual double deriv(double x0) const
Give the value of the derivative .
Definition: interp_krige.h:240
o2scl::vector_copy_jackknife
void vector_copy_jackknife(const vec_t &src, size_t iout, vec2_t &dest)
From a given vector, create a new vector by removing a specified element.
Definition: vector.h:2379
o2scl::interp_krige_optim::covar
double covar(double x1, double x2)
The covariance function.
Definition: interp_krige.h:301
o2scl::interp_krige
Interpolation by Kriging with a user-specified covariance function.
Definition: interp_krige.h:67
o2scl::interp_krige_optim::def_min
min_brent_gsl def_min
Default minimizer.
Definition: interp_krige.h:483
o2scl::interp_krige_optim::qual
double qual
The quality factor of the optimization.
Definition: interp_krige.h:298
o2scl::min_base
One-dimensional minimization [abstract base].
Definition: min.h:42
o2scl::interp_krige::f
covar_func_t * f
Pointer to user-specified covariance function.
Definition: interp_krige.h:87
O2SCL_ERR
#define O2SCL_ERR(d, n)
Set an error with message d and code n.
Definition: err_hnd.h:273
o2scl::funct
std::function< double(double)> funct
One-dimensional function typedef in src/base/funct.h.
Definition: funct.h:48
o2scl::interp_krige::df
covar_func_t * df
Pointer to user-specified derivative.
Definition: interp_krige.h:91
o2scl::min_base::min
virtual int min(double &x, double &fmin, func_t &func)=0
Calculate the minimum min of func w.r.t 'x'.
o2scl::interp_krige::set_covar
virtual int set_covar(size_t n_dim, const vec_t &x, const vec_t &y, covar_func_t &fcovar, bool err_on_fail=true)
Initialize interpolation routine.
Definition: interp_krige.h:222
o2scl::interp_krige_optim::qual_fun
double qual_fun(double x, double noise_var, int &success)
Function to optimize the covariance parameters.
Definition: interp_krige.h:329
o2scl::interp_base< vec_t, vec_t >::px
const vec_t * px
Independent vector.
Definition: interp.h:125
o2scl_linalg::diagonal_has_zero
int diagonal_has_zero(const size_t N, mat_t &A)
Return 1 if at least one of the elements in the diagonal is zero.
Definition: lu_base.h:54
o2scl::szttos
std::string szttos(size_t x)
Convert a size_t to a string.
o2scl::interp_base< vec_t, vec_t >::sz
size_t sz
Vector size.
Definition: interp.h:131
o2scl::interp_krige::Kinvf
ubvector Kinvf
Inverse covariance matrix times function vector.
Definition: interp_krige.h:83
o2scl::interp_base
Base low-level interpolation class [abstract base].
Definition: interp.h:107
o2scl::interp_krige::set_covar_noise
virtual int set_covar_noise(size_t n_dim, const vec_t &x, const vec_t &y, covar_func_t &fcovar, double noise_var, bool err_on_fail=true)
Initialize interpolation routine.
Definition: interp_krige.h:143
o2scl::interp_krige_optim::mode
size_t mode
Function to minimize (default mode_loo_cv)
Definition: interp_krige.h:471
o2scl::min_brent_gsl
One-dimensional minimization using Brent's method (GSL)
Definition: min_brent_gsl.h:95

Documentation generated with Doxygen. Provided under the GNU Free Documentation License (see License Information).