In [1]:
# Configure Jupyter so figures appear in the notebook
%matplotlib inline
# Configure Jupyter to display the assigned value after an assignment
%config InteractiveShell.ast_node_interactivity='last_expr_or_assign'
# import functions from the modsim library
from modsim import *
# set the random number generator
np.random.seed(7)
# If this cell runs successfully, it produces no output.
We'll start with a State object that represents the number of bikes at each station.
When you display a State object, it lists the state variables and their values:
In [2]:
bikeshare = State(olin=10, wellesley=2)
Out[2]:
We can access the state variables using dot notation.
In [3]:
bikeshare.olin
Out[3]:
In [4]:
bikeshare.wellesley
Out[4]:
Exercise: What happens if you spell the name of a state variable wrong? Edit the previous cell, change the spelling of wellesley, and run the cell again.
The error message uses the word "attribute", which is another name for what we are calling a state variable.
Exercise: Add a third attribute called babson with initial value 0, and display the state of bikeshare again.
In [5]:
bikeshare.olin -= 1
If we display bikeshare, we should see the change.
In [6]:
bikeshare
Out[6]:
Of course, if we subtract a bike from olin, we should add it to wellesley.
In [7]:
bikeshare.wellesley += 1
bikeshare
Out[7]:
In [8]:
def bike_to_wellesley():
bikeshare.olin -= 1
bikeshare.wellesley += 1
When you define a function, it doesn't run the statements inside the function, yet. When you call the function, it runs the statements inside.
In [9]:
bike_to_wellesley()
bikeshare
Out[9]:
One common error is to omit the parentheses, which has the effect of looking up the function, but not calling it.
In [10]:
bike_to_wellesley
Out[10]:
The output indicates that bike_to_wellesley is a function defined in a "namespace" called __main__, but you don't have to understand what that means.
Exercise: Define a function called bike_to_olin that moves a bike from Wellesley to Olin. Call the new function and display bikeshare to confirm that it works.
In [11]:
# Solution
def bike_to_olin():
bikeshare.wellesley -= 1
bikeshare.olin += 1
In [12]:
# Solution
bike_to_olin()
bikeshare
Out[12]:
modsim.py provides flip, which takes a probability and returns either True or False, which are special values defined by Python.
The Python function help looks up a function and displays its documentation.
In [13]:
help(flip)
In the following example, the probability is 0.7 or 70%. If you run this cell several times, you should get True about 70% of the time and False about 30%.
In [14]:
flip(0.7)
Out[14]:
In the following example, we use flip as part of an if statement. If the result from flip is True, we print heads; otherwise we do nothing.
In [15]:
if flip(0.7):
print('heads')
With an else clause, we can print heads or tails depending on whether flip returns True or False.
In [16]:
if flip(0.7):
print('heads')
else:
print('tails')
In [17]:
bikeshare = State(olin=10, wellesley=2)
Out[17]:
Suppose that in any given minute, there is a 50% chance that a student picks up a bike at Olin and rides to Wellesley. We can simulate that like this.
In [18]:
if flip(0.5):
bike_to_wellesley()
print('Moving a bike to Wellesley')
bikeshare
Out[18]:
And maybe at the same time, there is also a 40% chance that a student at Wellesley rides to Olin.
In [19]:
if flip(0.4):
bike_to_olin()
print('Moving a bike to Olin')
bikeshare
Out[19]:
We can wrap that code in a function called step that simulates one time step. In any given minute, a student might ride from Olin to Wellesley, from Wellesley to Olin, or both, or neither, depending on the results of flip.
In [20]:
def step():
if flip(0.5):
bike_to_wellesley()
print('Moving a bike to Wellesley')
if flip(0.4):
bike_to_olin()
print('Moving a bike to Olin')
Since this function takes no parameters, we call it like this:
In [21]:
step()
bikeshare
Out[21]:
In [22]:
def step(p1, p2):
if flip(p1):
bike_to_wellesley()
print('Moving a bike to Wellesley')
if flip(p2):
bike_to_olin()
print('Moving a bike to Olin')
Now we can call it like this:
In [23]:
step(0.5, 0.4)
bikeshare
Out[23]:
Exercise: At the beginning of step, add a print statement that displays the values of p1 and p2. Call it again with values 0.3, and 0.2, and confirm that the values of the parameters are what you expect.
In [24]:
# Solution
def step(p1, p2):
print(p1, p2)
if flip(p1):
bike_to_wellesley()
print('Moving a bike to Wellesley')
if flip(p2):
bike_to_olin()
print('Moving a bike to Olin')
step(0.3, 0.2)
Before we go on, I'll redefine step without the print statements.
In [25]:
def step(p1, p2):
if flip(p1):
bike_to_wellesley()
if flip(p2):
bike_to_olin()
And let's start again with a new State object:
In [26]:
bikeshare = State(olin=10, wellesley=2)
Out[26]:
We can use a for loop to move 4 bikes from Olin to Wellesley.
In [27]:
for i in range(4):
bike_to_wellesley()
bikeshare
Out[27]:
Or we can simulate 4 random time steps.
In [28]:
for i in range(4):
step(0.3, 0.2)
bikeshare
Out[28]:
If each step corresponds to a minute, we can simulate an entire hour like this.
In [29]:
for i in range(60):
step(0.3, 0.2)
bikeshare
Out[29]:
After 60 minutes, you might see that the number of bike at Olin is negative. We'll fix that problem in the next notebook.
But first, we want to plot the results.
In [30]:
results = TimeSeries()
Out[30]:
And we can add a value to the TimeSeries like this:
In [31]:
results[0] = bikeshare.olin
results
Out[31]:
The 0 in brackets is an index that indicates that this value is associated with time step 0.
Now we'll use a for loop to save the results of the simulation. I'll start one more time with a new State object.
In [32]:
bikeshare = State(olin=10, wellesley=2)
Out[32]:
Here's a for loop that runs 10 steps and stores the results.
In [33]:
for i in range(10):
step(0.3, 0.2)
results[i] = bikeshare.olin
Now we can display the results.
In [34]:
results
Out[34]:
A TimeSeries is a specialized version of a Pandas Series, so we can use any of the functions provided by Series, including several that compute summary statistics:
In [35]:
results.mean()
Out[35]:
In [36]:
results.describe()
Out[36]:
You can read the documentation of Series here.
In [37]:
plot(results, label='Olin')
decorate(title='Olin-Wellesley Bikeshare',
xlabel='Time step (min)',
ylabel='Number of bikes')
savefig('figs/chap02-fig01.pdf')
decorate, which is defined in the modsim library, adds a title and labels the axes.
In [38]:
help(decorate)
savefig() saves a figure in a file.
In [39]:
help(savefig)
The suffix of the filename indicates the format you want. This example saves the current figure in a PDF file named chap01-fig01.pdf.
Exercise: Wrap the code from this section in a function named run_simulation that takes three parameters, named p1, p2, and num_steps.
It should:
TimeSeries object to hold the results.step the number of times specified by num_steps, passing along the specified values of p1 and p2.TimeSeries.To test your function:
State object with the initial state of the system.run_simulation with appropriate parameters.Optional:
TimeSeries objects, keeps track of the number of bikes at Olin and at Wellesley, and plots both series at the end.
In [40]:
# Solution
def run_simulation(p1, p2, num_steps):
olin = TimeSeries()
wellesley = TimeSeries()
for i in range(num_steps):
step(p1, p2)
olin[i] = bikeshare.olin
wellesley[i] = bikeshare.wellesley
plot(olin, label='Olin')
plot(wellesley, label='Wellesley')
decorate(title='Olin-Wellesley Bikeshare',
xlabel='Time step (min)',
ylabel='Number of bikes')
In [41]:
# Solution
bikeshare = State(olin=10, wellesley=2)
run_simulation(0.3, 0.2, 60)
The functions in modsim.py are built on top of several widely-used Python libraries, especially NumPy, SciPy, and Pandas. These libraries are powerful but can be hard to use. The intent of modsim.py is to give you the power of these libraries while making it easy to get started.
In the future, you might want to use these libraries directly, rather than using modsim.py. So we will pause occasionally to open the hood and let you see how modsim.py works.
You don't need to know anything in these sections, so if you are already feeling overwhelmed, you might want to skip them. But if you are curious, read on.
This chapter introduces two objects, State and TimeSeries. Both are based on the Series object defined by Pandas, which is a library primarily used for data science.
You can read the documentation of the Series object here
The primary differences between TimeSeries and Series are:
I made it easier to create a new, empty Series while avoiding a confusing inconsistency.
I provide a function so the Series looks good when displayed in Jupyter.
I provide a function called set that we'll use later.
State has all of those capabilities; in addition, it provides an easier way to initialize state variables, and it provides functions called T and dt, which will help us avoid a confusing error later.
The plot function in modsim.py is based on the plot function in Pyplot, which is part of Matplotlib. You can read the documentation of plot here.
decorate provides a convenient way to call the pyplot functions title, xlabel, and ylabel, and legend. It also avoids an annoying warning message if you try to make a legend when you don't have any labelled lines.
In [42]:
help(decorate)
In [43]:
source_code(flip)
In [ ]: