Allows the user to specify the covariance structure for a Bayesian mass balance, simulates draws from reconciled masses and relevant covariance matrix, and approximates the log-marginal likelihood.
BMB(
X,
y,
cov.structure = c("indep", "component", "location"),
priors = "default",
BTE = c(500, 20000, 1),
lml = FALSE,
ybal = TRUE,
diagnostics = TRUE,
verb = 1
)
A matrix that maps constrained masses to observed masses. Can be built from the function constrainProcess
, see documentation for details.
A list of matrices of observed mass flow rates. Each matrix is a separate sample component. The rows of each matrix index the sampling location, and the columns index the sample set number. Can be specified using the importObservations
function.
Character string. "indep"
allows for no correlation. "component"
indicates correlation within an individual sample component. "location"
indicates correlation within an individual sampling location. Not specifying cov.structure
defaults to the "indep"
structure.
List or character string. When the default value priors = "default"
is used, the BMB
uses a set of default conjugate priors. When passing a list to the argument, the list must contain user specified hyperparameter values for each conjugate prior. To see the required list structure run BMB
with BTE = c(1,2,1)
and inspect the output. When priors = "Jeffreys"
the Jeffreys priors for \(\Sigma\) and \(\sigma^2\) with known mean given in priorlistBayesMassBal. When priors = "Jeffreys"
, the prior used for \(\beta\) is proportional to the indicator function \(I\lbrack \beta > 0 \rbrack\). See Details for more information.
Numeric vector giving c(Burn-in, Total-iterations, and Every)
for MCMC approximation of target distributions. The function BMB
produces a total number of samples of \((T - B)/E\). \(E\) specifies that only one of every \(E\) draws are saved. \(E > 1\) reduces autocorrelation between obtained samples at the expense of computation time.
Logical indicating if the log-marginal likelihood should be approximated. Default is FALSE
, which reduces computation time. Log-marginal likelihood is approximated using methods in chibBayesMassBal.
Logical indicating if the mass balanced samples for each \(y\) should be returned. Default is TRUE
. Setting ybal=FALSE
results in a savings in RAM and computation time.
Logical or list indicating if diagnostic functions geweke.diag
and effectiveSize
codaBayesMassBal should computed for the obtained samples. The default of TRUE
indicates diagnostics should be run with their default parameters. Alternatively, passing a list of the structure list(frac1 = 0.1, frac2 = 0.5)
will run both diagnostics and allow geweke.diag
to be run with parameters other than the default.
Numeric indicating verbosity of progress printed to R-console. The default of 1 prints messages and a progress bar to the console during all iterative methods. verb = 0
indicates no messages are printed.
Returns a list of outputs
beta
List of matrices of samples from the distribution of reconciled data. Each matrix in the list is a separate sample component. Each column of a matrix in beta
is a draw from the target distribution.
Sig
List of matrices containing draws from the distribution of each covariance matrix. If S.t
is the \(t^{th}\) draw from the distribution of covariance matrix S
and:
cov.structure = "component"
or "location"
, the \(t^{th}\) column of a matrix in Sig
is equal to S.t[upper.tri(S.t, diag = TRUE)]
.
priors
List of prior hyperparameters used in generating conditional posterior distributions and approximating log-marginal likelihood. The structure of the input argument priors
is required to be the same as the structure of this returned list slice. See Details.
cov.structure
Character string containing the covariance structure used.
y.cov
List of character matrices indicating details for the structure of each covariance matrix. Only returned when cov.structure = "location"
lml
Numeric of the log-marginal likelihood approximation. Returns NA
when lml = FALSE
diagnostics
List containing results from diagnostic functions geweke.diag
and effectiveSize
ybal
List of samples from the distribution of reconciled mass flow rates, in the same format as the function argument y
. Produced with argument ybal = TRUE
. Equivalent to lapply(BMB(...)$beta,function(X,x){x %*% X} , x = X)
. Viewing this output is more intuitive than viewing samples of beta
, at the expense of RAM and some computation time.
X
The function argument X
is passed to the output so that it can be used with other BayesMassBal
functions.
type
Character string used by plot.BayesMassBal
. type = "BMB"
for an object returned from the BMB
function.
See vignette("Two_Node_Process", package = "BayesMassBal")
for further details on how to use function outputs.
When the priors
argument is left unspecified, a set of default conjugate priors are used, which are chosen to allow BMB()
to work well in a general setting. In the current version of the BayesMassBal
package, only the conjugate priors stated below can be used, but hyperparameter values can be specified by the user using the priors
argument.
The prior distribution on beta
is a normal distribution truncated at 0. The mean of this distribution before truncation is the ordinary least squares (OLS) estimate of \(\beta\). OLS estimates less than 0, are changed to 0. The prior variance, before truncation, of each element of \(\beta\) is set to:
$$10^{\mathrm{number of integer digits of an element of } \beta + 6}$$
Currently, there is only support for diagonal prior covariance matrices for \(\beta\)
When cov.structure = "indep"
the error of all observations in a sample set are independent. An inverse gamma prior distribution, with \(\alpha_0 = 0.000001\) and \(\beta_0 = 0.000001\), is placed on the variance of the mass flow rate for each sample component at each sample location.
When cov.structure = "component"
or "location"
, the prior distribution on \(\Sigma_i\) is inverse Wishart \((\nu_0, \nu_0 \times S_0)\). The degrees of freedom parameter, \(\nu_0\), is equal to the dimension of \(\Sigma_i\). The scale matrix parameter is equal to a matrix, \(S_0\), with the sample variance of the relevant observation on the diagonal, multiplied by \(\nu_0\).
The user is able to specify the prior hyperparameters of the mean and variance of beta
, \(\alpha_0\) and \(\beta_0\) for each \(\sigma^2\), and the degrees of freedom and scale matrix for each \(\Sigma_i\) using the priors
argument. It is advisable for the user to specify their own prior hyperparameters for \(p(\sigma^2)\) if the variance of any element is well under 1, or \(p(\beta)\) if the there is a wide range in the magnitude of observations.
When priors = "Jeffreys"
Jeffreys priors are used for the prior distribution of the variance and covariance parameters. Priors used are \(p(\sigma^2) \propto \frac{1}{\sigma^2}\) and \(p(\Sigma) \propto |\Sigma|^{-(p+1)/2}\), as listed in priorlistBayesMassBal. The Jeffreys prior for a \(\beta\) with infinite support is \(p(\beta) \propto 1\). To preserve the prior information that \(\beta > 0\), \(p(\beta)\propto I\lbrack \beta > 0 \rbrack\) is chosen. It is not possible to calculate log-marginal likelihood using the methods in chibBayesMassBal with Jeffreys priors. Therefore, if priors = "Jeffreys"
and lml = TRUE
, the lml
argument will be ignored and a warning will be printed.
lml
is reported in base \(e\). See here for some guidance on how to interpret Bayes Factors, but note log base 10 is used on Wikipedia.
chibBayesMassBal gibbsexplBayesMassBal codaBayesMassBal priorlistBayesMassBal
# NOT RUN {
y <- importObservations(file = system.file("extdata", "twonode_example.csv",
package = "BayesMassBal"),
header = TRUE, csv.params = list(sep = ";"))
C <- matrix(c(1,-1,0,-1,0,0,1,-1,0,-1), byrow = TRUE, ncol = 5, nrow = 2)
X <- constrainProcess(C = C)
BMB_example <- BMB(X = X, y = y, cov.structure = "indep",
BTE = c(10,300,1), lml = FALSE, verb=0)
summary(BMB_example)
# }
Run the code above in your browser using DataLab