Simulate a random variable from given PDF

Question

Let $f(x) = \begin{cases} \frac{1}{a}\log(\frac{a}{x}), & \text{if } 0 < x < a \\ 0,& x \notin(0,a) \end{cases}$, $a \text{ fixed}$

I want to generate a random variable $X \sim f(x)$, where $f(x)$ is the PDF. I tried with inverse method but no result (transcedental equation). I want to try rejection method now, but suppose that I let $Y \sim \operatorname{Uniform}(0,a) \implies g(x) = \frac{1}{a}, x\in(0,a)$ is the density of $Y$. Now when I want to find the constant $c \ge \frac{f(x)}{g(x)} \forall x\in(0,1)$ I can't because $\frac{f(x)}{g(x)}$ is unbounded close to $0$. So ... how I can simulate a random variable $X \sim f(x)$? Thanks!

Answer 1

Let $Y \sim \operatorname{Gamma}(2,1)$ with PDF $$f_Y(y) = y e^{-y}, \quad y > 0.$$ Then $X = g(Y) = a e^{-Y}$ for a parameter $a > 0$ has PDF $$f_X(x) = f_Y(g^{-1}(x)) \left|\frac{dg^{-1}}{dx}\right| = f_Y\left(-\log \frac{x}{a}\right) \frac{1}{x} = -\frac{x}{ax} \log \frac{x}{a} = \frac{1}{a} \log \frac{a}{x}, \quad 0 < x < a.$$ Consequently, we can generate realizations from $X$ by first generating a gamma random variable $Y$, then compute $X = ae^{-Y}$.

In turn, we generate realizations from $Y$ by simply computing pairs of IID exponential variates with mean $1$, and taking their sum. And the exponential variates can be generated using the usual inverse transform method of a uniform random variate.

---

The theorem used for the transformation (restated here for the appropriate random variables) is:

Let $Y$ have PDF $f_Y(y)$ and let $X = g(Y)$, where $g$ is a monotone function. Let $\mathcal X$ and $\mathcal Y$ be defined by $\mathcal Y = \{y : f_Y(y) > 0\}$, $\mathcal X = \{x : x = g(y) \text{ for some } y \in \mathcal Y\}$. Suppose that $f_Y(y)$ is continuous on $Y$ and that $g^{-1}(x)$ has a continuous derivative on $\mathcal X$. Then the PDF of $X$ is given by $$f_X(x) = \begin{cases} f_Y(g^{-1}(x)) \left|\frac{d}{dx}\left[g^{-1}(x)\right]\right|, & x \in \mathcal X \\ 0, & \text{otherwise}. \end{cases}$$

Answer 2

@heropup has a nice, direct solution and gets my upvote. An alternative is to use acceptance-rejection. Without loss of generality assume $a = 1$ for the scale factor. The density $f(x) := -\log x$ on $(0, 1)$ is dominated by the density $g(x) := x^{-1/2}/2$ in the sense that $Cg(x) \geq f(x)$ for large enough $C$. A little computation shows the best choice for the constant, the smallest, is $C = 4/e$.

The distribution function associated with $g$ is $G(x) = x^{1/2}$, so $Y := U^2$ is distributed as $g$ when $U$ is uniform on $(0, 1)$. The acceptance-rejection algorithm becomes:

1. Generate $Y \sim g$.
2. Generate $V$ uniform on $(0, 1)$.
3. If $Cg(Y)V \leq f(Y)$ accept $Y$, otherwise repeat at Step 1.

---

Numerical tests show this algorithm has an efficiency of 86 percent.