ForK Library
Functions/Subroutines

krigingfuncvariance.f90 File Reference

Subroutine to Predict variance of function values from Kriging Surface (requires buildkriging to be called first) More...

Go to the source code of this file.

Functions/Subroutines

subroutine krigingfuncvariance (ndim, ntot, X, stot, H, hyper, mtot, Xm, Hm, S, covarflagi)
 Because of the stochastic nature of a Kriging model, the variance associated with mean function values can be predicted throughout the domain. The variance for function predictions are given as:

\[ V\left[y_*\right]= cov(\vec{x}_{*},\vec{x}_{*})-k_*^T K^{-1} k_*+ R(\vec{x}_*) A^{-1} R(\vec{x}_*)^{T} \]



Detailed Description

Subroutine to Predict variance of function values from Kriging Surface (requires buildkriging to be called first)

Definition in file krigingfuncvariance.f90.


Function Documentation

subroutine krigingfuncvariance ( integer,intent(in)  ndim,
integer,intent(in)  ntot,
real(8),dimension(ndim,ntot),intent(in)  X,
integer,intent(in)  stot,
real(8),dimension(stot,ntot),intent(in)  H,
real(8),dimension(ndim+2),intent(in)  hyper,
integer,intent(in)  mtot,
real(8),dimension(ndim,mtot),intent(in)  Xm,
real(8),dimension(stot,mtot),intent(in)  Hm,
real(8),dimension(mtot),intent(out)  S,
integer,intent(in)  covarflagi 
)

Because of the stochastic nature of a Kriging model, the variance associated with mean function values can be predicted throughout the domain. The variance for function predictions are given as:

\[ V\left[y_*\right]= cov(\vec{x}_{*},\vec{x}_{*})-k_*^T K^{-1} k_*+ R(\vec{x}_*) A^{-1} R(\vec{x}_*)^{T} \]

The first two terms, \(cov(\vec{x}_{*},\vec{x}_{*})-k_*^T K^{-1} k_*\), account for the variance of associated with a zero mean gaussian process while the last term \(R(\vec{x}_*) A^{-1} R(\vec{x}_*)^{T}\) accounts for the variance associated with the uncertainty associated with the regression parameters, where \( R(\vec{x}_*)\) is given as:

\[ R(\vec{x}_{*}) = h^{T} - k_*^T K^{-1} H^{T} \]


To compute the function variance, the hyperparameters within the covariance function are required (supplied by buildkriging). Hence, these subroutines must be used in conjunction with buildkriging called first.
Using the variance, a confidence interval can be established on the function predictions. For example, for a well training model, \( 95 \%\) of the true function values should lie within 2 standard deviations from the mean prediction, represented functionally as:

\begin{eqnarray*} y_{lower} &=& y_{*} - 2 \sqrt{V\left[y_*\right]} \\ y_{upper} &=& y_{*} + 2 \sqrt{V\left[y_*\right]} \end{eqnarray*}

Using this variance prediction, predictions from the model can be bounded. Additionally, the quality of the Kriging model can be tested using cross-validation procedures.

Author:
Brian Lockwood
Department of Mechanical Engineering
University of Wyoming
Date:
May 2, 2012
Parameters:
in)ndim : The dimension of the problem
in)ntot : The number of Training points
in)X : The location of the training points (size=[ndimxntot])
in)stot : Number of Terms in the regression
in)H: The collocation matrix for the regression (size=[stotxntot])
in)hyper: Hyperparameters for the Kriging Model (size=[ndim+2])
Supplied by buildkriging subroutine
in)mtot : The number of test points, the places where function prediction are desired
in)Xm : The location of the test points (size=[ndimxmtot])
in)Hm : The collocation matrix evaluated at the test points (size=[stotxmtot])
out)S: The variance of the function values (size=[mtot])
Predicting the variance values requires the construction of the covariance matrix for the training data so it is more expensive than function predictions
Best to predict the variance associated with multiple points using a single function call if possible.
in)covarflagi: Flag to govern which covariance function is used

  • covarflag==0 Uses Matern function with \(\nu=1/2\)
  • covarflag==1 Uses Matern function with \(\nu=3/2\)
  • covarflag==2 Uses Matern function with \(\nu=5/2\)


The parameter \(\nu\) governs the smoothness and differentiability of the covariance function. For function only Kriging, all three options are available.
Must supply the same covariance flag as used in buildkriging.

Definition at line 52 of file krigingfuncvariance.f90.

References covars::covarflag, covars::covarfunc(), choleskymod::invertsymmetric(), choleskymod::matrixmulttrans(), choleskymod::symmatrixmulttrans(), and choleskymod::symmatvec().

Referenced by funcmod::func_kriging(), and krigingwrapper().

 All Classes Namespaces Files Functions Variables