[Math] Pseudorandom Number Generator Using Uniform Random Variable

probabilityprobability distributionsrandomstatisticsuniform distribution

I am working out of Mathematical Statistics and Data Analysis by John Rice and ran into the following interesting problem I'm having trouble figuring out.

Ch 2 (#65)

How could random variables with the following density function be generated from a uniform random number generator?

$$f(x) = \frac{1 + \alpha x}{2}, \quad -1 \leq x \leq 1,\quad -1 \leq \alpha \leq 1$$

So I believe I'm suppose to use the following fact to solve the problem

Proposition D

Let U be uniform on [0, 1], and let X = $F^{-1}$(U). Then the cdf of X is F.

Proof

$$P(X \leq x) = P(F^{-1}(U) \leq x) = P(U \leq F(x)) = F(x)$$

That is, we can use uniform random variables to generate other random variables that will have cdf F

So my goal should then be to find a cdf and it's inverse then give as input to the inverse the uniform random variable. I've included my attempt.

Given $f(x) = \frac{1 + \alpha x}{2}$

$$F(X) = \int_{-1}^{x} \frac{1 + \alpha t}{2} dt \; = \; \frac{x}{2} + \frac{\alpha x}{4} + \frac{1}{2} – \frac{\alpha}{4}$$

$$4 \cdot F(X) – 2 + \alpha = 2x + \alpha x$$

$$F^{-1}(X) = \frac{4X – 2 + \alpha}{2 + \alpha}$$

So our random variable is, for example, T where

$$T = F^{-1}(U) = \frac{4U – 2 + \alpha}{2 + \alpha}$$

The answer in the back of the book is

$$X = [-1 + 2 \sqrt{1/4 – \alpha(1/2 – \alpha / 4 – U)}]/ \alpha$$

I'm not really sure where I went wrong. Any help?

Best Answer

The cdf appears to be wrong. When $-1\leq x\leq 1$, \begin{align*} F_X(x) &= \int_{-1}^{x} \frac{1 + \alpha t}{2} dt\\ &=\int_{-1}^x \frac{1}{2}+\frac{\alpha}{2}t\,dt\\ &=\frac{1}{2}[x+1]+\frac{\alpha}{4}[x^2-1]\\ \end{align*} Other than that, your approach seems fine.

Related Solutions

[Math] Random number generator from a piecewise PDF

You have

$$F(x) = \left\{ \begin{array}{lr} C\log_e(x/a) ,& \qquad x \in (a,b)\\ C\log_e(b/a)+\gamma C\log_e(x/b) ,& \qquad x \in (b,c) \end{array} \right.$$

$$F^{-1}(y) = \left\{ \begin{array}{lr} a\exp\left(\frac{y}{C}\right) ,& \qquad y \in (0,C\log_e(b/a))\\ b\exp\left(\frac{y-C\log_e(b/a)}{\gamma C}\right) ,& \qquad y \in (C\log_e(b/a),1) \end{array} \right.$$

and if $u$ is uniformly distributed, you want $F^{-1}(u)$, with the choice depending on whether it is more or less than $C\log_e(b/a)$.

Generating Standard Uniform random variable

Many statistical software programs use the 'Mersenne Twister' pseudorandom number generator. See this link for some information on this generator (the default in R) and other well vetted generators. More particularly, see Wikipedia on the Mersenne-Twister.

(There are no set rules for making a pseudorandom number generator, using congruential generators or otherwise. We know lots of things not to do, but no sure rules for success. Thus generators are tested using 'batteries' of problems that have proved difficult to simulate. A well-vetted generator is one that has passed many such tests.)

In R statistical software, the function runif samples the indicated number of observations from a uniform distribution. Thus the following R code samples $m = 10,000$ observations from $\mathsf{Exp}(\lambda = 3)$ according to the inverse CDF transformation shown by @peter 5 (+1).

set/seed(416)  # for reproducibility
m = 10^4;  u = runif(m)
x = -(1/3)*log(u)

Histograms of the vectors u and x, each containing $m$ elements, are shown below, along with the respective density functions.

par(mfrow=c(1,2))
hist(u, prob=T, col="skyblue2", main="Uniform Sample with Density of UNIF(0,1)")
 curve(dunif(x), add=T, n=10001, lwd=2, col="red")
hist(x, prob=T, col="skyblue2", main="Exponential Sample with Density of EXP(rate=3)")
 curve(dexp(x, 3), add=T, n=10001, lwd=2, col="red")
par(mfrow=c(1,1))

Notice that the plots show good agreement of the histograms of samples with the density functions of the respective distributions.

More formally, here are results of Kolmogorov-Smirnov goodness-of-fit tests based on the first 5000 observations in each sample (the samplie-size limit allowed by the implementation of this test in R). P-values far above 5% indicate that samples are consistent with the claimed distributions.

ks.test(u[1:5000], "punif")

        One-sample Kolmogorov-Smirnov test

data:  u[1:5000]
D = 0.013297, p-value = 0.3396
alternative hypothesis: two-sided


ks.test(x[1:5000], "pexp", 3)

        One-sample Kolmogorov-Smirnov test

data:  x[1:5000]
D = 0.013297, p-value = 0.3396
alternative hypothesis: two-sided

Best Answer

Related Solutions

[Math] Random number generator from a piecewise PDF

Generating Standard Uniform random variable

Related Question