Procedure: Multivariate lognormal approximation for correlation hyperparameters

Description and Background

The posterior distribution \(\pi^*(\delta)\) of the hyperparameter vector \(\delta\) is given in the procedure page for building a Gaussian process emulator for the core problem (ProcBuildCoreGP) in the case of the core problem (with linear mean function and weak prior information). The ideal way to compute the emulator in this case is to generate a sample of values of \(\delta\) from this distribution, but that is itself a complex computational task. We present here a simpler method based on a lognormal approximation to \(\pi^*(\delta)\).

Inputs

  • An emulator as defined in ProcBuildCoreGP, using a linear mean and a weak prior.
  • The mode of the posterior \(\hat{\delta}\) as defined in the discussion page on finding the posterior mode of correlation lengths (DiscPostModeDelta).
  • The \(p \times p\) Hessian matrix \(\displaystyle \frac{\partial ^2 g(\tau)}{\partial \tau^2}\) with (k,l)-th entry \(\displaystyle \frac{\partial ^2 g(\tau)}{\partial \tau_l \partial \tau_k}\), as defined in DiscPostModeDelta.

Outputs

  • A set of \(s\) samples for the correlation lengths \(\delta\), denoted as \(\tilde{\delta}\).
  • A posterior mean \(\tilde{m}^*(\cdot)\) and a covariance function \(\tilde{u}^*(\cdot,\cdot)\), conditioned on the samples \(\tilde{\delta}\)

Procedure

  • Define \(\displaystyle V = -\left(\frac{\partial ^2 g(\tau)}{\partial \tau^2}\right)^{-1}\). Draw \(s\) samples from the \(p\)-variate normal distribution \({\cal N} (\hat{\tau},V)\), call these samples \(\tilde{\tau}\).
  • Calculate the samples \(\tilde{\delta}\) as \(\tilde{\delta} = \exp(\tilde{\tau}/2)\).
  • Given the set of \(s\) samples \(\tilde{\delta}\), the posterior mean and variance \(\tilde{m}^*(\cdot)\), \(\tilde{u}^*(\cdot,\cdot)\) can be calculated with the same formulae given in the procedure page for sampling the posterior distribution of the correlation lengths (ProcMCMCDeltaCoreGP), or in more detail in the procedure page for predicting the simulator’s outputs using a GP emulator (ProcPredictGP).

Additional Comments

Most standard statistical computing packages have facilities for taking random samples from a multivariate normal distribution.

When an input is not particularly active, the posterior distribution of the correlation lengths \(\pi^*_{\delta}(\delta)\) can be very flat with respect to that input and obtain its maximum for a large value of \(\delta\). This can cause the respective entry of the matrix \(V\) to be very large, and the samples \(\tilde{\delta}\) that correspond to this input to have both unrealistically large and small values. An inspection of the samples that are returned by the above procedure is recommended, especially in high dimensional input problems, where less active inputs are likely to exist. If the sampled correlation lengths that correspond to one or more inputs are found to have very large (e.g. >50) and very small (e.g. \(<0.5\)) values at the same time, a potential remedy could be to fix the values of these samples to the respective entries of the \(\hat{\delta}\) vector.

References

This method is introduced in the following report.

  • Nagy B., Loeppky J.L. and Welch W.J. (2007). Fast Bayesian Inference for Gaussian Process Models. Technical Report 230, Department of Statistics, University of British Columbia.

Note however that the authors indicate also that the method works well when the correlation function has the Gaussian form but may not work so well in the case of the exponential power form (see the alternatives page on emulator prior correlation function (AltCorrelationFunction):

  • Nagy B., Loeppky J.L. and Welch W.J. (2007). Correlation parameterization in random function models to improve normal approximation of the likelihood or posterior. Technical Report 229, Department of Statistics, University of British Columbia.