# How to draw from conditional Weibull distribution

I would like to draw random numbers in R from a conditional Weibull distribution which is given by P(x<=t+T| x>T) = 1-(1-F(t+T))/(1-F(T)), where F is the Weibull-CDF and T is the time that has elapsed without event. I know the function rweibull for drawing random numbers, but is there a similar function for the conditional distribution?

• There is the truncdist package which implements the inversion method mentioned by @jcken Jun 3, 2020 at 8:42
• @Jarle Tufto, thanks. Will look into this. Jun 3, 2020 at 12:30

Assume $$T, \lambda, \kappa$$ are fixed. We can use the Probability Integral Transform (PIT). In essense, the PIT says that the CDF of a random variable follows a $$U(0,1)$$ distribution, so drawing $$U(0,1)$$ samples and running them through $$F^{-1}(x)$$ gives draws of $$X$$. If we can draw samples of $$x$$, we can then estimate the $$P(x < t + T | X > T)$$ via monte carlo. Nb I use the form of the Weibull that is used on wikipedia [ https://en.wikipedia.org/wiki/Weibull_distribution ].

To apply the inversion method, we want to sample values of $$x = t + T$$.

We have \begin{align} p = P(x < t+T | X>T) &= 1 - \frac{1 - F(t+T)}{1 - F(T)} \\ & = 1 - \frac{1 - \left( 1 - \exp\left\{ -\lambda^{-\kappa}(t+T)^{\kappa} \right\} \right)}{1 - \left( 1 - \exp\left\{ -\lambda^{-\kappa}(T)^{\kappa} \right\} \right)}\\ & = 1 - \exp \left\{ \lambda^{-\kappa}T^\kappa - \lambda^{-\kappa}(t + T)^\kappa \right\} \end{align}

Some algebra leads to $$t + T = \left\{ T^\kappa - \lambda^\kappa \log(1-p) \right\}^{1/\kappa}$$

The following R code can then be used to draw samples of $$X$$.

draw.x <- function(n, lambda, kappa, T0){
p <- runif(n)

( T0^kappa - lambda^kappa *log(1-p))^(1/kappa)

}

## estimate P(X < 5+1 = 6 | X > 5)
n.draws <- 10^4

samples <- draw.x(n.draws, 1.5, 1, 5)

prob.estimate <- sum(samples < 6)/n.draws
prob.estimate
## if kappa = 1=> exponential distribution rate = 1/lambda
## compare to exponential 1/lambda
pexp(6-5, 1/1.5); prob.estimate

hist(samples); abline(v = 6,col = 2, lwd = 2)

$$$$

• Amazing, thank you! I was aware of the Inverse Method (what you call PIT), but didn't realize that we can manually derive the inverse in this case. Helped me a lot. Jun 3, 2020 at 12:32

In Mathematica you can do it easily,

 f[a_, b_, T_] =
ProbabilityDistribution[
1 - (1 - CDF[WeibullDistribution[a, b], t + T])/
CDF[WeibullDistribution[a, b], T], {t, T, \[Infinity]}];
data = RandomVariate[f[3.5, 1, 2], 25]
`

{2.32932, 2.63388, 2.29273, 3.60003, 3.64601, 3.27736, 2.89225, \ 2.9748, 3.18756, 2.49238, 2.24663, 3.58997, 3.282, 3.28332, 2.54555, \ 3.33817, 2.6722, 2.94851, 3.07248, 2.23742, 2.49302, 3.13033, \ 3.23821, 3.47035, 3.50443}

I am not sure about t and T, you can adjust them yourself.

• Good to know, thanks for the answer! I was just looking for an R solution in this case. Jun 3, 2020 at 12:33