Statistics: How your intuition may not be correct!

I was playing Terraria with some friends, and one of my friends went for a fishing mini-game. Basically, he had a 10% chance to get the fish he wanted. After a dozen or so attempts, he was still unsuccessful. He then complained to me

"I have a 10% chance to get the fish I want, I should've gotten at least one of the damn thing in 10 attempts, so what gives? 10% over 10 attempts is 100% right?"

And this made me think a lot about it, and the math that goes behind it. Intuitively, a 10% chance over 10 or so attempts should be a guarantee that you could encounter what you want at some point, right? Well, it wasn't the case for my friend, so that means this intuition is wrong. So what are the odds?

In statistics, you use probabilities. My professor in college describes them as "the chances of an event occurring". So in this example, the probability that you get the result you want is 10%, meaning P(success) = 10%. The probability that you fail (harsh word) is 90%, so P(failure) = 90%. The chances that you get at least one success is a combination of 10 different outcomes, with each success getting pushed slot to the right. The first outcome is getting a success on the first attempt. The odds of this are P(success) (10%). Since we don't care about the rest of the attempts, the odds of the rest of the attempts is 100%:

S _ _ _ _ _ _ _ _ _

The second outcome is if we fail once, then succeed. The odds of this are 90% x 10%:

U S _ _ _ _ _ _ _ _

The third outcome is if we fail twice, then succeed. The odds of this are 90% x 90% x 10%:

U U S _ _ _ _ _ _ _

Seeing a pattern here? With every additional outcome, the odds of getting a success multiply the previous unsuccessful attempts by P(failure). It's equivalent to the following:

P(at_least_one) = 10% + 90% x 10% + 90% x 90% x 10% + 90% x 90% x 90% x 10% + ...

Which generalizes to

P(at_least_one) = 10% x (1 + 90% ^ 1 + 90% ^ 2) => 10% x Sum (0 to 9) of (90% ^ i)

This form has a closed expression, leading to

P(at_least_one) = 1 - (1 - P(fail)) ** 10 => 1 - (1 - 0.9) ** 10 ≈ 65%.

There it is. The likelihood that you have at least 1 success out of 10 attempts with a 10% success rate is 65%. That leaves a 35% that you don't succeed. This explains why my friend wasn't able to get what he wanted.

One of the neat things about this formula is that it doesn't work with just 10 attempts - it can be generalized to n attempts. Thinking about this in different way, you could say the chance you get at least one of something over a given number of attempts, is the change you get none of that thing for all attempts subtracted by 1. Now I'm writing this, I really wish I thought of this approach earlier, so I didn't have to do this though experiment...

Okay, so that's cool, but what if we take this one step further? Let's say that I want to get at least a certain number of successes out of a given number of attempts. How do I do that?

The black sheep: bringing out the coding and the dynamic programming

To get the cat out of the bag, I'd like to approach this problem in a programming language, and you'll hopefully see why as I'm explaining it.

It is easy to get overwhelmed by this problem - "now there's two, three successes I have to get!?", it may appear there are way too many variables and the program is impossible to solve. What if I told you however, that it isn't impossible, at least if you look at the problem from the right perspective?

Step 1

We need to think about this problem, and the easiest way to do that is to come back to the original problem of getting at least 1 successful attempt:

S _ _ _ _ _ _ _ _ _

U S _ _ _ _ _ _ _ _

...

Notice that for every outcome, we now need to figure out the odds of getting another success with fewer attempts, we're basically solving the same problems with reduced parameters. For the mathematically inclined folks out there, this is kind of like a reverse inductive reasoning, the solution is comprised on a some operation between a constant factor and the problem with reduced parameters.

Let's say for example we want at least 3 successes out of 10 attempts, where the chance of getting a success is 10%. The odds of this are the odds we get 1 success in 10 attempts multiplied by the odds we get the remaining 2 successes in however few attempts we have left, which is determined by where get this first success. Since there are 10 attempts, we need to calculate the odds of getting 2 successes out of 9 attempts given we get a success on the first attempt plus the odds we get 2 successes out of 8 attempts given we get a success on the second attempt and so on and so forth...

I took the liberty of outlining the algorithm in python. From the way I described it, this algorithm will use some kind of recursion:

