dynemu_mc {dynemu}R Documentation

Predictive posterior computation via Monte Carlo (MC) approximation for one-steap-ahead Gaussian process (GP) emulators

Description

The function computes the predictive posterior distribution (mean and variance) for GP emulators using MC method, given uncertain inputs.

Usage

dynemu_mc(mean, var, dynGP, mc.sample)

Arguments

mean

numeric vector or matrix specifying the mean of uncertain inputs. The number of columns must match the input dimension of 'dynGP'.

var

numeric vector or matrix specifying the variance of uncertain inputs. The number of columns must match the input dimension of 'dynGP'.

dynGP

list of trained GP emulators fitted by dynemu_GP, each corresponding to a state variable.

mc.sample

a number of mc samples generated for MC approximation. Default is 1000.

Details

Given a trained set of GP emulators 'dynGP' fitted using dynemu_GP, this function: 1. Computes the MC-approximated predictive posterior mean and variance for each state variable. 2. Incorporates input uncertainty by integrating over the input distribution via MC sampling.

Unlike dynemu_exact, which provides a closed-form solution, this function uses MC sampling to approximate the posterior.

Value

A list containing:

References

H. Mohammadi, P. Challenor, and M. Goodfellow (2019). Emulating dynamic non-linear simulators using Gaussian processes. Computational Statistics & Data Analysis, 139, 178-196; doi:10.1016/j.csda.2019.05.006

Examples

library(lhs)
### Lotka-Volterra equations ###
LVmod0D <- function(Time, State, Pars) {
  with(as.list(c(State, Pars)), {
    IngestC <- rI * P * C
    GrowthP <- rG * P * (1 - P/K)
    MortC <- rM * C

    dP <- GrowthP - IngestC
    dC <- IngestC * AE - MortC
    return(list(c(dP, dC)))
  })
}

### Define parameters ###
pars <- c(rI = 1.5, rG = 1.5, rM = 2, AE = 1, K = 10)

### Define time sequence ###
times <- seq(0, 30, by = 0.01)

### Initial conditions ###
set.seed(1)
d <- 2
n <- 12*d
Design <- maximinLHS(n, d) * 5 # Generate LHS samples in range [0,5]
colnames(Design) <- c("P", "C")

### Fit GP emulators ###
fit.dyn <- dynemu_GP(LVmod0D, times, pars, Design)

### Define uncertain inputs ###
xmean <- c(P = 1, C = 2)
xvar <- c(P = 1e-5, C = 1e-5)

### Define number of MC samples ###
mc.sample <- 1000

### Compute the next point ###
dynemu_mc(xmean, xvar, fit.dyn, mc.sample)


[Package dynemu version 1.0.0 Index]