Sampling

Uniform Disc

One way to sample points in a disc is to generate a random number between $0$ and $r$ for the radius and another one for the angle in $[0,2\pi]$ . The points won’t be uniformly distributed like that, since the area of segments on the inside is different from that of the outside. The sampling will look something like this:

function uniform_radial_circle(){

    const u_1 = Math.random();
    const u_2 = Math.random();
    
    const r = u_1;
    const alpha = u_2 * 2.0 * Math.PI;

    return [r,alpha];
}

The points are focused in the center! With uniform area sampling it looks like the following:

function uniform_circle(){

    const u_1 = Math.random();
    const u_2 = Math.random();
    
    const r = Math.sqrt(u_1);
    const alpha = u_2 * 2.0 * Math.PI;

    return [r,alpha];
}

Algorithm:

Generate points uniformly distributed over the area of a disk with radius $R$

Choose uniform numbers $u_1, u_2$ in $[0,1]$
Calculate $r = \sqrt{u_1}R$
Calculate $\alpha = 2\pi u_2$

Optionally convert to cartesian coordinates with:

\begin{aligned} x &= r\cos\alpha \\ y &= r\sin\alpha \end{aligned}

PDF Used In Monte Carlo Integration

\begin{aligned} p(\alpha) &= \frac{1}{\pi R^2} \\ p(r,\alpha) &= \frac{r}{\pi R^2} \end{aligned}

Derivation:

We want our distribution to be uniform in area over the disk. We will therefore integrate over the entire area $A$ .

The probability distribution function (pdf) has to integrate to $1$ over the complete region. Since we want it to be uniform, it will just be some constant for each area element => $p(a)=c$

\int_{A}p(a)da = \int_{A}cda = c \int_{A}da = 1

The easy way:

We know that the area of a disc is $\pi R^2$ .

c \int_{A}da = 1 \Rightarrow c\pi R^2 = 1 \Rightarrow c = \frac{1}{\pi R^2}

The complicated way:

An area element $da$ can be decomposed as $da=dxdy$ (the very small area is just a rectangle). We then define a mapping from polar to cartesian coordinates:

x = r\cos\alpha = f_x(r, \alpha)\frac{\partial f_x(r, \alpha)}{\partial r}

y = r\sin\alpha = f_y(r, \alpha)

To do a coordinate transformation we have to find the determinant of the Jacobian of this transformation.

\begin{aligned} J &= \begin{pmatrix} \frac{\partial f_x(r, \alpha)}{\partial r}& \frac{\partial f_x(r, \alpha)}{\partial \alpha} \\ \frac{\partial f_y(r, \alpha)}{\partial r} & \frac{\partial f_y(r, \alpha)}{\partial \alpha} \end{pmatrix} \\ &= \begin{pmatrix} \cos\alpha & -r\sin\alpha \\ \sin\alpha & r\cos\alpha \end{pmatrix} \end{aligned}

The determinant is thus:

\det{J} = r\cos\alpha \cos\alpha - (-r\sin\alpha \sin\alpha) = r\cos^2\alpha + r\sin^2\alpha = r(\cos^2\alpha + \sin^2\alpha) = r

Computing the integral from before:

\begin{aligned} c \int_{A}da &= c \int_{A}dxdy \\ &= c \int_0^{R}\int_0^{2\pi} rdrd\alpha \\ &= c \int_0^{R}rdr\int_0^{2\pi} d\alpha \\ &= c \int_0^{R}rdr [\alpha]_0^{2\pi} \\ &= c \int_0^{R}rdr \\ &= 2 c \pi \int_0^{R}rdr \\ &= 2 c \pi [\frac{r^2}{2}]_0^R \\ &= c\pi R^2 \end{aligned}

Setting this result equal to $1$ and solving for $c$ will give the desired result

Marginal And Conditional Densities:

We have $p(a)=\frac{1}{\pi R^2}$ . With the transformation from the last step we get $p(r,\alpha)=\frac{r}{\pi R^2}$ (Transform into polar angles by adding the factor $r$ )

First the marginal density:

\begin{aligned} p_r(r) &= \int_0^{2\pi} p(r,\alpha)d\alpha\\ &= \int_0^{2\pi} \frac{r}{\pi R^2}d\alpha \\ &= \frac{r}{\pi R^2}\int_0^{2\pi} d\alpha \\ &= \frac{r}{\pi R^2} 2\pi \\ &= \frac{2r}{R^2} \end{aligned}

The conditional probability:

p(\alpha \vert r) = \frac{p(r,\alpha)}{p_r(r)} = \frac{rR^2}{\pi R^2 2r} = \frac{1}{2\pi}

Compute CDFs

P_r(r) = \int_0^r p_r(r')dr' = \int_0^r \frac{2r'}{R^2}dr' = \frac{2}{R^2}[\frac{r^2}{2}]_0^r = \frac{r^2}{R^2}

P(\alpha \vert r) = \int_0^\alpha p(\alpha' \vert r)d\alpha' = \int_0^\alpha \frac{1}{2\pi}d\alpha' = \frac{\alpha}{2\pi}

Invert CDFs

\begin{aligned} u_1 &= P_r(r) \\ &= \frac{r^2}{R^2} \\ u_1 R^2 &= r^2 \\ \sqrt{u_1R^2} &= r \\ \sqrt{u_1}R &= r \end{aligned}

\begin{aligned} u_2 &= P(\alpha \vert r) \\ &= \frac{\alpha}{2\pi} \\ 2\pi u_2 &= \alpha \end{aligned}