def calculate_chance_of_at_least_n_occurring(n: int, attempts: int, p_success: float, cache: dict[tuple[int, int], float] = {}) -> float:
    # Cache is a dict that stores the chance of at least n occurring at the specified n and attempts
    if n <= 0:
        return 1
    elif attempts < n: # Number of attempts less than requirement, return 0
        return 0
    elif attempts == n: # we have the same number of attempts as the requirements - so we need all the attempts to succeed.
        return p_success ** n
    elif (n, attempts) in cache: # return if we already computed this result
        return cache[(n, attempts)]
     
    total_prob: float = 0.0
    for i in range(attempts + 1 - n):
        change_of_at_least_n_minus_1_occurring: float = calculate_chance_of_at_least_n_occurring(n - 1, attempts - i - 1, p_success, cache)
        cache[(n - 1, attempts - i - 1)] = change_of_at_least_n_minus_1_occurring # Set the cache
        total_prob = total_prob + (1 - p_success) ** i * p_success * change_of_at_least_n_minus_1_occurring
    
    return total_prob

In this algorithm, I implemented some base cases, which are easy to evaluate:

  1. If we need 0 (or less), the chances are automatically 100%, we don't need anything, so each outcome can be anything.
  2. If we need more than what attempts we have, there is no chance we can do this, so we return 0
  3. If we have the same number of attempts as needs, then we simply calculate the odds of every attempt succeeding
  4. If we already computed the call with the provided arguments, we can return that computed result (more on this later)

We then go through every attempt and stop at a point where we have enough residual attempts to meet the needs demand. So for instance, we need 3 out of 10 attempts to succeed, we would only go 8 attempts because we still need 2 attempts for 2 needs in the worst case.

Now it's time to explain the cache, why do we need it? Well the issue with this algorithm is that we call the same function multiple times with the same arguments - and for heavy cases we call this function with the same arguments a lot. This is wasted work, so we can store the computed result of this function call into a cache that remembers every result, and returns these results as necessary.

This algorithm takes ~9.5 seconds when computing the odds of getting at least 99 successes out of 1000 attempts given an 85% chance of succeed every attempt

Optional: Python's Built-In Caching

This is language specific, but Python has a built-in way of managing the cache for this specific problem. Using their built-in cache is really easy, all we need to do is add a @cache decorator to our algorithm at the top and remove the cache parameter and edge case:

# Recursive solution with lru cache
from functools import cache

@cache
def calculate_chance_of_at_least_n_occuring(n: int, attempts: int, p_success: float) -> float:
    # Cache is a dict that stores the chance of at least n occuring at the specified n and attempts
    if n <= 0:
        return 1
    elif attempts < n: # Number of attempts less than requirement, return 0
        return 0
    elif attempts == n: # we have the same number of attempts as the requirements - so we need all the attempts to succeed.
        return p_success ** n
     
    total_prob: float = 0.0
    for i in range(attempts + 1 - n):
        change_of_at_least_n_minus_1_occuring: float = calculate_chance_of_at_least_n_occuring(n - 1, attempts - i - 1, p_success)
        total_prob = total_prob + (1 - p_success) ** i * p_success * change_of_at_least_n_minus_1_occuring
    
    return total_prob

This algorithm takes ~4.8 seconds when computing the odds of getting at least 99 successes out of 1000 attempts given an 85% chance of succeed every attempt. This surprised me, but I assume it's because the underlying caching implementation is much more efficient at caching parameters than storing them as tuple-keys in a hashmap!

Step 2

Now for those of you experienced in programming, you may realize this is a dynamic programming problem, and you'd be right, it most definitely is. We have a problem composed up a finite number of subproblems, many of which are the same subproblem. But like all dynamic programming problems, there exists a solution for this problem that doesn't have to deal with recursion!

Think about it this way - instead of coming from the solution, and calculating the sub problems, why not start from the sub problems, and work your way up to the solution? This is what the matrix approach does:

def calculate_chance_of_at_least_n_occurring(n: int, attempts: int, p_success: float) -> float:
    if attempts < n:
        return 0.0
    if attempts == n:
        return p_success ** attempts
    if n == 1: # Closed form of at least 1 success is known beforehand
        return 1 - (1 - p_success) ** attempts
    # Our probability matrix. Rows are the the requirements (I need at least n of this), columns are the attempts (I have x attempts for this)
    prob_mat: list[list[float]] = [[0 for _ in range(attempts + 1)] for _ in range(n + 1)] 
    prob_mat[0] = [1 for _ in range(attempts + 1)] # We will always guarantee to get at least 0 regardless of attempts
    prob_mat[1] = [(1 - (1 - p_success) ** attempt) if attempt > 0 else 0 for attempt in range(attempts + 1)] # We can solve for at least one attempt via the closed formula

    for requirement in range(2, len(prob_mat)):
        for current_attempts in range(requirement, attempts + 1): # all numbers of attempts below our requirement are 0 - we can't achieve those cases
            if requirement == current_attempts:
                prob_mat[requirement][current_attempts] = p_success ** current_attempts
            else:
                # Need to go through every spot where I get a success (followed by prior failures),  then multiply this probability with the computed previous result of n - 1
                # requirements, with current_attempts - 1 - attempt (- 1 to account for the success slot) attempts.
                for attempt in range(current_attempts):
                    prob_mat[requirement][current_attempts] += (1 - p_success) ** attempt * p_success * prob_mat[requirement - 1][current_attempts - attempt - 1] 


    return prob_mat[n][attempts]

