rstatisticssimulationprobabilityexponential-distribution

Defining exponential distribution in R to estimate probabilities

I have a bunch of random variables (X1,....,Xn) which are i.i.d. Exp(1/2) and represent the duration of time of a certain event. So this distribution has obviously an expected value of 2, but I am having problems defining it in R. I did some research and found something about a so-called Monte-Carlo Stimulation, but I don't seem to find what I am looking for in it.

An example of what i want to estimate is: let's say we have 10 random variables (X1,..,X10) distributed as above, and we want to determine for example the probability P([X1+...+X10<=25]).

Thanks.

Solution

  • Are you aware of rexp() function in R? Have a look at documentation page by typing ?rexp in R console.

    A quick answer to your Monte Carlo estimation of desired probability:

    mean(rowSums(matrix(rexp(1000 * 10, rate = 0.5), 1000, 10)) <= 25)

    I have generated 1000 set of 10 exponential samples, putting them into a 1000 * 10 matrix. We take row sum and get a vector of 1000 entries. The proportion of values between 0 and 25 is an empirical estimate of the desired probability.

    Thanks, this was helpful! Can I use replicate with this code, to make it look like this: F <- function(n, B=1000) mean(replicate(B,(rexp(10, rate = 0.5)))) but I am unable to output the right result.

    replicate here generates a matrix, too, but it is an 10 * 1000 matrix (as opposed to a 1000* 10 one in my answer), so you now need to take colSums. Also, where did you put n?

    The correct function would be

    F <- function(n, B=1000) mean(colSums(replicate(B, rexp(10, rate = 0.5))) <= n)

    For non-Monte Carlo method to your given example, see the other answer. Exponential distribution is a special case of gamma distribution and the latter has additivity property.

    I am giving you Monte Carlo method because you name it in your question, and it is applicable beyond your example.