submit to programming reddit
 
(November 2010)

A nice puzzle... Read the code below...

#!/usr/bin/env python
import random

# A sickness with no symptoms has appeared in the general population:
# Chances of getting it, are one in ten thousand.
#
# Whoever has it, dies within a year.
#
# Naturally, people want to know if they are sick or not - to make sure
# their affairs are in order, fix their wills and testaments... say their
# farewells to their loved ones, and come to terms with their fate.
#
# Scientists have devised a test, which can show whether you have the 
# sickness or not, with an accuracy of 99%.
# 
# You just took the test - and you are told that it is positive.
#
# Should you start taking care of your post-mortem affairs?
# Or is it too soon?
#
# Don't read the rest, before you think about this.
# Be honest - no cheating!
#
#
#
# Well, let's assume that we have one million test subjects.
# One in ten thousand of them will be sick, which means that...
#
#     1.000.000 x 1 / 10.000 = 100 sick people
#
# ...100 of them will be sick - and 999.900 will be healthy.
# Now we perform the test on all one million of them: remember,
# the test is 99% accurate, so...
#
#                    1.000.000 people
#                     /           \
#                    /             \
#                   /               \
#              100 sick      999.900 healthy
#                 /\             /     \
#                /  \           /       \
#          99 pos   1 neg  9.999 pos  989.901 neg
#
# ...out of the 100 sick, it will give a positive result to 99,
# and miss only one, who will get a negative. Out of the 999.900
# healthy ones, 1% (that is, 9.999) will get a false positive,
# and the rest (999.900 - 9.999 = 989.901) will get a correct,
# negative result.
#
# So, grouping all who got positive results, you have:
#
# - 99 that are REALLY sick
# - and 9.999 that are healthy
#
# So, given the fact that the result is positive, your chance
# of ACTUALLY being sick, is:
#
#                                            99
#  P(I am sick if my test is positive) = ---------
#                                        99 + 9999
#
# ...which amounts to a little less than 1%.
# In other words, keep partying!
#
#
# The above is a simple example of the Bayes theorem, whose
# formal application is the following:
# ( the symbol '|' means 'given' )
#
#                          P(Sick)*P(TestPositive|Sick)
# P(Sick|TestPositive) = ------------------------------
#                                P(TestPositive)
#
#                     P(Sick)*P(TestPositive|Sick)
# =  -----------------------------------------------------------------
#    P(Sick)*P(TestPositive|Sick) + P(Healthy)*P(TestPositive|Healthy)
#
#               (1/10000)*0.99
# =  ------------------------------------- = 0.0098039215
#      (1/10000)*0.99 + (9999/10000)*0.01
#
# What follows below, is an experimental verification, 
# via a simple Python reproduction of the experiment.

# First, a helper function:
# When passed a number between 0 and 100, this function
# returns True with exactly the requested chance
# e.g. BoolProb(40) would return True with a 40% chance
#
def BoolProb(x):
    return random.random() < x/100.

# And now, the complete experiment:

# The experiment's sample size
SampleSize = 10000000
# How many positive results we got
PositiveTests = 0
# How many of the positive results were actually sick
PositiveTestAndSick = 0

for i in xrange(0, SampleSize):
    sick = BoolProb(0.01) # 1 in 10.000 => 0.01
    if sick:
	# If you are sick, you get a 99% chance of a positive
	testResultPositive = BoolProb(99) 
    else:
	# If you are healthy, you get a 1% chance of a positive
	testResultPositive = BoolProb(1)

    # Count the results
    if testResultPositive:
	PositiveTests += 1
        if sick:
	    PositiveTestAndSick += 1

print "Chance of being sick if test is positive:", \
    float(PositiveTestAndSick)/PositiveTests
Running the code will give you something like this:
bash$ python ./Am.I.sick.py
Chance of being sick if test is positive: 0.00955057847804
...that is, less than one per cent chance.

profile for ttsiodras at Stack Overflow, Q&A for professional and enthusiast programmers
GitHub member ttsiodras
Back to index  My CV  About meLast update on: Fri Aug 9 22:48:30 2013 (Valid HTML)

The comments on this website require the use of JavaScript. Perhaps your browser isn't JavaScript capable or the script is not being run for another reason. If you're interested in reading the comments or leaving a comment behind please try again with a different browser or from a different connection.