In this solution, we no longer need a cache, because we use a matrix, which has a fast lookup time. The rows are the number of needs, and the columns are the attempts that we have. As an example row 2 and column 2 contains the odds that we get 2 successes with 2 attempts. We can fill out row 0 and row 1 for 0 and 1 needs respectfully (those cases should be trivial by now) and then calculate attempts 2 to n. We also assign every row column value where the row is greater than the column 0 since the number of attempts can't match the number of needs, and for every needs, we start on the attempts that are at least "needs" given by current_attempts in range(requirement, attempts + 1).

This algorithm takes ~3.4 seconds when computing the odds of getting at least 99 successes out of 1000 attempts given an 85% chance of succeed every attempt. Not a significant improvement over recursion, but that's to be expected considering the amount of space being used

Step 3

The above approach is optimal for time, but you can reduce the amount of space used if you look at what previous sub problems we are accessing. When you are trying to solve at least n successes given a attempts, you really only need to access the probability outcomes for n - 1 success given all possible attempts. This sounds fancy, but all it means is that we only need the previous row for a matrix and not all previous rows. The final solution is my magnum opus, and only needs 2 rows to produce the same number:

def calculate_chance_of_at_least_n_occurring(n: int, attempts: int, p_success: float) -> float:
    if attempts - n < 0:
        return 0.0
    if attempts - n == 0:
        return p_success ** attempts
    if attempts - n == 1: # Change of perspective, calculate the chance of at most one failure
        return ((p_success) ** (attempts - 1) * (1 - p_success) * attempts) + p_success ** attempts
    if n == 1: # Closed form of at least 1 success is known beforehand
        return 1 - (1 - p_success) ** attempts
    
    reading_row = [(1 - (1 - p_success) ** attempt) if attempt > 0 else 0 for attempt in range(attempts + 1)] # We can solve for at least one attempt via the closed formula
    writing_row = [0 for _ in range(attempts + 1)]

    for requirement in range(2, n + 1):
        for current_attempts in range(requirement, attempts + 1): # all numbers of attempts below our requirement are 0 - we can't achieve those cases
            if requirement == current_attempts:
                writing_row[current_attempts] = p_success ** current_attempts
            else:
                # Need to go through every spot where I get a success (followed by prior failures),  then multiply this probability with the computed previous result of n - 1 requirements, with current_attempts - 1 - attempt (- 1 to account for the success slot) attempts.
                for attempt in range(current_attempts):
                    writing_row[current_attempts] += (1 - p_success) ** attempt * p_success * reading_row[current_attempts - attempt - 1] 

        reading_row = writing_row
        writing_row = [0 for _ in range(attempts + 1)]

    return reading_row[-1]

In this solution, we added one more edge case. If the number of attempts outpaces our needs by 1, that means we can only have at most one failure. This can be interpreted as the chance we get exactly one failure + the chance we get all successes. This is what the line ((p_success) ** (attempts - 1) * (1 - p_success) * attempts) + p_success ** attempts calculates, since multiplication is commutative, all outcomes in this case can be represented the same way, multiplied by attempts. Then, we simply need to add the odds of every attempt being successful to get the odds of no failures.

This algorithm takes ~3.0 seconds when computing the odds of getting at least 99 successes out of 1000 attempts given an 85% chance of succeed every attempt. Again, not a significant improvement, but it is a tiny bit faster for much more efficient space usage. Instead of O(a * n) space usage for a attempts and n successes, you have a O(a) space usage. We can trim the space usage down by a constant factor, by doing fancy number magic and removing the 0 columns, reducing subsequent array sizes, but I don't think that would save much more space, and it would make the code much less readable.

Epilogue (Sort Of)

Try these algorithms out and measure their execution time! It was really fun to calculate these because not only does it give me more experience as a developer, it also teaches me that what's intuitive, is not always correct!


Back to Table of Contents