# Generation of random variables via composition and inversion

## What are the main pros and cons of each method and when to use each one?

Law [2007] mentions that:

"Again, the reader is encouraged to develop the inverse-transform method for generating a random variate from the right-trapezoidal distribution in Example 8.4. Note that especially if a is large, the composition method will be faster than the inverse transform, since the latter always requires that a square root be taken, while it is quite likely (with probability a) that the former will simply return $$X= U_2 \sim U(0, 1)$$. This increase in speed must be played off by the analyst against the possible disadvantage of having to generate two or three random numbers to obtain one value of X."

So basically, for the right-trapezoidal distribution (mentioned in the quote) it might be useful to use composition since it's likely to be faster. However, I am not entirely sure why it so?

Say I am generating a random variable using the inverse method. Have a CDF and need to find it's inverse. If it's impossible I am likely to use the inverse transform sampling using the fact that the CDFs are weakly monotonic and right-continuous and will generalize the inverse CDF to a form of $$F^{-1}_X (u) = \inf\{ x|F_X(x)\geq u\}$$ for $$0. It can be a little time-consuming but say I will get the inverse CDF. Then I am creating uniform variable $$U$$ and use it to generate the random variable using the inverse CDF.

I do not understand how the composition method is quicker. Since for composition, I first need to identify that I can in fact use composition from its PDF. Then I create the first uniform variable $$U_1$$ use it to decide on the distribution (case of composition I assume there are at least two). And then I create another uniform variable $$U_2$$ and use it to generate random variable. But I still need to know the inverse so that I know what to output? It just seems to me that in the composition method I still need to calculate the inverse in order to know what the random variable will be.

LAW, A. M. [2007], Simulation Modeling and Analysis, 4th ed., McGraw-Hill, NewYork.

The minimum working example that I can think of is $$f(x) = 0.5 f_A(x) + 0.5 f_B(x)$$, where A is exponential with parameter $$0.5$$ and B is exponential with parameter $$2$$. By looking at the PDF $$f(x)$$ I can see that we can here use the composition method. And the generation method will be something like:

1. Generate $$U_1,U_2 \sim U(0,1)$$.
2. If $$U_1 < 0.5$$, then output $$X = -\frac{\log(U_2)}{0.5}$$
3. Else output $$X = -\frac{\log(U_2)}{2}$$

But if I didn't know what the inverse of exponential distribution is, I would not know what to output. And the composition method would also involve the inversion part. Am I right? I don't see how it would be faster in that scenario.

• Hi: You get the inverse by calculating the cdf and then setting it to the uniform and solving for $X$ so how would not know what it was ? I may not be folllowing and not intending to sound rude. Mar 1 at 19:22
• The composition method lets you model your distribution as a mixture of other densities. So as you say, you start by choosing one of them appropriately, then generating a sample from it. When you get to that stage, however, you don't necessarily need to invert the CDF or numerically solve $F_k(X) = U$ for $X$. You can also use rejection-acceptance or maybe an off-the-shelf generator. Mar 1 at 20:14
• @Aruralreader You could post that comment as an answer Mar 1 at 21:03
• Inverting $0.5(1-e^{-0.5x})+0.5(1-e^{-2x})=u$ requires numerical steps hence is more costly than the composition approach. If this is not clear enough, imagine a mixture of 10³ exponential densities... Mar 1 at 21:03
• @Xi'an, I kind of get your point but not entirely. The way I understand it is: the inverse method in this case requires one $U$, one compare and one logarithm. And composition, requires two $U$, one compare and one logarithm. Therefore, I would argue that it is more efficient to use inverse method. Mar 1 at 21:16

In general, you should not be concerned about efficiency (performance) between random variate generation methods, unless you have written an implementation of them, compared their running time, used them in your application, and found the running time to be unacceptable in your application. This is a general issue in programming that is obviously not limited to random variate generation, an issue commonly known as "premature optimization".

The bigger concern is ease of implementation or ease of sampling, or rather the availability of algorithms that sample the distribution in question. (Accuracy of sampling the distribution is also often important.) Although numerous algorithms are now available for sampling "standard" distributions, such as normal, gamma, and beta, many other distributions that occur in practice have non-trivial probability functions (such as PDF, CDF, and/or inverse CDF), and many of them have PDFs whose normalizing constant is intractable, for example. A notable example is certain Bayesian posteriors whose PDFs cannot be evaluated pointwise.

The following question shows one example involving the beta distribution:

Here is a crude R code comparing composition versus inversion for a mixture of $$K$$ exponential distributions:

library(rbenchmark)
#target creation
K=100
we=sort(runif(K),d=T) #mixture weight
we=we/sum(we)
wes=cumsum(we) #cumulated sum
la=rexp(K)     #expomential rates
lah=median(la)
F=function(x)sum(we*exp(-la*x)) #tail cdf

benchmark("compo"={
#composition
x=rexp(1,la[1+sum(runif(1)>wes)])
},
"inve"={
#Newton inversion
u=runif(e<-1)
x0=-log(u)/lah
#precision of 1e-3
while(abs((f<-u-F(x0)))>1e-3)x0=x0-f/lah
},
replication=1e3,
columns = c("test","elapsed"))


returning the times as

         test elapsed
1       compo   0.001
2        inve   0.231


Even using a much more efficient inversion algorithm preverses the comparison: with

Q = function(u) uniroot((function(x) F(x) - u), lower = 0,
upper = qexp(.99,rate=min(la)))[1] #numerical tail quantile
x=Q(runif(1))


the benchmarked times are

         test elapsed
1       compo   0.001
2        ïnve   0.235
3     uniroot   0.017


Moving to K=99999 components in the mixture leads to

         test elapsed
1       compo   0.057
2        inve  45.736
3     uniroot   5.814


which makes the case in favour of composition versus inversion...

• What is the uniroot method? Mar 31 at 6:13
• ?uniroot provides the R documentation on this function: The function ‘uniroot’ searches the interval from ‘lower’ to ‘upper’ for a root (i.e., zero) of the function ‘f’ with respect to its first argument. Mar 31 at 6:15
• ah so 'uniroot' in your table denotes the inversion method using the uniroot root finding algorithm. Which performs much better than inversion method with Newton's method Mar 31 at 6:17
• I wonder, can we find a case where inversion method outperforms composition by a big margin? I can't see how we could: the composition step will always be very cheap, right? unless somehow the whole expression was easier to invert than the individual parts. hmm Mar 31 at 6:20
• While checking Reuven Rubinstein's book, I noticed that he classified stratified sampling as a composition method. In which case finding the probability weights of each stratum may prove costly. Another potentially costly setting may be the one of infinite mixtures if selecting the component proves dear. Mar 31 at 7:53