The number of successes in a certain number of trials

In risk analysis, we often attempt to predict the results of a set of random trials, where the trials can result in either a 'success' (the outcome we are most interested in) or a 'failure'. For example:

• How many airplane flights will result in crashes;

• How many people eating hamburgers will get E.coli infections;

• How many sales pitches will result in a sale;

• How many people entering a shop will make a purchase;

• How many women receiving fertility treatment will fall pregnant;

• How many cars will need a replacement engine within the guarantee period; etc.

Independent trials

If each of these n trials is independent (meaning the result of each trial is not influenced by the result of any previous trial), and if all trials have the same probability of success, the outcome of these trials conforms to a binomial process. Moreover, we can model the number of successes s using a Binomial distribution:

s = VoseBinomial(n,p)

A binomial process is a random counting system where there are n independent identical trials, each one of which has the same probability of success p, which produces s successes from those n trials (where 0≤s≤n and n > 0 obviously). There are thus three parameters {n, p, s} that between them completely describe a binomial process.

To model the number of successes in a certain number of trials we will use:

s = VoseBinomial(n,p)

where n is the number of trials and p is the probability of the trial becoming a success.

The simplest example of a binomial process is the toss of a coin. If I toss a fair coin (a coin with 50% probability of returning either head or tail) 10 times, what is the distribution of the number of heads I will get? The answer to this question can be modelled using just one formula:

s = VoseBinomial(10,0.5) , which will produce the following outcome distribution:

As expected, the most likely value and the mean number of successes equal 5.

Problem: What is the probability that I will have no heads at all?

As we see from the graph above, this probability is very low, so is better not determined by simulation as we would need very many iterations to get an accurate answer. However, the probability could be calculated using the VoseBinomialProb function:

P(s=0) = VoseBinomialProb(0,10,0.5,0) = 0.000976563

Problem: What is the distribution of the maximum number of heads I can get in a row by tossing a fair coin 10 times?

The solution to this problem is provided in the example model - Coins.

Hypergeometric process

The hypergeometric process occurs when one is sampling randomly without replacement from some population, and where one is counting the number in that sample that have some particular characteristic. In this situation we have four parameters: M, the population size; D, the number of successes in the population; n, the sample size; and s, the number of successes in the sample. We do not have a probability of success parameter (p) here since it is changing because the proportion of successes in the remaining population is changing with every sample we take.

Let's imagine a pack of playing cards (52 cards, without jokers), and 13 of them are hearts. If we consider a card of hearts to be a success, then we have: M = 52 and D = 13. If we are to pick a single card from the pack at random, the probability of picking a heart is equal to 13/52 = 1/4.

What is the probability of picking two hearts in a row? Well, unlike in the binomial process, we cannot just multiply the probability of 1/4 by itself, since the probability of the second card being a heart depends on the suit of the first card:

1. if the first card was a heart (a success), then there are only 12 hearts remaining in the pack and the probability of picking another heart reduces to 12/51

2. if the first card was not a heart (a failure), then all 13 hearts are still in the pack and the probability of picking a heart increases to 13/51

Excel's HYPGEOMDIST(s,n,D,M) function calculates the probability of picking 2 hearts from two trials:

P(2 hearts) = HYPGEOMDIST(2,2,13,52) = 0.058823529

In general, direct calculation of a probability is to be preferred over simulation because it is faster of more accurate, however, probability problems quickly become far too complex for us to calculate and we then resort to simulation.

Let's now consider the following problem: I have three full packs of cards, and I draw 10 cards from each pack. What is the probability that I draw at least 10 hearts in total? This problem is more complicated as there are many combinations of samples from each pack that would give the required result, so we resort to simulation.

The solution to this problem is provided in the example model   cards

Clustered trials

Sometimes, the trials are not independent of each other but grouped together in fixed, or variable, sized groups, or clusters. For example, the number of airline passengers that might die in a year from a plane crash is strongly grouped because they are very likely to suffer the same fate if they are in an aircraft together that crashes. Other examples include:

• How many people in a village get divorced (they are paired) in a year

• How many infected blood samples are mis-diagnosed by a laboratory (if, for example, a lab tests samples in batches and makes a mistake with a batch);

• How many manufactured items in a consignment fail to meet the required tolerance (if, for example, a machine is not set up correctly for a production run).

In such situations, one can often model the group using a Binomial distribution. Then, if the number of individuals in the group is constant (say takes a value k), the number of successes is =k*VoseBinomial(n,p). So, for example, if we had 290 married couples in the village, and believed a married couple had a 3% probability of divorcing in a year, we would estimate seeing 2*VoseBinomial(290,3%) divorcees next year.

Alternatively, if the number in a group is variable, we need to create a model that sums a variable number of random variables. For example, imagine that we consider that there is a 5% chance of incorrectly setting up a machine to produce widgets resulting in a bad batch of out of tolerance widgets, and that we set up 10 machines on a day's production run, but that each machine will produce Poisson(250) widgets, what fraction of my production will be out of tolerance?

Model widgets provides the answer.

Probability randomly varies for the set of trials as a whole

If the probability of a trial becoming a success is a random variable itself, the resultant distribution of the number of successes is wider than a Binomial distribution. For example, the probability that the sheep in a flock will survive the winter depends on whether the winter is particularly harsh. All of the sheep will endure the same conditions, so if the probability is higher for one it is higher for all, so you could argue that the probability is variable, but the same for each sheep. A convenient way of modelling a probability that varies is to use the Beta distribution. In ModelRisk, we could then model the number of successes as:

s = VoseBinomial(n,VoseBeta(a,b))

where a and b are two parameters used to create the required shape for the Beta distribution. In fact, this distribution for s is known in probability theory as the Beta-Binomial distribution. The Beta-Binomial distribution always has a greater spread than its most closely matching Binomial distribution.

But how about scenarios where the probability of a trial's success is a random variable but each trial's probability is independent of the others? You might think that we would need to model the probability for each trial separately as a random variable, for example as shown in the following spreadsheet: indepprob.

But, in fact, we just need to use the following formula:

s = VoseBinomial(n, a/(a+b))

where a/(a+b) is the mean of a Beta(a,b) distribution. The above model runs the calculation both ways for you to compare, but can you work out why this is true?