Experimental Mathematics

Sometimes you need to prove something for a project mathematically. However what if there might be an error in your reasoning that means that your proof is incorrect. That could lead to all sorts of problems when you rely on your incorrect proof.

One way to guard against this is to do some experimental mathematics to validate your proof.

 

Example

Lets look at an example. What is the expected number of rolls to roll a six on a six sided die. Well we know that the probability of rolling a six is independent of previous rolls so we can create a formula as follows

E = 1/6 * 1 + 5/6 * (E+1)

Where E is the expected number of rolls. This is because if we do not roll a six, the expectation after the first roll will be the same expectation E plus the one failed roll. Rearranging this we get:

E – 5/6 * E = 1/6 + 5/6

1/6 E = 1

E = 6

There are other proofs involving infinite series but personally I find this one easier to explain.

How can we check that what we have asserted here is sensible? Well lets run a quick experiment.

import random

def turns_to_roll_six():
    count = 0
    while True:
        roll = random.randint(1,6)
        count += 1
        if roll == 6:
            return count

valid_results = []
# run 10000 simulations
for turn in range(0, 10000):
    length = turns_to_roll_six()
    valid_results.append(length)

valid_paths = len(valid_results)
total_rolls = sum(valid_results)
print("number of valid paths: " + str(valid_paths))
print("total rolls in valid paths: " + str(total_rolls))
print("average expectation: " + str(total_rolls / valid_paths))

Running this code gives us a result like

This is in close agreement with our calculation of an expectation value of six. Thus it seems likely that our logic is sound, our experiment has supported our mathematical conclusion. You can test this for yourself by running the code on Trinket

Trickier Example

So applying the same logic we could ask for the expected number of rolls in the following situation posed by Gil Kalai. “You throw a dice until you get 6. What is the expected number of throws given that all throws gave even numbers”

Now at first glance it might seem that this is simple. No evens were thrown so the dice must behave like a 3 sided dice giving an expectation value of 3 rolls. In fact we have made a subtle mathematical error, but that is not necessarily apparent. We might think we have a valid proof.

It is at this point that doing an experiment can come to the rescue. Rather than just taking our “proof” and moving on lets write ourselves a second python script to do a Monte Carlo simulation of the expectation number for this situation.

import random

def turns_to_roll_six():
    count = 0
    while True:
        roll = random.randint(1,6)
        count += 1
        if roll == 6:
            return count
            # discount odds!
        if roll % 2 == 1:
            return 0

valid_results = []
# run 10000 simulations
for turn in range(0, 10000):
    length = turns_to_roll_six()
    # only count the paths that did not result in odds
    if length > 0 :
        valid_results.append(length)

valid_paths = len(valid_results)
total_rolls = sum(valid_results)
print("number of valid paths: " + str(valid_paths))
print("total rolls in valid paths: " + str(total_rolls)) 
print("average expectation: " + str(total_rolls / valid_paths))

Running this code tends to result in results in the range 1.4 to 1.6!

Once again I have put a copy on Trinket so that you can try this for yourself. Clearly this is a fair way away from our predicted result of 3 so we know that we have made a mistake in our proof.

What went wrong?

The error lies in the fact that odd rolls are removed from our sample pool in such a way that they remove entire expectation chains. This means that this system does not behave like a 3 sided die. We have made an incorrect assumption. In fact, 4/6 of the time (6s and odds) the chain terminates. Only in 2/6 cases does it proceed to a subsequent roll. While we eventually discount all the chains that terminated due to odds being rolled (3/4 of all chains), they still affect the distribution. This means that the equation we are looking for is as follows

E = 4/6 * 1  + 2/6 * (E + 1)

Which resolves to:

E – 2/6 * E = 4/6 + 2/6

4/6 E = 1

E = 6/4

E = 1.5

Therefore the result we are actually looking for is 1.5. This is is in accordance with what we found from our experiment.

Conclusion

I hope I have convinced you that performing experiments can be useful when attempting to validate your mathematical reasoning. Obviously you can’t guarantee to catch a proof error this way. However by covering any edge or corner cases you can think of with practical tests you can reduce the likelihood of error. I have found this approach especially useful when reasoning about probability. It helps to avoid subtle errors in reasoning that will be hard to track down when embedded in larger systems.

Leave a Reply

Your email address will not be published. Required fields are marked *

This site uses Akismet to reduce spam. Learn how your comment data is processed.