# Generating a Random Value Vector from an Exponential Distribution using R

Given a standard PDF of the form $$f(x)=ae^{-ax}$$ with domain $$[0,+\infty)$$, its CDF being $$F(x)=1-e^{-ax}$$, and a mutated CDF that takes $$p \in [0,1]$$ as a probability and returns the corresponding $$x$$ value at that percentile: $$x=-\frac{ln(1-p)}{a}$$ (quantile equation)

I am trying to write an R program that returns a randomized vector of size n of from this PDF distribution. The function below represents the quantile equation from above.

makeExpDist <- function(a, p) {
return(log(1-p)/(-a))
}


I then use this last function to map a randomized vector (size n) from a uniformly distributed variable $$[0,1]$$ into the above function, expecting that the returned vector will represent a sample from the original PDF.

giveSample <- function(cdf, a, n){
random <- runif(n)
le <- c()
for(i in random){
le <- append(le, cdf(a,i))
}
return(data.frame(id = 1:n, le_prob = random, le_num = le))
}


I would expect the mean to be $$86.6434$$ (rounded) since $$-\frac{ln(.5)}{.008}=86.6434$$ (rounded), which can be supported by $$\int_{0}^{86.64}.008e^{-.008x}dx=.5$$

To demonstrate:

makeExpDist(.008,.5)


returns

 86.6434


I would also expect that to be the mean of a sufficiently large sample from the PDF. I try this below, but the mean of the resulting vector is always $$~125$$:

testframe <- giveSample(makeExpDist, .008, 100000) #test
mean(testframe$le_num) mean(testframe$le_prob)


... which on my last run printed:

 125.7898
 0.5016412


The value of the second mean as expected is around 0$$.5$$ since that simply is the uniform variable I used.

So my question is why do I keep receiving $$~125$$ as the mean of randomized values when the true mean of the PDF is $$~86$$? I have tried many fixes in R and get the same exact mean of $$~125$$. Also, when I run the same process in excel of mapping a randomized uniform vector onto this distribution, I get the same output of $$~125$$. Is it wrong to use a random uniform variable $$[0,1]$$ to randomize a value in the PDF? Given the expected value of runif() is 0.5 why is the expected value of the return vector not equal to makeExpDist(.008,.5)?

Thanks a bunch!

• Where did the "$ln(.5)$" come from?? The mean of your distribution is $1/a.$
– whuber
Jul 8, 2021 at 19:37
• @whuber I think OP just confused 50th percentile (plugging in $p=0.5$ into their quantile function) for the mean. Jul 8, 2021 at 19:49
• @bdeonovic Good interpretation! (+1)
– whuber
Jul 8, 2021 at 20:34
• You're right, I did confuse the two terms. The mean is $1/.008=125$ which is what I was looking for. Used the wrong integration technique for mean. Jul 9, 2021 at 20:11

The quantity you calculated (86.6434) is the 50th percentile, ie the median, not the mean. Since the exponential is a skewed distribution the mean does not equal the median.

Try

> median(giveSample(makeExpDist, 0.008, 1000)\$le_num)
 87.15622


The median for an $$\text{Exponential}(a)$$ is $$log(2)/a$$ while the mean is $$1/a$$.

> log(2)/0.008
 86.6434
> 1/0.008
 125


This is not the place for debugging R code, so I will focus on the the theory. It seems you are trying to use the quantile function (inverse CDF) of a sample from $$\mathsf{Unif}(0,1)$$ to get a random sample from $$\mathsf{Exp}(\mathrm{rate}=\lambda=0.008)$$

As you say, the quantile function is $$X = -\log(1-U)/\lambda.$$ Notice that $$U^\prime = 1 - U \sim\mathsf{Exp}(\lambda)$$ also, so it is a little more efficient to use $$X = -\log(U)/\lambda.$$ (But either method works.) Here is simple R code, which gives the mean of a sample of a million observations from $$\mathsf{Exp}(\lambda=0.008)$$ as $$124.93 \pm 0.25,$$ as anticipated.

set.seed(2021)
u = runif(10^6)
x = -log(u)/.008
mean(x)
 124.9272     ## aprx E(X) = 125
2*sd(x)/1000
 0.250065     ## aprx 95% margin of simulation error

mean(-log(1-u)/0.008)
 125.1309     ## alternate method


In the histograms below images and objects are the same color. R code for figure:

frb = rainbow(8)
cutp.u=seq(0,1, by=.2)
par(mfrow=c(1,2))
hist(u, prob=T, br=cutp.u, col=frb[c(1,2,3,5,6)])