glsalgo: Generalized Least Squares fitting

In gschnabel/nucdataBaynet: Bayesian networks for nuclear data evaluation

Generalized Least Squares fitting

Description

Estimate the variables in a Bayesian network via the Generalized Least Squares method.

Usage

glsalgo(
  map,
  zprior,
  U,
  obs,
  zref = zprior,
  adjust_idcs = NULL,
  damp = 0,
  ret.list = FALSE
)

Arguments

map	Mapping object. Usually a compound map, see create_compound_map
zprior	Vector of prior estimates of the independent variables (i.e., associated with nodes without parent nodes)
U	Prior covariance matrix of the independent variables
obs	Vector with observed values of dependent nodes. Must be of same size as zprior. An NA value in this vector means that the corresponding variable was not observed.
zref	Values of the independent variable that should be used as reference point to calculate the Taylor approximation (using the jacobian function of map)
adjust_idcs	Indices of variables that should be adjusted, all other independent variables are assumed to be fixed.
damp	A value indicating the magnitude of the damping term applied to the inverse posterior covariance matrix. For a pure GLS update, leave it zero. This parameter is used by the in the case of a non-linear mapping.
ret.list	If FALSE, return a vector with the estimated variables of the independent variables. Otherwise, return a list with more information. Default is FALSE

Value

Return a vector with the posterior estimates of the independent variables if ret.list=FALSE. Otherwise return a list with the following fields:

zpost	Vector of posterior estimates of the independent variables
ypost	Vector of posterior estimates of the dependent variables

See Also

LMalgo

Examples

library(Matrix)
params <- list(
  mapname = "mymap",
  maptype = "linearinterpol_map",
  src_idx = 1:10,
  tar_idx = 11:15,
  src_x = 1:10,
  tar_x = 3:7
)
mymap <- create_linearinterpol_map()
mymap$setup(params)

U <- Diagonal(n=15, x=c(rep(1e3, 10), rep(1, 5)))
zprior <- rep(0, 15)
zref <- rep(0, 15)
obs <- c(rep(NA,10), 5:9)

# glsalgo only works for linear relationships
# LMalgo can also deal with non-linear relationships
zpost <- glsalgo(mymap, zprior, U, obs)
optres <- LMalgo(mymap, zprior, U, obs)
zpost2 <- optres$zpost

# posterior estimates of values on computational grid
zpost[1:10]
# posterior estimates of error variables associated with observations
zpost[11:15]
# get posterior covariance block of independent variables
get_posterior_cov(mymap, zpost, U, obs, 1:5, 5:10)
# draw samples of independent variables from posterior distribution
get_posterior_sample(mymap, zpost, U, obs, 10)

gschnabel/nucdataBaynet documentation built on Feb. 3, 2023, 4:13 a.m.