Question

Suppose we want to simulate a large number \(n\) of independent exponentials with rate \(1-\) call them \(X_{1}, X_{2}, \ldots, X_{n} .\) If we were to employ the inverse transform technique we would require one logarithmic computation for each exponential generated. One way to avoid this is to first simulate \(S_{n}\), a gamma random variable with parameters \((n, 1)\) (say, by the method of Section 11.3.3). Now interpret \(S_{n}\) as the time of the \(n\) th event of a Poisson process with rate 1 and use the result that given \(S_{n}\) the set of the first \(n-1\) event times is distributed as the set of \(n-1\) independent uniform \(\left(0, S_{n}\right)\) random variables. Based on this, explain why the following algorithm simulates \(n\) independent exponentials: Step 1: Generate \(S_{n}\), a gamma random variable with parameters \((n, 1)\). Step 2: Generate \(n-1\) random numbers \(U_{1}, U_{2}, \ldots, U_{n-1}\). Step 3: Order the \(U_{i}, i=1, \ldots, n-1\) to obtain \(U_{(1)}

Short answer

The provided algorithm simulates \(n\) independent exponentials with rates equal to 1 by following these steps: 1) Generate a Gamma random variable, \(S_n\), which consists of the sum of \(n\) independent exponentials; 2) Generate \(n-1\) uniform random numbers, and order them in ascending order to represent the first \(n-1\) event times; 3) Calculate the difference between successive uniform random variables to simulate the exponentials with the given rate. The combination of Gamma distribution and ordered uniform random numbers confirms the algorithm's effectiveness in simulating the desired exponentials.

Step-by-step solution

Generate a Gamma random variable with parameters (n, 1)

In this step, we generate a random variable \(S_n\) that follows a gamma distribution with parameters \((n, 1)\). The parameter \((n, 1)\) corresponds to \(n\) independent exponentially distributed random variables \(X_1, X_2, ... , X_n\), each with a rate of 1. You can think of \(S_n\) as the sum of the \(n\) exponential random variables: \(S_n = X_1 + X_2 + ... + X_n\).

Generate n-1 uniform random numbers

In this step, we generate \(n-1\) independent random numbers \(U_1, U_2, ..., U_{n-1}\), each uniformly distributed between 0 and 1. These random numbers will be used to calculate the difference between successive uniform random variables, which will later help us in obtaining exponential random variables.

Order the uniform random numbers

After generating the \(n-1\) uniform random numbers, we order them in ascending order. Let the ordered uniform random numbers be denoted by \(U_{(1)} < U_{(2)} < ... < U_{(n-1)}\).

Calculate the n exponential random variables

In this step, we first set \(U_{(0)}=0\) and \(U_{(n)}=1\). Next, we calculate the n exponential random variables \(X_i\) using the equation: \(X_i = S_n\left(U_{(i)} - U_{(i-1)}\right), i=1, \ldots, n\). The above statement shows that given \(S_n\), the time of the \(n\)-th event of a Poisson process with rate 1, the set of the first \(n-1\) event times is distributed as the set of \(n-1\) independent uniform \((0, S_n)\) random variables. Thus, when we calculate the \(X_i = S_n\left(U_{(i)} - U_{(i-1)}\right)\) values, what we are doing is essentially simulating \(n\) independent exponentials with rates equal to 1. This proves the effectiveness of the provided algorithm in simulating the given conditions.