pkern_sim: Random draw from multivariate normal distribution for grids

In deankoch/pkern: Fast Kriging with Product Kernels

Random draw from multivariate normal distribution for grids

Description

Generates a random draw from the multivariate Gaussian distribution for the covariance model pars on grid g, with mean zero.

Usage

pkern_sim(g, pars = pkern_pars(g), fac = NULL)

Arguments

g	any object accepted or returned by pkern_grid
pars	list, covariance parameters in form returned by pkern_pars
fac	list, optional pre-computed factorization of component correlation matrices

Details

pars and g define the model's covariance matrix V. This function uses base::rnorm to get a vector of independent standard normal variates, which it multiplies by the square root of the covariance matrix, V, for the desired model (as defined by pars and g). The result has a multivariate normal distribution with mean zero and covariance V.

Multiple independent draws can be computed more efficiently by reusing the factorization of V. This can be pre-computed with pkern_var and supplied in fac, or a multi-layer g can be supplied (see examples).

Value

numeric vector, the vectorized grid data

Examples

# example grid and covariance parameters
gdim = c(100, 200)
g = pkern_grid(gdim)
pars_gau = pkern_pars(g)

# this example has a large nugget effect
gval = pkern_sim(g, pars=pars_gau)
pkern_plot(matrix(gval, gdim))

# plot with yx coordinates
g_sim = modifyList(g, list(gval=gval))
pkern_plot(g_sim)

# repeat with smaller nugget effect for less noisy data
pars_smooth = modifyList(pars_gau, list(eps=1e-2))
gval_smooth = pkern_sim(g, pars_smooth)
g_sim_smooth = modifyList(g, list(gval=gval_smooth))
pkern_plot(g_sim_smooth)

# the nugget effect can be very small, but users should avoid eps=0
pars_smoother = modifyList(pars_gau, list(eps=1e-12))
gval_smoother = pkern_sim(g, pars_smoother)
g_sim_smoother = modifyList(g, list(gval=gval_smoother))
pkern_plot(g_sim_smoother)

# multi-layer example
n_pt = prod(gdim)
n_layer = 3
g_multi = pkern_grid(list(gdim=gdim, gval=matrix(NA, n_pt, n_layer)))
gval_multi = pkern_sim(g_multi, pars_smoother)
g_sim_multi = modifyList(g, list(gval=gval_multi))
pkern_plot(g_sim_multi, layer=1)
pkern_plot(g_sim_multi, layer=2)
pkern_plot(g_sim_multi, layer=3)

deankoch/pkern documentation built on Oct. 26, 2023, 8:54 p.m.