Lab Exercise II: Generating Continuous Random Variables

The Inverse Transform Method

Let \(U\) be a uniform \((0, 1)\) random variable. For any continuous distribution function \(F\) the random variable \(X\) defined by

\(X = F^{−1}(U)\) has distribution \(F\). [\(F^{−1}(u)\) is defined to be that value \(x\) such that \(F(x) = u\).]

Algorithm: Main Steps

  • STEP 1: Generate \(U\) having a uniform density.

  • STEP 2: Set \(X=F^{−1}(U)\)

Exercise 1:

Consider a random variable \(X\) with distribution function: \(F(x)= x^n\), where \(0 <x< 1\)

  1. Give the inverse transform algorithm for generating a random variable \(X\).

  2. Use the inverse transform algorthm to write an R code for generating \(10\) random values of \(X\)

  3. Re-execute the R code in (part 2) with \(1000\) values of \(X\) and illustrate the distribution of the random values graphically.

  4. Consider \(1000\) random values of \(X\), estimate \(E(X)\)

Answers:

  1. List the steps/algorithm ….

    - compute X using inverse transform: \(X=F^{−1}(U)= U^{1/n}\)

    - generate U values from uniform distribution(0,1)

    - apply the inverse function for U,n to obtain X

    -plot the distribution of the generated X values

    - estimate E(X) by finding the mean of the generated x

  2. Write your R code in the code chunk below:

set.seed(111)
### Write your R code here
### part 1
inverseu= function(u,n){
  return(u^(1/n))
}
### part 2
U1= runif(10)
n= 4
x= inverseu(U1,n)
hist(x,main = "10 Random values of X", col="#557D99")

### part 3
U2= runif(1000)
x2= inverseu(U2,n)
hist(x2,main = "10 Random values of X", col="#557D99" )

### part 4
estimated = mean(x2)
paste("E(X)= ",estimated)
## [1] "E(X)=  0.801461879482378"

The Rejection Method

Suppose we have a method for generating a random variable having density function \(g(x)\). We can use this as the basis for generating from the continuous = distribution having density function of \(f (x)\) by generating \(Y\) from \(g\) and then accepting this generated value with aprobability proportional to \(\frac{f(Y)}{g(Y)}\). Specifically, let \(c\) be a constant such that \(\frac{f (y)}{g(y)} \le c\) for all \(y\).

Algorithm: Main Steps

  • STEP 1: Generate \(Y\) having density \(g\).

  • STEP 2: Generate a random number \(U\).

  • STEP 3: If \(U\le \frac{f(Y)}{cg(Y)}\), set \(X = Y\). Otherwise, return to Step 1

Exercise 2:

Consider the density function \(f (x) = 20x(1 − x)^3\), where \(0 <x< 1\)

  1. What is your \(g(x)\).[ Hint: see example 5f, page 74]

  2. Give the acceptance-rejection algorithm for generating a random variable $X$.

  3. Use the acceptance-rejection algorithm to write an R code for generating \(10\) random values of \(X\)

  4. Re-execute the R code in (part 3) with \(1000\) values of \(X\) and illustrate the distribution of the random values graphically.

Answers:

  1. Proposal function:

we choose g(x) to be uniform dist over(0,1) => g(x)= 1

  1. List the steps/algorithm

    - define f function

    - define g function

    - define Accept- reject function:

    - generate y from uniform (0,1)

    -and the same for u

    - compute f(y)/g(y) = 20y(1 − y)^3 / 1 = 20y(1 − y)^3

    - if u<= that value , then accept the value of y and store it in X, otherwise go to the next value

    - repeat the process and use the function on n= 10 and n =1000

    - find the histogram of these generated values of X

  2. Write your R code in the code chunk below:

set.seed(111)
### part 1 
f= function(x){
  return(20 * x * (1 - x)^3)
}
g= function(x){
  return(1)
}
### part 2
accrej= function(n){
  x= numeric(n)
  i=1
  while (i<=n) {
    y= runif(1)
    u= runif(1)
    if (u<=20 * y * (1 - y)^3 ) {
      x[i] = y
      i = i+1
    }
  }
  return(x)
}
### part 3
xv= accrej(10)
hist(xv,main = "10 Random values of X", col="#557D99")

### part 4
xw= accrej(1000)
hist(xw,main = "10 Random values of X", col="#557D99")

N.B Submit your answers as a html file.