Title: | Uncertainty Propagation for R Vectors |
---|---|
Description: | Support for measurement errors in R vectors, matrices and arrays: automatic uncertainty propagation and reporting. Documentation about 'errors' is provided in the paper by Ucar, Pebesma & Azcorra (2018, <doi:10.32614/RJ-2018-075>), included in this package as a vignette; see 'citation("errors")' for details. |
Authors: | Iñaki Ucar [aut, cph, cre] , Lionel Henry [ctb], RStudio [cph] (Copyright for code written by RStudio employees.) |
Maintainer: | Iñaki Ucar <[email protected]> |
License: | MIT + file LICENSE |
Version: | 0.4.2.1 |
Built: | 2024-11-21 05:05:44 UTC |
Source: | https://github.com/r-quantities/errors |
Support for measurement errors in R vectors, matrices and arrays: automatic uncertainty propagation and reporting.
Every measurement has an unknown error associated. Uncertainty is the
acknowledgement of that error: we are aware that our representation of reality
may differ from reality itself. This package provides support for measurement
errors in R vectors, matrices and arrays. Uncertainty metadata is associated
to quantity values (see errors
), and this uncertainty is
automatically propagated when you operate with errors
objects (see
groupGeneric.errors
), or with errors
and numeric objects
(then numeric values are automatically coerced to errors
objects with
no uncertainty).
Correlations between measurements are also supported. In particular, any
operation (e.g., z <- x + y
) results in a correlation between output
and input variables (i.e., z
is correlated to x
and y
,
even if there was no correlation between x
and y
). And in
general, the user can establish correlations between any pair of variables
(see correl
).
This package treats uncertainty as coming from Gaussian and linear sources (note that, even for non-Gaussian non-linear sources, this is a reasonable assumption for averages of many measurements), and propagates them using the first-order Taylor series method for propagation of uncertainty. Although the above assumptions are valid in a wide range of applications in science and engineering, the practitioner should evaluate whether they apply for each particular case.
Iñaki Ucar
Iñaki Ucar, Edzer Pebesma and Arturo Azcorra (2018). Measurement Errors in R. The R Journal, 10(2), 549-557. doi:10.32614/RJ-2018-075
datasets
for a description of the datasets used in the
examples below.
## Simultaneous resistance and reactance measurements # Obtain mean values and uncertainty from measured values V <- mean(set_errors(GUM.H.2$V)) I <- mean(set_errors(GUM.H.2$I)) phi <- mean(set_errors(GUM.H.2$phi)) # Set correlations between variables correl(V, I) <- with(GUM.H.2, cor(V, I)) correl(V, phi) <- with(GUM.H.2, cor(V, phi)) correl(I, phi) <- with(GUM.H.2, cor(I, phi)) # Computation of resistance, reactance and impedance values (R <- (V / I) * cos(phi)) (X <- (V / I) * sin(phi)) (Z <- (V / I)) # Correlations between derived quantities correl(R, X) correl(R, Z) correl(X, Z) ## Calibration of a thermometer # Least-squares fit for a reference temperature of 20 degC fit <- lm(bk ~ I(tk - 20), data = GUM.H.3) # Extract coefficients and set correlation using the covariance matrix y1 <- set_errors(coef(fit)[1], sqrt(vcov(fit)[1, 1])) y2 <- set_errors(coef(fit)[2], sqrt(vcov(fit)[2, 2])) covar(y1, y2) <- vcov(fit)[1, 2] # Predicted correction for 30 degC (b.30 <- y1 + y2 * set_errors(30 - 20))
## Simultaneous resistance and reactance measurements # Obtain mean values and uncertainty from measured values V <- mean(set_errors(GUM.H.2$V)) I <- mean(set_errors(GUM.H.2$I)) phi <- mean(set_errors(GUM.H.2$phi)) # Set correlations between variables correl(V, I) <- with(GUM.H.2, cor(V, I)) correl(V, phi) <- with(GUM.H.2, cor(V, phi)) correl(I, phi) <- with(GUM.H.2, cor(I, phi)) # Computation of resistance, reactance and impedance values (R <- (V / I) * cos(phi)) (X <- (V / I) * sin(phi)) (Z <- (V / I)) # Correlations between derived quantities correl(R, X) correl(R, Z) correl(X, Z) ## Calibration of a thermometer # Least-squares fit for a reference temperature of 20 degC fit <- lm(bk ~ I(tk - 20), data = GUM.H.3) # Extract coefficients and set correlation using the covariance matrix y1 <- set_errors(coef(fit)[1], sqrt(vcov(fit)[1, 1])) y2 <- set_errors(coef(fit)[2], sqrt(vcov(fit)[2, 2])) covar(y1, y2) <- vcov(fit)[1, 2] # Predicted correction for 30 degC (b.30 <- y1 + y2 * set_errors(30 - 20))
S3 method for errors
objects (see as.data.frame
).
## S3 method for class 'errors' as.data.frame(x, row.names = NULL, optional = FALSE, ...)
## S3 method for class 'errors' as.data.frame(x, row.names = NULL, optional = FALSE, ...)
x |
any R object. |
row.names |
|
optional |
logical. If |
... |
additional arguments to be passed to or from methods. |
x <- set_errors(1:3, 0.1) y <- set_errors(4:6, 0.2) (z <- cbind(x, y)) as.data.frame(z)
x <- set_errors(1:3, 0.1) y <- set_errors(4:6, 0.2) (z <- cbind(x, y)) as.data.frame(z)
S3 method for errors
objects (see as.list
).
## S3 method for class 'errors' as.list(x, ...)
## S3 method for class 'errors' as.list(x, ...)
x |
object to be coerced or tested. |
... |
objects, possibly named. |
x <- set_errors(1:3, 0.1) as.list(x)
x <- set_errors(1:3, 0.1) as.list(x)
S3 method for errors
objects (see as.matrix
).
## S3 method for class 'errors' as.matrix(x, ...)
## S3 method for class 'errors' as.matrix(x, ...)
x |
an R object. |
... |
additional arguments to be passed to or from methods. |
as.matrix(set_errors(1:3, 0.1))
as.matrix(set_errors(1:3, 0.1))
S3 method for errors
objects (see c
).
## S3 method for class 'errors' c(..., recursive = FALSE)
## S3 method for class 'errors' c(..., recursive = FALSE)
... |
objects to be concatenated. All |
recursive |
logical. If |
c(set_errors(1, 0.2), set_errors(7:9, 0.1), 3)
c(set_errors(1, 0.2), set_errors(7:9, 0.1), 3)
S3 methods for errors
objects (see cbind
).
## S3 method for class 'errors' cbind(..., deparse.level = 1) ## S3 method for class 'errors' rbind(..., deparse.level = 1)
## S3 method for class 'errors' cbind(..., deparse.level = 1) ## S3 method for class 'errors' rbind(..., deparse.level = 1)
... |
(generalized) vectors or matrices. These can be given as named
arguments. Other R objects may be coerced as appropriate, or S4
methods may be used: see sections ‘Details’ and
‘Value’. (For the |
deparse.level |
integer controlling the construction of labels in
the case of non-matrix-like arguments (for the default method): |
x <- set_errors(1, 0.1) y <- set_errors(1:3, 0.2) (m <- cbind(x, y)) # the '1' (= shorter vector) is recycled (m <- cbind(m, 8:10)[, c(1, 3, 2)]) # insert a column cbind(y, diag(3)) # vector is subset -> warning cbind(0, rbind(x, y))
x <- set_errors(1, 0.1) y <- set_errors(1:3, 0.2) (m <- cbind(x, y)) # the '1' (= shorter vector) is recycled (m <- cbind(m, 8:10)[, c(1, 3, 2)]) # insert a column cbind(y, diag(3)) # vector is subset -> warning cbind(0, rbind(x, y))
errors
ObjectsSet or retrieve correlations or covariances between errors
objects.
See the details section below.
correl(x, y) correl(x, y) <- value set_correl(x, y, value) covar(x, y) covar(x, y) <- value set_covar(x, y, value)
correl(x, y) correl(x, y) <- value set_correl(x, y, value) covar(x, y) covar(x, y) <- value set_covar(x, y, value)
x |
an object of class |
y |
an object of class |
value |
a numeric vector of length 1 or the same length as |
The uncertainties associated to errors
objects are supposed
to be independent by default. If there is some known correlation, it can be
defined using these methods, and it will be used for the propagation of the
uncertainty by the mathematical and arithmetic operations.
The correl
method sets or retrieves correlations, i.e., a value (or
vector of values) between -1
and 1
(see base cor
on how to compute correlations). A covariance is just a correlation value
multiplied by the standard deviations (i.e., the standard uncertainty) of
both variables. It can be defined using the covar
method (see base
cov
on how to compute covariances). These methods are
equivalent; in fact, correl
calls covar
internally.
Every errors
object has a unique ID, and pairwise correlations are
stored in an internal hash table. All the functions or methods that modify
somehow the dimensions of errors
objects (i.e., subsets, binds,
concatenations, summaries...) generate new objects with new IDs, and
correlations are not, and cannot be, propagated. Only mathematical and
arithmetic operations propagate correlations, where appropriate, following
the Taylor series method.
correl
and covar
return a vector of correlations and
covariances respectively (or NULL
).
set_correl
and set_covar
, which are pipe-friendly versions of
the setters, return the x
object.
x <- set_errors(1:5, 0.1) y <- set_errors(1:5, 0.1) # Self-correlation is of course 1, and cannot be changed correl(x, x) ## Not run: correl(x, x) <- 0.5 ## End(Not run) # Cross-correlation can be set, but must be a value between -1 and 1 correl(x, y) ## Not run: correl(x, y) <- 1.5 ## End(Not run) correl(x, y) <- runif(length(x)) correl(x, y) covar(x, y)
x <- set_errors(1:5, 0.1) y <- set_errors(1:5, 0.1) # Self-correlation is of course 1, and cannot be changed correl(x, x) ## Not run: correl(x, x) <- 0.5 ## End(Not run) # Cross-correlation can be set, but must be a value between -1 and 1 correl(x, y) ## Not run: correl(x, y) <- 1.5 ## End(Not run) correl(x, y) <- runif(length(x)) correl(x, y) covar(x, y)
Datasets found in Annex H of the GUM (see reference below).
GUM.H.2 GUM.H.3
GUM.H.2 GUM.H.3
GUM.H.2
, from Section 2 of Annex H (Table H.2), provides
simultaneous resistance and reactance measurements. It is a data frame with 5
rows and 3 variables:
Voltage amplitude, in Volts.
Current amplitude, in Amperes.
Phase-shift angle of the voltage relative to the current, in radians.
GUM.H.3
, from Section 3 of Annex H (Table H.6), provides
thermometer readings and observed corrections to obtain a linear calibration
curve for some reference temperature. It is a data frame with 11 rows and
2 variables:
Thermometer reading, in Celsius degrees.
Observed correction, in Celsius degrees.
An object of class data.frame
with 11 rows and 2 columns.
BIPM, IEC, IFCC, ILAC, IUPAC, IUPAP, ISO, and OIML (2008). Evaluation of Measurement Data – Guide to the Expression of Uncertainty in Measurement, 1st edn. JCGM 100:2008. Joint Committee for Guides in Metrology. https://www.bipm.org/en/committees/jc/jcgm/publications
See errors-package
for examples.
S3 method for errors
objects (see diff
).
## S3 method for class 'errors' diff(x, lag = 1L, differences = 1L, ...)
## S3 method for class 'errors' diff(x, lag = 1L, differences = 1L, ...)
x |
a numeric vector or matrix containing the values to be differenced. |
lag |
an integer indicating which lag to use. |
differences |
an integer indicating the order of the difference. |
... |
further arguments to be passed to or from methods. |
diff(set_errors(1:10, 0.1), 2) diff(set_errors(1:10, 0.1), 2, 2) x <- cumsum(cumsum(set_errors(1:10, 0.1))) diff(x, lag = 2) diff(x, differences = 2)
diff(set_errors(1:10, 0.1), 2) diff(set_errors(1:10, 0.1), 2, 2) x <- cumsum(cumsum(set_errors(1:10, 0.1))) diff(x, lag = 2) diff(x, differences = 2)
Drop Uncertainty
drop_errors(x) ## S3 method for class 'data.frame' drop_errors(x)
drop_errors(x) ## S3 method for class 'data.frame' drop_errors(x)
x |
an |
the numeric without any errors
attributes, while preserving other
attributes like dimensions or other classes.
Equivalent to errors(x) <- NULL
or set_errors(x, NULL)
.
Set or retrieve uncertainty to/from numeric vectors.
errors(x) errors_max(x) errors_min(x) errors(x) <- value set_errors(x, value = 0) as.errors(x, value = 0)
errors(x) errors_max(x) errors_min(x) errors(x) <- value set_errors(x, value = 0) as.errors(x, value = 0)
x |
a numeric object, or object of class |
value |
a numeric vector of length 1 or the same length as |
`errors<-`
sets the uncertainty values (and converts x
into an object of class errors
). set_errors
is a pipe-friendly
version of `errors<-`
and returns an object of class errors
.
as.errors
is an alias for set_errors
.
See correl
on how to handle correlations between pairs of variables.
errors
returns a vector of uncertainty. errors_max
(errors_min
) returns the values plus (minus) the uncertainty.
groupGeneric.errors
, mean.errors
,
Extract.errors
, c
, rep
, cbind.errors
,
format.errors
, print.errors
, plot.errors
,
as.data.frame.errors
, as.matrix.errors
, t
.
x = 1:3 class(x) x errors(x) <- 0.1 class(x) x (x <- set_errors(x, seq(0.1, 0.3, 0.1))) errors_max(x) errors_min(x)
x = 1:3 class(x) x errors(x) <- 0.1 class(x) x (x <- set_errors(x, seq(0.1, 0.3, 0.1))) errors_max(x) errors_min(x)
S3 operators to extract or replace parts of errors
objects.
## S3 method for class 'errors' x[...] ## S3 method for class 'errors' x[[...]] ## S3 replacement method for class 'errors' x[...] <- value ## S3 replacement method for class 'errors' x[[...]] <- value
## S3 method for class 'errors' x[...] ## S3 method for class 'errors' x[[...]] ## S3 replacement method for class 'errors' x[...] <- value ## S3 replacement method for class 'errors' x[[...]] <- value
x |
object from which to extract element(s) or in which to replace element(s). |
... |
additional arguments to be passed to base methods
(see |
value |
typically an array-like R object of a similar class as |
x <- set_errors(1:3, 0.1) y <- set_errors(4:6, 0.2) (z <- rbind(x, y)) z[2, 2] z[2, 2] <- -1 errors(z[[1, 2]]) <- 0.8 z[, 2]
x <- set_errors(1:3, 0.1) y <- set_errors(4:6, 0.2) (z <- rbind(x, y)) z[2, 2] z[2, 2] <- -1 errors(z[[1, 2]]) <- 0.8 z[, 2]
errors
Format an errors
object for pretty printing.
## S3 method for class 'errors' format(x, digits = NULL, scientific = FALSE, notation = getOption("errors.notation", "parenthesis"), decimals = getOption("errors.decimals", FALSE), ...)
## S3 method for class 'errors' format(x, digits = NULL, scientific = FALSE, notation = getOption("errors.notation", "parenthesis"), decimals = getOption("errors.decimals", FALSE), ...)
x |
an |
digits |
how many significant digits are to be used for uncertainties.
The default, |
scientific |
logical specifying whether the elements should be encoded in scientific format. |
notation |
error notation; |
decimals |
logical specifying whether the uncertainty should be formatted
with a decimal point even when the |
... |
ignored. |
K. Nakamura et al. (Particle Data Group), J. Phys. G 37, 075021 (2010)
x <- set_errors(1:3*100, 1:3*100 * 0.05) format(x) format(x, digits=2) format(x, digits=2, decimals=TRUE) format(x, scientific=TRUE) format(x, notation="plus-minus") x <- set_errors(c(0.827, 0.827), c(0.119, 0.367)) format(x, notation="plus-minus", digits="pdg")
x <- set_errors(1:3*100, 1:3*100 * 0.05) format(x) format(x, digits=2) format(x, digits=2, decimals=TRUE) format(x, scientific=TRUE) format(x, notation="plus-minus") x <- set_errors(c(0.827, 0.827), c(0.119, 0.367)) format(x, notation="plus-minus", digits="pdg")
errors
objectsAutomatic errorbars for variables with uncertainty.
geom_errors(mapping = NULL, data = NULL, stat = "identity", position = "identity", ..., na.rm = FALSE, orientation = NA, show.legend = NA, inherit.aes = TRUE)
geom_errors(mapping = NULL, data = NULL, stat = "identity", position = "identity", ..., na.rm = FALSE, orientation = NA, show.legend = NA, inherit.aes = TRUE)
mapping |
Set of aesthetic mappings created by |
data |
The data to be displayed in this layer. There are three options: If A A |
stat |
The statistical transformation to use on the data for this layer.
When using a
|
position |
A position adjustment to use on the data for this layer. This
can be used in various ways, including to prevent overplotting and
improving the display. The
|
... |
Other arguments passed on to
|
na.rm |
If |
orientation |
The orientation of the layer. The default ( |
show.legend |
logical. Should this layer be included in the legends?
|
inherit.aes |
If |
geom_errors()
understands the following aesthetics
(required aesthetics are in bold):
x
or y
alpha
colour
group
linetype
linewidth
if (requireNamespace("ggplot2", quietly=TRUE)) { iris.e <- iris iris.e[1:4] <- lapply(iris.e[1:4], function(x) set_errors(x, x*0.02)) library(ggplot2) ggplot(iris.e) + aes(Sepal.Length, Sepal.Width, color=Species) + geom_point() + geom_errors() ggplot(iris.e) + aes(Sepal.Length, Sepal.Width, color=Species) + geom_point() + geom_errors(width=0.05, height=0.05, linewidth=0.2) }
if (requireNamespace("ggplot2", quietly=TRUE)) { iris.e <- iris iris.e[1:4] <- lapply(iris.e[1:4], function(x) set_errors(x, x*0.02)) library(ggplot2) ggplot(iris.e) + aes(Sepal.Length, Sepal.Width, color=Species) + geom_point() + geom_errors() ggplot(iris.e) + aes(Sepal.Length, Sepal.Width, color=Species) + geom_point() + geom_errors(width=0.05, height=0.05, linewidth=0.2) }
Math
, Ops
and Summary
group generic methods for
errors
objects with support for automatic uncertainty propagation (see
groupGeneric
for a comprehensive list of available methods).
## S3 method for class 'errors' Math(x, ...) ## S3 method for class 'errors' Ops(e1, e2) ## S3 method for class 'errors' Summary(..., na.rm = FALSE)
## S3 method for class 'errors' Math(x, ...) ## S3 method for class 'errors' Ops(e1, e2) ## S3 method for class 'errors' Summary(..., na.rm = FALSE)
x , e1 , e2
|
objects. |
... |
further arguments passed to methods. |
na.rm |
logical: should missing values be removed? |
Math
The sign
method returns a numeric value without uncertainty. floor
,
ceiling
, trunc
, round
and signif
add the rounding
error to the original uncertainty. lgamma
, gamma
, digamma
and
trigamma
are not implemented. The rest of the methods propagate the
uncertainty as expected from the first-order Taylor series method.
Ops
Boolean operators drop the uncertainty (showing a warning once) and operate on the
numeric values. The rest of the operators propagate the uncertainty as expected from
the first-order Taylor series method. Any numeric operand is automatically
coerced to errors
(showing a warning once) with no uncertainty.
Summary
The methods all
and any
are not supported for errors
objects and fail with an informative message. min
, max
(and
range
) return the minimum or (and) maximum value minus/plus its uncertainty.
sum
and prod
propagate the uncertainty as expected from the first-order
Taylor series method.
x <- set_errors(1:3, 0.1) exp(x) log(x) cumsum(x) cumprod(x) y <- set_errors(4:6, 0.2) x / sqrt(y) + y * sin(x) # numeric values are automatically coerced to errors x^2 # boolean operators drop uncertainty y > x c(min(x), max(x)) range(x) sum(y) prod(y)
x <- set_errors(1:3, 0.1) exp(x) log(x) cumsum(x) cumprod(x) y <- set_errors(4:6, 0.2) x / sqrt(y) + y * sin(x) # numeric values are automatically coerced to errors x^2 # boolean operators drop uncertainty y > x c(min(x), max(x)) range(x) sum(y) prod(y)
S3 methods for errors
objects.
## S3 method for class 'errors' mean(x, trim = 0, na.rm = FALSE, ...) ## S3 method for class 'errors' weighted.mean(x, ..., na.rm = FALSE) ## S3 method for class 'errors' median(x, na.rm = FALSE, ...)
## S3 method for class 'errors' mean(x, trim = 0, na.rm = FALSE, ...) ## S3 method for class 'errors' weighted.mean(x, ..., na.rm = FALSE) ## S3 method for class 'errors' median(x, na.rm = FALSE, ...)
x |
an |
trim |
the fraction (0 to 0.5) of observations to be
trimmed from each end of |
na.rm |
a logical evaluating to |
... |
further arguments passed to of from other methods. |
The mean
and weighted.mean
methods set the uncertainty as
the maximum of the standard deviation of the mean and the (weighted) mean of the uncertainty.
The median
method sets the uncertainty as 1.253 * errors(mean(x))
,
which is derived from the asymptotic variance formula of the median. Note that
this value is valid only if the sample is big enough.
An errors
object.
S3 method for errors
objects which automatically prints the error bars.
## S3 method for class 'errors' plot(x, y, ...)
## S3 method for class 'errors' plot(x, y, ...)
x , y
|
the |
... |
additional arguments (see |
cars <- as.matrix(cars) cars <- as.data.frame(set_errors(cars, cars * 0.05)) plot(cars$speed) plot(cars)
cars <- as.matrix(cars) cars <- as.data.frame(set_errors(cars, cars * 0.05)) plot(cars$speed) plot(cars)
S3 method for errors
objects.
## S3 method for class 'errors' print(x, ...)
## S3 method for class 'errors' print(x, ...)
x |
an |
... |
further arguments passed to or from other methods. |
x <- set_errors(1:10, 1:10 * 0.05) print(x) print(x[1:3]) print(x[1]) print(x[1], digits=2) print(x[1], notation="plus-minus")
x <- set_errors(1:10, 1:10 * 0.05) print(x) print(x[1:3]) print(x[1]) print(x[1], digits=2) print(x[1], notation="plus-minus")
S3 method for errors
objects (see rep
).
## S3 method for class 'errors' rep(x, ...)
## S3 method for class 'errors' rep(x, ...)
x |
a vector (of any mode including a |
... |
further arguments to be passed to or from other methods. For the internal default method these can include:
|
rep(set_errors(1, 0.1), 4)
rep(set_errors(1, 0.1), 4)
S3 method for errors
objects (see t
).
## S3 method for class 'errors' t(x)
## S3 method for class 'errors' t(x)
x |
a matrix or data frame, typically. |
a <- matrix(1:30, 5, 6) errors(a) <- 1:30 t(a)
a <- matrix(1:30, 5, 6) errors(a) <- 1:30 t(a)