Compute numerically expected values and variances of renewal count processes.

evCount_conv_bi(
  xmax,
  distPars,
  dist = c("weibull", "gamma", "gengamma", "burr"),
  method = c("dePril", "direct", "naive"),
  nsteps = 100,
  time = 1,
  extrap = TRUE
)

evCount_conv_user(
  xmax,
  distPars,
  extrapolPars,
  survR,
  method = c("dePril", "direct", "naive"),
  nsteps = 100,
  time = 1,
  extrap = TRUE
)

Arguments

xmax

unsigned integer maximum count to be used.

distPars

TODO

dist

TODO

method

TODO

nsteps

unsiged integer, number of steps used to compute the integral.

time

double, time at wich to compute the probabilities. Set to 1 by default.

extrap

logical, if TRUE, Richardson extrapolation will be applied to improve accuracy.

extrapolPars

ma::vec of length 2. The extrapolation values.

survR

function, user supplied survival function; should have signature function(t, distPars), where t is a positive real number (the time where the survival function is evaluated) and distPars is a list of distribution parameters. It should return a double value.

Value

a named list with components "ExpectedValue" and "Variance".

Details

evCount_conv_bi computes the expected value and variance of renewal count processes for the builtin distirbutions of inter-arrival times.

evCount_conv_user computes the expected value and variance for a user specified distirbution of the inter-arrival times.

Examples

pwei_user <- function(tt, distP) {
    alpha <- exp(-log(distP[["scale"]]) / distP[["shape"]])
    pweibull(q = tt, scale = alpha, shape = distP[["shape"]],
             lower.tail = FALSE)
}

## ev convolution Poisson count
lambda <- 2.56
beta <- 1
distPars <- list(scale = lambda, shape = beta)

evbi <- evCount_conv_bi(20, distPars, dist = "weibull")
evu <- evCount_conv_user(20, distPars, c(2, 2), pwei_user, "dePril")

c(evbi[["ExpectedValue"]], lambda)
#> [1] 2.560001 2.560000
c(evu[["ExpectedValue"]], lambda )
#> [1] 2.560001 2.560000
c(evbi[["Variance"]], lambda     )
#> [1] 2.560003 2.560000
c(evu[["Variance"]], lambda      )
#> [1] 2.56 2.56

## ev convolution weibull count
lambda <- 2.56
beta <- 1.35
distPars <- list(scale = lambda, shape = beta)

evbi <- evCount_conv_bi(20, distPars, dist = "weibull")
evu <- evCount_conv_user(20, distPars, c(2.35, 2), pwei_user, "dePril")

x <- 1:20
px <- dCount_conv_bi(x, distPars, "weibull", "dePril",
                     nsteps = 100)
ev <- sum(x * px)
var <- sum(x^2 * px) - ev^2

c(evbi[["ExpectedValue"]], ev)
#> [1] 1.968046 1.968046
c(evu[["ExpectedValue"]], ev )
#> [1] 1.968046 1.968046
c(evbi[["Variance"]], var    )
#> [1] 1.348213 1.348213
c(evu[["Variance"]], var     )
#> [1] 1.348214 1.348213