In [ ]:
import matplotlib.pyplot as plt
%matplotlib inline
import random
import math

Generating Random Numbers

Python has several functions built-in to the random library for generating random numbers. For example, the function random.uniform() may be useful for generating random floating point numbers:

In [ ]:
# random.uniform will generate a random number from a uniform range from low to high (including high!)
low = 0
high = 10
random_number = random.uniform(low, high)
print(random_number)

# run me!

You can also use random.randint() for random integers instead of floats:

In [ ]:
print(random.randint(low, high))

Activity 1) Game of craps

In a game of craps, you roll two dice and find their sum. Let's call this sum as sum_dice.

Here are the rules to find out if you win or lose:

1) If sum_dice is 7 or 11, you win

2) If sum_dice is 2, 3, or 12, you lose

3) If sum_dice is 4, 5, 6, 8, 9 or 10, you roll the dice again until you get either:

  • the new sum equal to sum_dice and you win

  • the new sum equal to 7 and you lose

Let's write a code snippet that simulates one game of craps:

Write the function sum_two_dice() which returns the sum of the roll of two dice

def sum_two_dice():
    # Write some code here
    return the_sum

Write the function check_craps_rule(sum_dice) that takes as argument the sum of the two dice, and returns 0 if the game corresponds to a loss and 1 if it corresponds to a win.

def check_craps_rule(sum_dice):
    # Write some code here
    # if win, won = 1
    # if loss, won = 0
    return won

Use your defined functions to simulate one game of craps:

check_craps_rule(sum_two_dice())

Did you win?

We can run the above cell many times and count by handhow many times we win to estimate the probability of winning. But that seems a little inefficient...

What is the probability of winning craps?

The exact value is

$$P(\text{win}) = \frac{244}{495} \approx 0.492929 $$

But we are not going to show you the proof here. Instead, we will perform many simulations of one game of crap (rolling two dice) and we will tally how many times we win (let us call this nwin). If the number of simulations N is large enough, we will get a good approximation for the probability of win as

$$P(\text{win}) = \frac{\text{nwin}}{N} $$

Write a code snippet that runs N=1000 simulations of a game of craps to get Pwin, the probability of win. Estimate also the relative error and save that in the variable error?

Repeat the simulations above, but instead use N=10000. What happens to your approximation? Compare your values of Pwin and error. Then try N=100000.

You just used Monte Carlo Methods to approximate the probability of winning of a game of craps. You will be learning a lot more about this in CS 357!

Activity 2) Can we help this person?

After having a "couple" of extra drinks, Mariastenchia was in urgent need to use the restroom. Luckly, she saw Murphy's pub open and decided to go in for a quick release. She is in a very delicate situation, and after a quick evaluation, she decided that if she cannot make it to the bathroom in less than 100 steps, she must be in serious trouble.

Do you think she can make a successful trip to the restroom? Let's help her estimating her odds.

Let's take a look at the floor plan of Murphy's pub

We will simplify the floor plan with a square domain of size room_size = 10. The bottom left corner of the room will be used as the origin, i.e. the position (x,y) = (0,0).

The bathroom location is indicated by a square, with the left bottom corner located at (8,8) and dimension h = 1. These quantities are stored as a tuple bathroom = (8,8,1).

The door is located at the bottom half, indicated by the blue line. Mariastenchia initial position is given by initial_position = (1,1), marked with a red half-dot.

In [ ]:
room_size = 10
bathroom = (8,8,1) 
initial_position = (1,1)

