interval {propagate}R Documentation

Uncertainty propagation based on interval arithmetics

Description

Calculates the uncertainty of a model by using interval arithmetics based on a "combinatorial sequence grid evaluation" approach, thereby avoiding the classical dependency problem that inflates the result interval.

Usage

interval(df, expr, seq = 10, plot = TRUE)

Arguments

df

a 2-row dataframe/matrix with lower border values AiA_i in the first row and upper border values BiB_i in the second row. Column names must correspond to the variable names in expr.

expr

an expression, such as expression(x/y).

seq

the sequence length from AiA_i to BiB_i in [Ai,Bi][A_i, B_i].

plot

logical. If TRUE, plots the evaluations and min/max values as blue border lines.

Details

For two variables x,y{\color{red}x}, {\color{blue}y} with intervals [x1,x2][{\color{red}x_1}, {\color{red}x_2}] and [y1,y2][{\color{blue}y_1}, {\color{blue}y_2}], the four basic arithmetic operations op{+,,,/}\langle \mathrm{op} \rangle \in \{+, -, \cdot, /\} are

[x1,x2] ⁣op ⁣[y1,y2]=[{\color{red}x_1}, {\color{red}x_2}] \,\langle\!\mathrm{op}\!\rangle\, [{\color{blue}y_1}, {\color{blue}y_2}] =

[min(x1 ⁣op ⁣y1,x1 ⁣op ⁣y2,x2 ⁣op ⁣y1,x2 ⁣op ⁣y2),max(x1 ⁣op ⁣y1,x1 ⁣op ⁣y2,x2 ⁣op ⁣y1,x2 ⁣op ⁣y2)]\left[ \min({\color{red}x_1} {\langle\!\mathrm{op}\!\rangle} {\color{blue}y_1}, {\color{red}x_1} \langle\!\mathrm{op}\!\rangle {\color{blue}y_2}, {\color{red}x_2} \langle\!\mathrm{op}\!\rangle {\color{blue}y_1}, {\color{red}x_2} \langle\!\mathrm{op}\!\rangle {\color{blue}y_2}), \max({\color{red}x_1} {\langle\!\mathrm{op}\!\rangle} {\color{blue}y_1}, {\color{red}x_1} {\langle\!\mathrm{op}\!\rangle} {\color{blue}y_2}, {\color{red}x_2} {\langle\!\mathrm{op}\!\rangle} {\color{blue}y_1}, {\color{red}x_2} {\langle\!\mathrm{op}\!\rangle} {\color{blue}y_2})\right]

So for a function f([x1,x2],[y1,y2],[z1,z2],...)f([{\color{red}x_1}, {\color{red}x_2}], [{\color{blue}y_1}, {\color{blue}y_2}], [{\color{green}z_1}, {\color{green}z_2}], ...) with kk variables, we have to create all combinations Ci=({{x1,x2},{y1,y2},{z1,z2},...}k)C_i = {\{\{{\color{red}x_1}, {\color{red}x_2}\}, \{{\color{blue}y_1}, {\color{blue}y2}\}, \{{\color{green}z_1}, {\color{green}z_2}\}, ...\} \choose k}, evaluate their function values Ri=f(Ci)R_i = f(C_i) and select R=[minRi,maxRi]R = [\min R_i, \max R_i].
The so-called dependency problem is a major obstacle to the application of interval arithmetic and arises when the same variable exists in several terms of a complicated and often nonlinear function. In these cases, over-estimation can cover a range that is significantly larger, i.e. minRiminf(x,y,z,...),maxRimaxf(x,y,z,...)\min R_i \ll \min f(x, y, z, ...) , \max R_i \gg \max f(x, y, z, ...). For an example, see http://en.wikipedia.org/w/index.php?title=Interval_arithmetic under "Dependency problem". A partial solution to this problem is to refine RiR_i by dividing [x1,x2][{\color{red}x_1}, {\color{red}x_2}] into ii smaller subranges to obtain sequence (x1,x1.1,x1.2,x1.i,x2)({\color{red}x_1}, x_{1.1}, x_{1.2}, x_{1.i}, {\color{red}x_2}). Again, all combinations are evaluated as described above, resulting in a larger number of RiR_i in which minRi\min R_i and maxRi\max R_i may be closer to minf(x,y,z,...)\min f(x, y, z, ...) and maxf(x,y,z,...)\max f(x, y, z, ...), respectively. This is the "combinatorial sequence grid evaluation" approach which works quite well in scenarios where monotonicity changes direction (see 'Examples'), obviating the need to create multivariate derivatives (Hessians) or use some multivariate minimization algorithm.
If intervals are of type [x1<0,x2>0][{\color{red}x_1} < 0, {\color{red}x_2} > 0], a zero is included into the middle of the sequence to avoid wrong results in case of even powers, i.e. [1,1]2=[1,1][1,1]=[1,1][-1, 1]^2 = [-1, 1][-1, 1] = [-1, 1] when actually the right interval is [0,1][0, 1], see curve(x^2, -1, 1).

