uncertMC estimates measurement uncertainty from a function,
expression or formula by Monte Carlo simulation.
uncertMC(expr, x, u, method = "MC", df, cor, cov, distrib, distrib.pars,
B = 200, keep.x = TRUE, vectorized=TRUE, ...)An object of class uncertMC. See uncertMC-class for details.
An expression, function, or formula with no left-hand side (e.g.
~a*x+b*x^2) which can be evaluated in the environment x to
provide a numeric value.
A named list or vector of parameters supplied to expr.
A named list or named vector of length length(x) of standard uncertainties.
Method of uncertainty evaluation. The only method currently supported
by uncertMC is "MC". If any other method is specified, control is
passed to uncert.
A named list or named vector of degrees of freedom. df can be
a partial named list if not all distributions (see below) use degrees of
freedom.
Optional (square, symmetric) correlation or covariance matrices, respectively.
If neither is specified, uncertMC assumes independent variables.
A character vector of length length(x) or a named list
of names of distribution functions associated with u. See Details
for defaults.
A named list of lists of parameters describing the distributions
associated with u to be passed to the relevant distribution function.
If distrib is present but distrib.pars is not, an attempt is made
to set defaults based on other parameters; see Details.
Number of Monte Carlo replicates.
If TRUE, the simulated replicates of x are included in
the return object.
If TRUE, expr is assumed to take vector arguments. If FALSE,
expr is treated as if it takes scalar arguments. See Details for the difference.
Additional parameters to be passed to a function (for the function method) or used in an expression (for expression or formula method).
S. L. R. Ellison s.ellison@lgcgroup.com
Although most likely to be called by uncert, uncertMC may be called directly.
If any of x, u, df, distrib or distrib.pars are not lists,
they are coerced to lists. If x is not named, arbitrary names of the form 'Xn'
are applied. If u, df, distrib or distrib.pars do not have
names, the names will be set to names(x) if they are of length exactly
length(x); if not, an error is returned.
For Monte Carlo evaluation, distributions and distribution parameters are needed but
defaults are used if some or all are absent. If distrib is missing, or
if it is a list with some members missing, the distribution is assumed Normal
and any missing member of distrib is set to "norm".
Distributions are usually identified by the root of the distribution function name; for example
to specify the Normal, distrib$name="norm". At present, only the random value
generator (e.g. rnorm) is used. Names of user-specified distributions functions can also be
used, provided they have a random value generator named r<dist> where <dist>
is the abbreviated distribution. Parameters are passed to distribution functions using
do.call, so the function must accept the parameters supplied in distrib.pars.
If distrib.pars or members of it are missing, an attempt is made to deduce
appropriate distribution parameters from x, u, df and distrib.
In doing so, the following assumptions and values apply for the respective distributions:
mean=x$name, sd=u$name.
min=x-sqrt(3)*u, max=x+sqrt(3)*u.
min=x-sqrt(6)*u, max=x+sqrt(6)*u, mode=x.
df=df, mean=x, sd=u.
If either cor or cov are present, a test is made to see if off-diagonal
elements are significant. If not, uncertMC treats the values as independent.
The test simply checks whether the sum of off-diagonal elements of cor (calculated
from cov if cov is present) is bigger than
.Machine.double.eps*nrow^2.
Correlation is supported as long as all correlated variables are normally distributed.
If correlation is present, uncertMC follows a two-stage simulation procedure.
First, variables showing correlation are identified. Following a check that
their associated distrib values are all "norm", mvrnorm from
the MASS library is called to generate the simulated x values for those variables.
Second, any remaining (i.e. independent) variables are simulated from their respective
distrib and distrib.pars.
Vectorisation makes a difference to execution speed. If vectorize=TRUE, MC evaluation
uses eval using the simulated data as the evaluation environment; if not, apply
is used row-wise on the simulated input matrix. This makes an appreciable difference to
execution speed (typically eval is faster by a factor of 5 or more) so the default
assumes vectorised expressions. However, not all functions and expressions take vector arguments,
especially user functions involving complicated arithmetic or numerical solutions. Use vectorize=FALSE
for functions or expressions that do not take vector arguments.
Note: One common symptom of an expression that does not take vector arguments is
an R warning indicating that only the first element (typically of a parameter in x) is used.
uncertMC may also return NA for u on attempting to take the sd of a single simulated point.
JCGM 100 (2008) Evaluation of measurement data - Guide to the expression of uncertainty in measurement. tools:::Rd_expr_doi("10.59161/JCGM100-2008E"). (JCGM 100:2008 is a public domain copy of ISO/IEC Guide to the expression of uncertainty in measurement (1995) ).
Kragten, J. (1994) Calculating standard deviations and confidence intervals with a universally applicable spreadsheet technique, Analyst, 119, 2161-2166.
Ellison, S. L. R. (2005) Including correlation effects in an improved spreadsheet calculation of combined standard uncertainties, Accred. Qual. Assur. 10, 338-343.
uncert, uncert-class, uncertMC-class
expr <- expression(a+b*2+c*3+d/2)
x <- list(a=1, b=3, c=2, d=11)
u <- lapply(x, function(x) x/10)
u.MC<-uncertMC(expr, x, u, distrib=rep("norm", 4), method="MC")
print(u.MC, simplify=FALSE)
#An example with correlation
u.cor<-diag(1,4)
u.cor[3,4]<-u.cor[4,3]<-0.5
u.formc.MC<-uncertMC(~a+b*2+c*3+d/2, x, u, cor=u.cor, keep.x=TRUE)
u.formc.MC
#A non-linear example
expr <- expression(a/(b-c))
x <- list(a=1, b=3, c=2)
u <- lapply(x, function(x) x/20)
set.seed(403)
u.invexpr<-uncertMC(expr, x, u, distrib=rep("norm", 3), B=999, keep.x=TRUE )
u.invexpr
#Look at effect of vectorize
system.time(uncertMC(expr, x, u, distrib=rep("norm", 3), B=9999, keep.x=TRUE ))
system.time(uncertMC(expr, x, u, distrib=rep("norm", 3), B=9999, keep.x=TRUE, vectorize=FALSE))
Run the code above in your browser using DataLab