Find Random Number Generator following the density $f (x) = \frac{1 + \alpha x}{2}$, $ −1 ≤ x ≤ 1$, $−1 ≤\alpha ≤ 1$

211 Views Asked by At

How could random variables with the following density function be generated from a uniform random number generator? $f (x) = \frac{1 + \alpha x}{2}$, $ −1 ≤ x ≤ 1$, $−1 ≤\alpha ≤ 1$.

To find the Random number generator one first needs to find the distribution function and then the quantile function. However, I have troubles solving this on my own:

\begin{align} F(x) & = \int_{-1}^x \frac{1+\alpha u}{2} du\\ & = \frac 12 \Big( (x+\frac{\alpha x^2}{2})-(-1+\frac{\alpha}{2})\Big)\\ & = \alpha x^2 + 2x + 2 - \alpha \end{align}

I am not sure how to handle the $\alpha$. For getting the quantile(=inverse) function i need to solve this last expression for x.

1

There are 1 best solutions below

0
On

Not sure how you get from your second line to the third, there is clearly error somewhere

\begin{equation} F(x) = \frac{\alpha}{4}(x^2 - 1) + \frac{1}{2}(x+1) \end{equation}

with obvious properties $F(x=-1)=0$ and $F(x=1)=1$. Now assign $F(x)$ to RV $\xi$ in the [0...1] range and find inverse:

\begin{equation} x = F^{-1}(\xi) \end{equation}

Simple quadratic equation with solution:

\begin{equation} x = \frac{2 \; \sqrt{\frac{1}{4} + \alpha(\frac{\alpha}{4} + \xi - \frac{1}{2})} - 1}{\alpha} \end{equation}

For the case when $\alpha=0$, $x=2\xi-1$

Easy to check:

\begin{equation} \xi=0, \;then\; x=-1 \end{equation} \begin{equation} \xi=1, \;then\; x=+1 \end{equation}

Python code

import numpy as np
import matplotlib.pyplot as plt

def PDF(x, alpha):
    return (1.0 + alpha*x) / 2

def sampleAlpha(n, alpha):
    U01 = np.random.uniform(size=n)
    if alpha == 0.0:
        return 2.0*U01 - 1.0

    D  = 0.25 + alpha*(U01 + 0.25*alpha - 0.5)
    return (2.0*np.sqrt(D) - 1.0) / alpha

nbins = 100
xmin = -1.0
xmax =  1.0
alpha = 0.3

t = np.linspace(xmin, xmax, nbins + 1)

q = sampleAlpha(100000, alpha)
h, bedges = np.histogram(q, bins = t, density = True)

p = PDF(t, alpha)

plt.bar(bedges[:-1], h, width = (xmax - xmin)/float(nbins),  align='edge')
plt.plot(t, p, 'r') # plotting PDF
plt.xlim(xmin, xmax)
plt.show()

generates graphs like one below for sampled vs PDF

enter image description here