Value

A 2-element vector with the resulting interval and an (optional) plot of all evaluations.

Author(s)

Andrej-Nikolai Spiess

References

Wikipedia entry is quite good, especially the section on the 'dependency problem':
http://en.wikipedia.org/wiki/Interval_arithmetic

Comparison to Monte Carlo and error propagation:
Interval Arithmetic in Power Flow Analysis.
Wang Z & Alvarado FL.
Power Industry Computer Application Conference (1991): 156-162.

Computer implementation
Interval arithmetic: From principles to implementation.
Hickey T, Ju Q, Van Emden MH.
JACM (2001), 48: 1038-1068.

Complete Interval Arithmetic and its Implementation on the Computer.
Kulisch UW.
In: Numerical Validation in Current Hardware Architectures. Lecture Notes in Computer Science 5492 (2009): 7-26.

Examples

## Example 1: even squaring of negative interval.
EXPR1 <- expression(x^2)
DAT1 <- data.frame(x = c(-1, 1))
interval(DAT1, EXPR1)

## Example 2: A complicated nonlinear model.
## Reduce sequence length to 2 => original interval
## for quicker evaluation.
EXPR2 <- expression(C * sqrt((520 * H * P)/(M *(t + 460))))
H <- c(64, 65)
M <- c(16, 16.2)
P <- c(361, 365)
t <- c(165, 170)
C <- c(38.4, 38.5)
DAT2 <- makeDat(EXPR2)
interval(DAT2, EXPR2, seq = 2)

## Example 3: Body Mass Index taken from
## http://en.wikipedia.org/w/index.php?title=Interval_arithmetic
EXPR3 <- expression(m/h^2)
m <- c(79.5, 80.5)
h <- c(1.795, 1.805)
DAT3 <- makeDat(EXPR3)
interval(DAT3, EXPR3)

## Example 4: Linear model.
EXPR4 <- expression(a * x + b)
a <- c(1, 2)
b <- c(5, 7)
x <- c(2, 3)
DAT4 <- makeDat(EXPR4)
interval(DAT4, EXPR4)

## Example 5: Overestimation from dependency problem.
# Original interval with seq = 2 => [1, 7]
EXPR5 <- expression(x^2 - x + 1)
x <- c(-2, 1)
DAT5 <- makeDat(EXPR5)
interval(DAT5, EXPR5, seq = 2)

# Refine with large sequence => [0.75, 7]
interval(DAT5, EXPR5, seq = 100)
# Tallies with curve function.
curve(x^2 - x + 1, -2, 1)

## Example 6: Underestimation from dependency problem.
# Original interval with seq = 2 => [0, 0]
EXPR6 <- expression(x - x^2)
x <- c(0, 1)
DAT6 <- makeDat(EXPR6)
interval(DAT6, EXPR6, seq = 2)

# Refine with large sequence => [0, 0.25]
interval(DAT6, EXPR6, seq = 100)
# Tallies with curve function.
curve(x - x^2, 0, 1)

[Package propagate version 1.0-6 Index]