How to do a number of random trials efficiently?

Question

Let say an event have the probability P to succeed. (0 < P < 1 ) and I have to make N tests to see if this happens and I want the total number of successes:

I could go

int countSuccesses = 0;
while(N-- > 0)
{
   if(Random.NextDouble()<P) countSuccesses++; // NextDouble is from 0.0 to 1.0
}

But is there not a more efficient way to do this? I want to have a single formula so I just can use ONE draw random number to determine the total number of successes. (EDIT The idea of using only one draw was to get below O(n))

I want to be able call a method

GetSuccesses( n, P)

and it to be O(1)

UPDATE
I will try to go with the MathNet.Numerics.Distributions.Binomial.Sample(P, n) even if it might be using more then only one random number I'll guess it will be faster than O(n) even if its not O(1). I'll benchmark that. Big thanks to David and Rici.

UPDATE
The binomial sample above was O(n) so it did not help me. But thanks to a comment done by Fred I just switched to
MathNet.Numerics.Distributions.Normal.Sample(mean, stddev) where
mean = n * P
stddev = Math.Sqrt(n * P * (1 - P)); and now it is O(1) !

Answer 1

Per @rici for small N you can use the CDF or PMF of the Binomial Distribution, and simply compare the random input with the probabibilities for 0,1,2..N successes.

Something like:

  static void Main(string[] args)
    {

        var trials = 10;
        var trialProbability = 0.25;
        for (double p = 0;  p <= 1; p += 0.01)
        {
            var i = GetSuccesses(trials, trialProbability, p);
            Console.WriteLine($"{i} Successes out of {trials} with P={trialProbability} at {p}");
        }
        Console.ReadKey();

    }
    static int GetSuccesses(int N, double P, double rand)
    {
        for (int i = 0; i <= N; i++)
        {
            var p_of_i_successes = MathNet.Numerics.Distributions.Binomial.PMF(P, N, i);
            if (p_of_i_successes >= rand)
                return i;

            rand -= p_of_i_successes;

        }
        return N;


    }

Answer 2

I'm not going to write formula here, as it's already in wiki, and I don't really know good formatting here for such things.

Probability for each outcome can be determined by Bernulli formula https://en.wikipedia.org/wiki/Bernoulli_trial

What you need to do is to calculate binominal coefficient, then probability calculation becomes quite easy - multiply binominal coefficient by p and q in appropriate powers. Fill in array P[0..n] that contains probability for each outcome - number of exactly i successes.

After set up go from 0 up to n and calculate rolling sum of probabilities. Check lower/upper bounds against random value and once it's inside current interval, return the outcome.

So, deciding part will be like this:

sum=0;
for (int i = 0; i <= n; i++)
  if (sum-eps < R && sum+P[i]+eps > R)
    return i;
  else
    sum+=P[i];

Here eps is small floating point value to overcome floating point rounding issues, R is saved random value, P is an array of probabilities I mentioned before.

Unfortunately, this method is not practical for big N (20 or 100+):

• you'll get quite big impact of rounding errors

• random numbers generator can be not determinitive enough to cover every possible outcome with proper probabilities distribution

Answer 3

Based on how you've phrased the question, this is impossible to do.

Your are essentially asking how to ensure that a single coin flip (i.e. one randomized outcome) is exactly 50% heads and 50% tails, which is impossible.

Even if you were to use two random numbers, where you expect one heads and one tails; this test would fail in 50% of all cases (because you might get two heads or two tails).

---

Probablility is founded on the law of large numbers. This explicitly states that a small sampling cannot be expected to accurately reflect the expected outcome.

The LLN is important because it guarantees stable long-term results for the averages of some random events. For example, while a casino may lose money in a single spin of the roulette wheel, its earnings will tend towards a predictable percentage over a large number of spins. Any winning streak by a player will eventually be overcome by the parameters of the game. It is important to remember that the law only applies (as the name indicates) when a large number of observations is considered. There is no principle that a small number of observations will coincide with the expected value or that a streak of one value will immediately be "balanced" by the others (see the gambler's fallacy).

---

When I asked this as a comment; you replied:

@Flater No, I am making N actual draws but with only one random number.

But this doesn't make sense. If you only use one random value, and keep using that same value, then every draw is obviously going to give you the exact same outcome (that same number).

---

The closest I can interpret your question in a way that is not impossible would be that you were mistakenly referring to a single random seed as a single random number.

A random seed (or seed state, or just seed) is a number (or vector) used to initialize a pseudorandom number generator.

For a seed to be used in a pseudorandom number generator, it does not need to be random. Because of the nature of number generating algorithms, so long as the original seed is ignored, the rest of the values that the algorithm generates will follow probability distribution in a pseudorandom manner.

However, your explicitly mentioned expectations seem to disprove that assumption. You want to do something like:

GetSuccesses( n, P, Random.NextDouble())

and you're also expecting to get a O(1) operation, which flies in the face of the law of large numbers.

If you're actually talking about having a single random seed; then your expectation are not correct.

• If you make N draws, the operation is still of O(N) complexity. Whether you randomize the seed after every draw or not is irrelevant, it's always O(N).
• GetSuccesses( n, P, Random.NextDouble()) is going to give you one draw, not one seed. Regardless of terminology used, your expectation of the code is not related to using the same seed for several draws.

---

As the question is currently phrased; what you want is impossible. Repeated comments for clarification by several commenter have not yet yielded a clearer picture.

As a sidenote, I find it very weird that you've answered to every comment except when directly asked if you are talking about a seed instead of a number (twice now).

Answer 4

@David and @rici pointed me to the
MathNet.Numerics.Distributions.Binomial.Sample(P, n)
A benchmark told me that it also O(n) and in par with my original

int countSuccesses = 0;
while(N-- > 0)
{
   if(Random.NextDouble()<P) countSuccesses++; // NextDouble is from 0.0 to 1.0
}

But thanks to a comment done by Fred:

You could turn that random number into a gaussian sample with mean N*P, which would have the same distribution as your initial function

I just switched to
MathNet.Numerics.Distributions.Normal.Sample(mean, stddev) where
mean = n * P
stddev = Math.Sqrt(n * P * (1 - P)); and now it is O(1) !

and the function I wanted got to be:

    private int GetSuccesses(double p, int n)
    {
        double mean = n * p;
        double stddev = Math.Sqrt(n * p * (1 - p));
        double hits = MathNet.Numerics.Distributions.Normal.Sample(Random, mean, stddev);
        return (int)Math.Round(hits, 0);
    }

As Paul pointed out, this is an approximation, but one that I gladly accept.