ForK Library
Functions/Subroutines

kriginggradvariance.f90 File Reference

Subroutine to compute the variance of derivative values from Kriging Surface (requires buildkriging to be called first) More...

Go to the source code of this file.

Functions/Subroutines

subroutine kriginggradvariance (ndim, ntot, X, stot, H, hyper, mtot, Xm, gtotm, ptsm, dimsm, Gm, Var, covarflagi)
 In order to provide a confidence measure for the derivative predictions associated with the Kriging surface, the variance of the derivative prediction can be calculated based on the covariance structure of the Kriging model. To predict the variance of derivative values, the variance between derivative and function values as well as derivative and derivative values must be determined. Due to the linearity of differentiation, these covariance values are determined by differentiation of the covariance function. As seen previously, the covariance between a derivative value and function value is given as:

Detailed Description

Subroutine to compute the variance of derivative values from Kriging Surface (requires buildkriging to be called first)

Definition in file kriginggradvariance.f90.


Function Documentation

subroutine kriginggradvariance ( 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,
integer,intent(in)  gtotm,
integer,dimension(gtotm),intent(in)  ptsm,
integer,dimension(gtotm),intent(in)  dimsm,
real(8),dimension(stot,gtotm),intent(in)  Gm,
real(8),dimension(gtotm),intent(out)  Var,
integer,intent(in)  covarflagi 
)

In order to provide a confidence measure for the derivative predictions associated with the Kriging surface, the variance of the derivative prediction can be calculated based on the covariance structure of the Kriging model. To predict the variance of derivative values, the variance between derivative and function values as well as derivative and derivative values must be determined. Due to the linearity of differentiation, these covariance values are determined by differentiation of the covariance function. As seen previously, the covariance between a derivative value and function value is given as:

\[ cov(\frac{\partial y_{i}}{\partial x_{k}}, y_{j}) = \frac{\partial}{\partial x_{k}} k(\vec{x}_{i},\vec{x}_{j}) \]

while the correlation between derivative values is given by:

\[ cov(\frac{\partial y_{i}}{\partial x_{k}}, \frac{\partial y_{j}}{\partial x_{l}}) = \frac{\partial^{2}}{\partial x_{k} \partial x_{l}} k(\vec{x}_{i},\vec{x}_{j}) \]

. Using these definitions, the variance between points can be calculated as:

\[ V\left[\frac{\partial y(\vec{x}_{*})}{\partial x_{k}} \right]= \frac{\partial^{2}}{\partial x_{k}^{2}} k(\vec{x}_{*},\vec{x}_{*})-l_{*,k}^T K^{-1} l_{*,k}+ S_{k}(\vec{x}_*) A^{-1} S_{k}(\vec{x}_*)^{T} \]

where \(l_{*,k}^{T}\) is a vector with elements representing the covariance between the derivative value with respect to \(x_{k}\) and the training point function values, given as:

\[ l_{*,k}^{T} = \left[\frac{\partial}{\partial x_{k}} k(\vec{x}_{*},\vec{X}_{1}) \dots \frac{\partial}{\partial x_{k}} k(\vec{x}_{*},\vec{X}_{N}) \right] \]

and \(S_{k}(\vec{x}_{*})\) represents the derivative of the mean function defined as:

\[ S_{k}(\vec{x}_{*}) = g_{k}(\vec{x}_{*}) - l_{*,k}^T K^{-1} H \]

where \(g_{k}\) is a vector containing the derivative of the regression basis functions with respect to \(x_{k}\) evaluated at the test point \(\vec{x}_{*}\) and \(G_{k}\) is a matrix containing the derivative of the regression basis evaluated at the training points \(\vec{X}\). Using this variance prediction, confidence bounds can be placed on the derivative predictions from the Kriging Model.
Because the variance of the derivative prediction requires the covariance function to be twice differentiable, only the Matern functions with \( \nu \geq 3/2 \) can be used.

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) beta: Regression coefficients based on the optimal estimate for the Kriging model (size=[stot])
Supplied by buildkriging subroutine
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)gtotm: The total number of derivatives to be predicted (ndim*mtot if the gradient at every test point is desired)
in) ptsm: A vector identifying the test point where a particular derivative is specified (size=[gtotm] with values ranging from 1 to mtot)
in) dimsm: A vector identifying the dimension the derivative is taken with respect to (size=[gtotm] with values ranging from 1 to ndim)
in)Gm : The derivative of the regression basis evaluated at the test points (size=[stotxgtotm])
out)Var: The variance of the derivative values (size=[gtotm])
Because the variance prediction requires the covariance matrix of the training data, it is more expensive than the mean prediction
If possible, predict the variance associated with multiple points using a single subroutine call
in)covarflagi: Flag to govern which covariance function is used

  • 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. When the variance of the derivative predictions is desired, a twice differentiable covariance function is required. This requires \(\nu \geq 3/2 \). Must supply the same covariance flag as used in buildkriging.

Definition at line 67 of file kriginggradvariance.f90.

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

Referenced by krigingwrapper().

 All Classes Namespaces Files Functions Variables