3
$\begingroup$

Is there any build in function or shorter version of this? Assume I have 3 reactions time and each reaction occurs randomly. Generate rand=RandomReal[] if rand<0.2 then result is 1, if rand<0.5 then result is 2 and if rand<1 then result is 3. Thanks..

 reaction = {0.2, 0.5, 1};
 First@Flatten@Position[reaction, SelectFirst[reaction, RandomReal[] < # &]]
$\endgroup$
6
  • $\begingroup$ Do you mean With[{q = RandomReal[]}, First@Flatten@Position[reaction, SelectFirst[reaction, q < # &]]] ` $\endgroup$
    – Coolwater
    Dec 27, 2017 at 18:36
  • $\begingroup$ Not really. Maybe MapIndex or Ordering $\endgroup$ Dec 27, 2017 at 18:41
  • $\begingroup$ wd=WeightedData[Range[3],Differences[{0, ##&@@reaction}]]; First@RandomSample[wd,1]? $\endgroup$
    – kglr
    Dec 27, 2017 at 18:56
  • $\begingroup$ shorter: First@FirstPosition[reaction,_?(#>RandomReal[]&)]? $\endgroup$
    – kglr
    Dec 27, 2017 at 19:03
  • $\begingroup$ @kglr I can accept your solution if there is no shorter one. Please write as an answer. Thanks. $\endgroup$ Dec 27, 2017 at 19:07

4 Answers 4

3
$\begingroup$
wd = WeightedData[Range[3], Differences[{0, ## & @@ reaction}]]; 

SeedRandom[1]
RandomChoice[wd] 

3

RandomChoice[wd, 5] 

{1,3,1,2,1}

Update: a shorter alternative using FirstPosition

SeedRandom[123]
FirstPosition[reaction - RandomReal[], _?Positive][[1]]

2

Note: This function uses a single random number while the alternative

FirstPosition[reaction, _?(#>RandomReal[]&)][[1]]

creates a new random number for every comparison until the condition is satisfied (as does the code in OP), and hence, it can give a different result

SeedRandom[123]
FirstPosition[reaction, _?(#>RandomReal[]&)][[1]]

3

$\endgroup$
2
  • $\begingroup$ I'll go with this since it is more readable for me. First@FirstPosition[reaction,_?(RandomReal[]<#&)] $\endgroup$ Dec 27, 2017 at 21:22
  • $\begingroup$ First@FirstPosition[reaction,_?(RandomReal[]<#&)] does not work in Do loop but this works perfectly FirstPosition[reaction - RandomReal[], _?Positive][[1]] see mathematica.stackexchange.com/questions/165073/… $\endgroup$ Apr 28, 2018 at 7:46
7
$\begingroup$

Assuming you meant to only use one RandomReal call (per @Coolwater's comment), you can use RandomChoice:

RandomChoice[{.2, .3 ,.5} -> {1, 2, 3}]
RandomChoice[{.2, .3 ,.5} -> {1, 2, 3}, 10]

3

{3, 1, 3, 2, 3, 2, 3, 2, 3, 1}

Addendum

To mimic the OP code with RandomChoice, by changing the random number for each comparison, you just need to adjust the probabilities. 20% of the time, the first random number is less than .2. Of the other 80%, 50% of the time the second random number is less than .5, and 50% of the time it is greater than .5. So, the following should produce the same distribution as the OP function:

RandomChoice[{.2, .4, .4} -> {1, 2, 3}]

We can check by running the OP code multiple times:

reaction = {0.2, 0.5, 1};
data = Table[First@Flatten@Position[reaction, SelectFirst[reaction, RandomReal[] < # &]], {10^5}];

Counts[data]/10.^5

<|1 -> 0.20025, 2 -> 0.4018, 3 -> 0.39795|>

$\endgroup$
11
  • $\begingroup$ What I want is generate random number rand=RandomReal[] and look at its location in reaction and extract index of reaction. $\endgroup$ Dec 27, 2017 at 18:45
  • $\begingroup$ @OkkesDulgerci I think that's exactly what my function does. 20% of the time, the random number is less than .2, 30% of the time the random number is between .2 and .5 and 50% of the time it is more than .5. $\endgroup$
    – Carl Woll
    Dec 27, 2017 at 18:53
  • $\begingroup$ RandomChoice[{.2, .3 ,.5} -> {1, 2, 3}] might be equivalent what I am doing. I will think about it. $\endgroup$ Dec 27, 2017 at 18:54
  • $\begingroup$ When I replaced your code with mine it does not give desired result. $\endgroup$ Dec 27, 2017 at 19:06
  • 1
    $\begingroup$ @OkkesDulgerci - your code (generating a random number for each test done by SelectFirst) is equivalent to RandomChoice[{0.2, 0.4, 0.4} -> Range[3]] since the first test extracts ~20%, if the second test is needed it extracts ~50% of the remaining ~80% (i.e., ~40%), and the last test will extract all of the remaining. Using a single random number gives the result that Carl provides RandomChoice[{0.2, 0.3, 0.5} -> Range[3]] $\endgroup$
    – Bob Hanlon
    Dec 27, 2017 at 19:12
1
$\begingroup$

Another possibility is to make this a distribution and then use all the machinery that goes with that:

SeedRandom["December 27, 2017"];
distReaction = With[
    {
       reaction = {0.2, 0.5, 1}
    },
    EmpiricalDistribution[ Differences[{0} ~ Join ~ reaction ] -> Range[3] ]
]


RandomVariate[ distReaction, 5 ]

{2, 1, 3, 3, 3}

We could then ask for the probability of this specific result assuming that the results are i.i.d.:

Probability[ 
    {x1, x2, x3, x4, x5} == {2, 1, 3, 3, 3}, 
    Thread [ {x1, x2, x3, x4, x5} \[Distributed] distReaction ]
]

0.0075

$\endgroup$
1
$\begingroup$

The first response that came to mind, after reading the question, was

BlockRandom[
 With[{rand = RandomReal[]},

  Which[
   rand <= 0.2, reaction[[1]],
   0.2 < rand <= 0.5, reaction[[2]],
   0.5 < rand <= 1, reaction[[3]]]], RandomSeeding -> 123456798]

where, obviously, reaction = {0.2, 0.5, 1}. By replacing reaction[[i]] with i it is easy to obtain the index of the reaction (not the time it takes).

I think that, this is as literal an implementation as it gets.

A couple of other approaches are

BlockRandom[
 With[{rand = RandomReal[]},

  Extract[
   reaction,
   First[Position[Thread[rand <= reaction], True]]]
  ], RandomSeeding -> 123456798]

and

BlockRandom[
 With[{rand = RandomReal[]},

  Last[Flatten[Reap[
     Scan[
      If[rand <= #, Sow[#]; Return[]] &, reaction]]]]
  ], RandomSeeding -> 123456798]

In the former of the last two approaches, the index of the reaction instead of it's time can be obtained by replacing Extract[<>] with First[Flatten[Position[Thread[rand <= reaction], True]]].

The later code segment cannot be readily adjusted to obtain the index of the reaction instead of the time (at least not without doing something like First@Position[reaction, x] where x stands for the output of the last approach)

$\endgroup$

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.