1
$\begingroup$

I want to generate numbers from a 2D "cylindrical" distribution, whose PDF would be:

$f(x,y) = \begin{cases}1/(45^2 \pi) & x^2+y^2 \le 45\\0 & \text{otherwise}\end{cases}$

I don't know if it is still called a uniform distribution since it is in 2D with the domain restricted to $x^2+y^2 \le 45$, so it is not rectangular.

I wrote the function dist[x_, y_] := If[Sqrt[(0 - x)^2 + (0 - y)^2] <= 45, 1./(Pi*45^2), 0]

and I confirmed that Integrate[dist[x, y], {x, -Infinity, Infinity}, {y, -Infinity, Infinity}] == 1.

but I have no idea how to generate random $(x,y)$ couples from that point.

I'd like to use RandomVariate but it cannot be inputted with user-made PDFs.

Thanks !

$\endgroup$
  • 2
    $\begingroup$ Why not RandomPoint@Disk[{0, 0}, 45]? $\endgroup$ – corey979 Jan 15 '17 at 7:13
  • $\begingroup$ Seems I don't have access to RandomPoint in Home Edition v.10.0.2.0 $\endgroup$ – Domiinic Jan 15 '17 at 7:28
6
$\begingroup$

You can take a random angle and a uniform random number between 0 and $45^2$ in the following manner:

n = 1000;
θ = RandomVariate[UniformDistribution[{0, 2 π}], n];
r = (RandomVariate[UniformDistribution[{0, 45^2}], n]^0.5);

x = r  Cos[θ];
y = r  Sin[θ];

ListPlot[Transpose[{x, y}], AspectRatio -> 1]

Uniform samples on a circle

$\endgroup$
  • $\begingroup$ It works. Thanks ! $\endgroup$ – Domiinic Jan 15 '17 at 8:02
3
$\begingroup$

One can always do rejection sampling:

With[{n = 1000, r = 45},
     k = 0;
     pts = Reap[While[cand = RandomReal[{-r, r}, 2]; 
                      If[cand.cand < r^2, Sow[cand]; k++]; k < n]][[-1, 1]]];

but if you are unwilling or unable to do this, use the results from here:

With[{n = 1000, r = 45},
     pts = Table[With[{v = RandomVariate[NormalDistribution[], 2], 
                       y = -Log[RandomReal[]]}, r v/Sqrt[v.v + y]], {n}];]

Both approaches easily extend to general $n$-dimensional ball sampling.

$\endgroup$
3
$\begingroup$

Here's another way:

ran = RandomPoint[Disk[{0, 0}, Sqrt[45]], 10^3];

We can check the condition:

Max[(Norm /@ ran)^2] <= 45

which returns True.

$\endgroup$
  • $\begingroup$ The OP noted in the comments that his version does not yet have RandomPoint[], but it is nevertheless good to have it here as an explicit answer. $\endgroup$ – J. M. is away Jan 27 '17 at 3:41
  • $\begingroup$ @J.M. I apologize for not reading carefully. And I'm glad my answer is useful nevertheless. $\endgroup$ – mef Jan 27 '17 at 12:30

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.