We will simplify Mariastenchia's challenge and remove all the tables and the bar area inside Murphy's. Here is how the room looks like in this simplified floor plan (you will be using this helper function later to vizualize Mariastenchia's path):

In [ ]:
# Helper function to draw the floor plan
# You should not modify anything here
def draw_murphy(wc,person,room_size):
    fig = plt.figure(figsize=(6,6))
    ax = fig.add_subplot(111)
    plt.xlim(0,room_size)
    plt.ylim(0,room_size)
    plt.xlabel("X Position")
    plt.ylabel("Y Position");
    ax.set_aspect('equal')
    
    rect = plt.Rectangle(wc[:2], wc[-1], wc[-1], color=(0.6,0.2,0.1) )
    ax.add_patch(rect)
    plt.text(wc[0],wc[1]+wc[2]+0.2,"WC")

    rect = plt.Rectangle((0,0),2,0.1, color=(0,0,1) )
    ax.add_patch(rect)
    plt.text(0.5,0.2,"door")
   
    circle = plt.Circle(person[:2], 0.3, color=(1,0,0))
    ax.add_patch(circle)
    
draw_murphy(bathroom,initial_position,room_size)

How are we going to model Mariastenchia's walk?

  • Since Mariastenchia had too many drinks, we will model her steps as a random walk. Generate a random angle $\alpha$ between $0$ and $2\pi$ to decide which direction she will move.
  • Each step has magnitude of 1. Once you have the angle, you can find the components of her step as:

$$ step = [\cos(\alpha), \sin(\alpha)]$$

Write the function random_step that takes as argument the current position, and returns the new position based on a random step with orientation $\alpha$.

def random_step(current_position):
    # current_position is a list: [x,y] 
    # x is the position along x-direction, y is the position along y-direction
    # new_position is also a list: [xnew,ynew]

    # Write some code here 
    return new_position

Let's make Mariastenchia give her 100 steps, using the function random_step. Complete the code snippet below, and plot her path from the door (given by the variable initial_position above) to her final location. Did she reach the bathroom?

Code snippet A

position = #create a list to store all the 100 positions

# Write code to simulate Mariastenchia 100 steps

draw_murphy(bathroom,initial_position,room_size)
x,y = zip(*position)
plt.plot(x,y)

You probably noticed Mariastenchia hitting walls, or even walking through them! Let's make sure we impose this constraints to the random walk, so we can get some feasible solution. Here are the rules for Mariastenchia's walk:

  • If Mariastenchia runs into a wall, the simulation stops, and we call it a "failure" (ops!)
  • If the sequence of steps takes Mariastenchia to the restroom, the simulation stops, and we call it a success (yay!).
  • To simplify the problem, the "restroom" does not have a door, and Mariastenchia can "enter" the square through any of its sides.
  • Mariastenchia needs to reach the restroom in less than 100 steps, otherwise the simulation stops (at this point, it is a little too late for her...). This is also called a failure.

Write the function check_rules that checks if the new_position is a valid one, according to the rules above. The function should return 0 if the new_position is a valid step (still inside Murphy's and searching for the restroom), 1 if new_position is inside the restroom (sucess), and -1 if new_position is a failure (Mariastenchia became a super woman and crossed a wall)

def check_rules(room_size,wc,current_position):
    # write some code here
    # return 0, 1 or -1
In [ ]:
# Try your function check_rules with the following conditions:
print(check_rules(room_size,bathroom,[-1,0]))
print(check_rules(room_size,bathroom,[5,5]))
print(check_rules(room_size,bathroom,[8.5,8.5]))

Modify your code snippet A above, so that for every step, you check for the constraints using check_rules. Instead of giving all the 100 steps, you should stop earlier in case she reaches the restroom (check_rules == 1) or she hits a wall ( check_rules == -1)

Code snippet B

In [ ]:

Let's estimate the probability Mariastenchia reaches the restroom:

You should now run many simulations of code snippet B (one attempt to reach the restroom), and tally how many attempts are successful. Start with N = 1000 simulations, and count number of success. The probability to reach the bathroom is n_success/N. Then you increase the number of simulations, try N = 10000 and N = 100000.

Code snippet C

n_success = 0
N = 1000
for i in range(N):
    position = # initial position
    # write some code here (use code snippet B)
print("probability is ", n_success/N)

It looks like this random walk does not give her much of a chance to get to the restroom. She may need more steps, or we can modify her walk to be a little less random. All you need to do is to modify your function random_step. What about we make her move forwards (in the positive direction of y) with a 70% probability (meaning she would move backwards with 30% probability). Then you just run code snippet C again.

Do you want to keep having fun with this? Add the bar in the middle of the room. Now your check_rules function needs to check for additional constraints. I think Mariastenchia should be able to take extra steps now that she has an extra challenge. Can she go around the bar and make it to the bathroom with "reasonable" probabilities if she can give 300 steps?