Suppose I have 8 possible independent things to happen this weekend. what are the odds that five (or more) of them happening? and how do you calculate it?
Est. Success
Scenario Rate
A … 95%
B … 40%
C … 30%
D … 10%
E … 5%
F … 15%
G … 6%
H … 8%
Have to ask, is this a homework question?
No, I am trying to figure out odds of a Tiger Woods falling out of the top 50 in the World Golf rankings this week. at least 5 of 8 possible scenarios have to happen and these are the percentages that I have assigned to the success rate of any one thing happening. I know how to do it if all of them happen and none of them happen. but not in between.
its actually a little more complicated than what I have summarized above. Because Scenario F and G can happen two or three times if the stars align and each event “counts” as a success. but if I have formula, then I can calculate it.
There are 8!/(5! * 3!) + 8!/(6! * 2!) + 8!/(7!*1!) + 8!/(8! + 0!) = 93 different ways for 5 or more independent events to occur out of 8 trials. Should be relatively simple to write a program to enumerate those 93, calculate the probability of each, and sum up the probabilities.
Since they joined in 2000, it’s very doubtful!
Whether it’s homework or not that looks like a messy problem. Given the number of combinations I’d say you would end up with a sum of 93 terms. Good luck with the calculator…
ETA: What P just mentioned above.
Would it be simpler to ask “what are the odds of 4 of them **not **happening”?
I don’t think it’s quite that simple, because you could end up with a probability greater than one if you just summed the probabilities of all 93 cases. I think you have to find 1 - [ (1 - p[sub]1[/sub]) x (1 - p[sub]2[/sub]) x … x (1 - p[sub]93[/sub]) ]
[edit] Actually, that’s probably rubbish. They can’t sum to more than one.
It’s not all that difficult to compute numerically. I get that it’s about 9.64 x 10^-3.
rm(list = ls())
probs <- c(.95, .4, .3, .1, .05, .15, .06, .08)
d <- c(0, 1)
events <- expand.grid(d, d, d, d, d, d, d, d)
a <- apply(events, 1, sum) > 4
b <- rep(0, nrow(events))
for (i in 1:nrow(events))
{
b* <- prod(probs^events[i, ] * (1 - probs)^(1 - events[i, ]))
}
p <- t(a) %*% b
Well, I could have a kid (but I don’t). thanks for the help.
0.964% chance? Less than 1%? That sounds low to me. I thought the percentage would be somewhere around 5-10%.
Out of curiosity, what language is that, ultrafilter?
I agree with ultrafilter’s answer. There’s another more efficient approach (quadratic in the number of events, instead of exponential) to reach this answer, based on generating functions. The idea is to represent a random Bernoulli event (biased coin flip) with success probability p by a polynomial
(1-p) + p*x .
Now multiplying two of these polynomials together gives a quadratic polynomial; it’s easy to see that the constant term is the probability that both events fail, the linear coefficient is the probability that exactly one of the two events succeeds, and the quadratic coefficient is the probability that both events succeed. More generally, if you multiply a bunch of these polynomials together, the coefficient of the x[sup]k[/sup] term is just the probability that exactly k of the events succeed. So all you have to do is multiply out the eighth-order polynomial
(0.05+0.95x)(0.60+0.40x)(0.70+0.30x)(0.90+0.10x)(0.95+0.05x)(0.85+0.15x)(0.94+0.06x)(0.92+0.08x)
and sum the coefficients of the x[sup]5[/sup], …, x[sup]8[/sup] terms.
Oh, that’s awesome, Omphaloskeptic.
Just gut feel.
I ran a simulation in Excel (10000x) and and it appears your number is correct.
Many thanks. Appreciate it.
It looks like Omphaloskeptic has come up with the elegant formula I was looking for. But I am having a major brain fart. I don’t see how it comes up with 0.964% chance. Can someone step me through it?
Here’s a Perl program that performs Omphaloskeptic’s method that I’ll step you through. I put in a print statement in the middle so that you can see how the polynomial builds up. At each step, the first number, index 0, is the chance that no successes happen. So the initial value there is 1, since with no trials yet, there’s a 100% chance that there are no successes. As more and more values get put into @poly, values at index 1 are the probability that exactly one success will have been seen, at index 2 is the probability that exactly two successes have been seen, etc. Then we sum up all of the coefficients with an index greater than or equal to 5.
#!/usr/bin/perl -w
use strict;
my @poly = ( 1 );
foreach my $probability ( 0.95, 0.4, 0.3, 0.1, 0.05, 0.15, 0.06, 0.08 )
{
my @trialFail = (map { $_ * (1 - $probability) } @poly, 0);
my @trialSuccess = (0, map { $_ * $probability } @poly);
@poly = map { $_ + shift @trialSuccess } @trialFail;
print "\@poly = (@poly)
";
}
# sum up the coefficients for index >= 5
my $sum = 0;
foreach my $probability (@poly[5 .. $#poly])
{
$sum += $probability;
}
print "$sum
"
So, first @poly = ( 1 )
The first time through the loop, @trialFail = ( 0.05, 0 ) and @trialSuccess = (0, 0.95). The third line sums them up, setting @poly = (0.05, 0.95).
The second time through the loop, @trialFail = (0.03, 0.57, 0), @trialSuccess = (0, 0.02, 0.38), so their sum is (0.03, 0.59, 0.38).
The third time through the loop, @trialFail = (0.021, 0.413, 0.266, 0), @trialSuccess = (0, 0.009, 0.177, 0.114), so their sum is (0.021, 0.422, 0.443, 0.114).
For the following loops, we see, in order:
(0.0189 0.3819 0.4409 0.1469 0.0114)
(0.017955 0.36375 0.43795 0.1616 0.018175 0.00057)
(0.01526175 0.31188075 0.42682 0.2030525 0.03968875 0.00321075 8.55e-05)
(0.014346045 0.29408361 0.419923645 0.21647855 0.049490575 0.00539943 0.000273015 5.13e-06)
(0.0131983614 0.2717046048 0.4098564422 0.2327541576 0.062849613 0.0089267216 0.0006831282 2.65608e-05 4.104e-07)
We then chop off the entries from index 0 to index 4 to get:
(0.0089267216 0.0006831282 2.65608e-05 4.104e-07)
…which then sum to 0.009636821, i.e., 0.9637%
Thanks again. I need to wrap my head around this post when I am not so tired.