6
$\begingroup$

Random numbers in a range[pi,pf] can be generated using

p = RandomReal[{pi, pf}];

How can I generate a random number from this range excluding the range[a,b] where a>pi and b<pf?

$\endgroup$

5 Answers 5

6
$\begingroup$

Generate a real number uniformly from a range with the total width of each possible range (that is, [pi, a] + [b, pf]). Then map the value from this range to the split range of interest:

f[pi_, pf_, a_, b_] := 
  If[# > a - pi, # + b - a, #] &[RandomReal[{0, pf - b + a - pi}]] + 
   pi;

The random number is generated in RandomReal[{0, pf - b + a - pi}].

The re-mapping occurs in If[# > a - pi, # + b - a, #] &[...] + pi. This could also be written as If[# > a - pi, # + b - a + pi, # + pi] &[...] for clarity. The original re-mapping can be read as "if the generated value is greater than a - pi (placing it after the split), then increase its position by the width of the split (b - a), otherwise keep it the same, then add pi to it.

$\endgroup$
0
14
$\begingroup$

You can use RandomPoint with Interval:

First @ RandomPoint[Interval[{1, 3}, {5, 7}]]

5.23755

If you want to generate multiple points:

RandomPoint[Interval[{1, 3}, {5, 7}], 10][[All, 1]]

{6.29151, 5.69004, 1.61121, 6.25677, 5.87869, 2.44657, 6.98119, 2.19049,
2.63382, 1.08334}

$\endgroup$
1
  • $\begingroup$ (+1) This is a nice solution! Thanks for sharing that. $\endgroup$
    – Shasa
    Commented Jun 9, 2022 at 17:43
7
$\begingroup$

Another approach:

(* Set parameters *)
pi = 1;
pf = 7;
a = 3;
b = 5;

(* Set sample size *)
n = 1000000;

(* Generate a random indicator variable that decides which section to sample from *)
(* 1 -> from pi to a and 0 -> from b to pf *)
α = RandomVariate[BernoulliDistribution[(a - pi)/(pf - pi - (b - a))], n];

(* Generate a set of uniform (0,1) numbers *)
x = RandomVariate[UniformDistribution[{0, 1}], n];

(* Combine α and x to get a random sample from the two segments *)
y = α (pi + x (a - pi)) + (1 - α) (b + x (pf - b));

Histogram[y, "FreedmanDiaconis", "PDF"]

Histogram

$\endgroup$
3
  • 3
    $\begingroup$ (Of course, RandomVariate[UniformDistribution[{0, 1}], n] is functionally equivalent to RandomReal[1, n].) $\endgroup$ Commented Jun 9, 2022 at 20:34
  • 3
    $\begingroup$ @J.M. One of my many quirks is to avoid RandomReal because (to me) that doesn't explicitly imply a uniform distribution as explicitly as stating Uniform[{0,1}]. Part of this comes from cringing when folks add "random noise" to a regression using RandomReal when that is rarely appropriate (using RandomReal, that is, not the cringing). If there was some great added efficiency in using RandomReal, then that would be different. $\endgroup$
    – JimB
    Commented Jun 10, 2022 at 5:31
  • $\begingroup$ Oh, you're certainly right on the abuse of RandomReal[] in that regard. From a readability viewpoint, explicitly invoking a uniform distribution makes sense; testing the overhead looks like an interesting experiment... $\endgroup$ Commented Jun 10, 2022 at 5:48
7
$\begingroup$

I include this answer for completeness, not because it's better than any of the other answers.

You could use "rejection sampling." For the stated problem, rejection sampling amounts to the following simple algorithm: Draw from the full domain until you get a draw from the restricted domain and return that draw.

Here's an implementation:

ff = 
Function[{pi, pf, a, b},
  Module[{ran},
    While[
      ran = RandomVariate[UniformDistribution[{pi, pf}]];
      a < ran < b
    ];
    ran
  ]
]

However, if $(\text{b}-\text{a})$ is a large fraction of $(\text{pf}-\text{pi})$ then the sampler will not be very efficient in terms of the number of draws it must make from the full domain to generate a draw from the restricted domain.

Of course, one can make the function run faster by compiling it (about 30 times faster):

cf = 
Compile[{pi, pf, a, b},
  Module[{ran = 1.},
    While[
      ran = RandomVariate[UniformDistribution[{pi, pf}]];
      a < ran < b
    ];
    ran
  ]
]

(Compile complains if ran is not initialized.)

$\endgroup$
1
  • 3
    $\begingroup$ We were lucky here that there were (clever!) ways to assemble the required distribution from built-in distributions, as featured in other answers; otherwise, rejection sampling is really one of the most generally applicable methods for getting variates, so it is useful to have this here for completeness. $\endgroup$ Commented Jun 10, 2022 at 17:48
5
+50
$\begingroup$

Here is an approach when there are multiple subintervals within uniform distribution.

(* Set intervals *)
intervals = {{1, 3}, {5, 7}, {8, 11}};

(* Get proportions of distribution associated with each interval *)
p = ((#[[2]] - #[[1]]) & /@ intervals)/(Total[intervals[[All, 2]] - intervals[[All, 1]]])
(* {2/7,2/7,3/7} *)

(* Sample size *)
n = 1000000;

(* Select a intervals randomly *)
SeedRandom[12345];
α = RandomVariate[MultinomialDistribution[n, p], 1][[1]]
(* {285630,285547,428823} *)

(* Combine interval information and number of samples from each interval *)
sampling = {#[[1, 1]], #[[1, 2]], #[[2]]} & /@ Transpose[{intervals, α}]
(* {{1,3,285630},{5,7,285547},{8,11,428823}} *)

(* Generate observations *)
x = RandomVariate[UniformDistribution[{#[[1]], #[[2]]}], #[[3]]] & /@ sampling // Flatten;

(* Show histogram of observations *)
Histogram[x, "FreedmanDiaconis", "PDF"]

Histogram of observations

$\endgroup$
1
  • $\begingroup$ (+1) This is a nice answer. Thanks for sharing that. $\endgroup$
    – Shasa
    Commented Jun 10, 2022 at 5:43

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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