Probability Density Functions for Probabilistic Uncertainty Analysis

Define a distribution for PSA parameters.

Usage

normal(mean, sd)

lognormal(mean, sd, meanlog, sdlog)

gamma(mean, sd)

binomial(prob, size)

multinomial(...)

logitnormal(mu, sigma)

beta(shape1, shape2)

triangle(lower, upper, peak = (lower + upper)/2)

poisson(mean)

define_distribution(x)

beta(shape1, shape2)

triangle(lower, upper, peak = (lower + upper)/2)

use_distribution(distribution, smooth = TRUE)

Arguments

mean: Distribution mean.
sd: Distribution standard deviation.
meanlog: Mean on the log scale.
sdlog: SD on the log scale.
prob: Proportion.
size: Size of sample used to estimate proportion.
...: Dirichlet distribution parameters.
mu: Mean on the logit scale.
sigma: SD on the logit scale.
shape1: for beta distribution
shape2: for beta distribution
lower: lower bound of triangular distribution.
upper: upper bound of triangular distribution.
peak: peak of triangular distribution.
x: A distribution function, see details.
distribution: A numeric vector of observations defining a distribution, usually the output from an MCMC fit.
smooth: Use gaussian kernel smoothing?

Details

These functions are not exported, but only used in define_psa(). To specify a user-made function use define_distribution().

use_distribution() uses gaussian kernel smoothing with a bandwidth parameter calculated by stats::density(). Values for unobserved quantiles are calculated by linear interpolation.

define_distribution() takes as argument a function with a single argument, x, corresponding to a vector of quantiles. It returns the distribution values for the given quantiles. See examples.

Examples

define_distribution(
  function(x) stats::qexp(p = x, rate = 0.5)
)
#> [[1]]
#> function (x) 
#> stats::qexp(p = x, rate = 0.5)
#> <environment: 0x559b7134bfb0>
#> 

# a mixture of 2 gaussians
x <- c(rnorm(100), rnorm(100, 6))
plot(density(x))


use_distribution(x)
#> [[1]]
#> function (x) 
#> {
#>     if (smooth) {
#>         noise <- stats::rnorm(n = length(x), mean = 0, sd = stats::density(distribution)$bw)
#>     }
#>     else {
#>         noise <- 0
#>     }
#>     (stats::approxfun(x = seq(0, 1, length = length(distribution)), 
#>         y = distribution))(x) + noise
#> }
#> <bytecode: 0x559b711a49c8>
#> <environment: 0x559b711a2bf8>
#>