Introductory Tutorial#

The Boltzmann Wealth Model#

Important:

  • If you are just exploring Mesa and want the fastest way to execute the code we recommend executing this tutorial online in a Colab notebook. Colab

  • If you have installed mesa and are running locally, please ensure that your Mesa version is up-to-date in order to run this tutorial.

Tutorial Description#

Mesa is a Python framework for agent-based modeling. This tutorial will assist you in getting started. Working through the tutorial will help you discover the core features of Mesa. Through the tutorial, you are walked through creating a starter-level model. Functionality is added progressively as the process unfolds. Should anyone find any errors, bugs, have a suggestion, or just are looking for clarification, let us know!

The premise of this tutorial is to create a starter-level model representing agents exchanging money. This exchange of money affects wealth.

Next, space is added to allow agents to move based on the change in wealth as time progresses.

Two of Mesa’s analytic tools: the data collector and batch runner are then used to examine the dynamics of this simple model.

More Tutorials:#

Visualization: There is a separate visualization tutorial that will take users through building a visualization for this model (aka Boltzmann Wealth Model).

Advanced Visualization (legacy): There is also an advanced visualization tutorial that will show users how to use the JavaScript based visualization option, which also uses this model as its base.

Model Description#

This is a starter-level simulated agent-based economy. In an agent-based economy, the behavior of an individual economic agent, such as a consumer or producer, is studied in a market environment. This model is drawn from the field econophysics, specifically a paper prepared by Drăgulescu et al. for additional information on the modeling assumptions used in this model. [Drăgulescu, 2002].

The assumption that govern this model are:

  1. There are some number of agents.

  2. All agents begin with 1 unit of money.

  3. At every step of the model, an agent gives 1 unit of money (if they have it) to some other agent.

Even as a starter-level model the yielded results are both interesting and unexpected to individuals unfamiliar with it the specific topic. As such, this model is a good starting point to examine Mesa’s core features.

Tutorial Setup#

Create and activate a virtual environment. Python version 3.9 or higher is required.

Install Mesa:

pip install --upgrade mesa

Install Jupyter Notebook (optional):

pip install jupyter

Install Seaborn (which is used for data visualization):

pip install seaborn

If running in Google Colab run the below cell to install Mesa. (This will also work in a locally installed version of Jupyter.)

# SKIP THIS CELL unless running in colab

%pip install --quiet mesa
# The exclamation points tell jupyter to do the command via the command line
Note: you may need to restart the kernel to use updated packages.

Building the Sample Model#

After Mesa is installed a model can be built. A jupyter notebook is recommended for this tutorial, this allows for small segments of codes to be examined one at a time. As an option this can be created using python script files.

Good Practice: Place a model in its own folder/directory. This is not specifically required for the starter_model, but as other models become more complicated and expand multiple python scripts, documentation, discussions and notebooks may be added.

Create New Folder/Directory#

  • Using operating system commands create a new folder/directory named ‘starter_model’.

  • Change into the new folder/directory.

Creating Model With Jupyter Notebook#

Write the model interactively in Jupyter Notebook cells.

Start Jupyter Notebook:

jupyter notebook

Create a new Notebook named money_model.ipynb (or whatever you want to call it).

Creating Model With Script File (IDE, Text Editor, Colab, etc.)#

Create a new file called money_model.py (or whatever you want to call it)

Code will be added as the tutorial progresses.

Import Dependencies#

This includes importing of dependencies needed for the tutorial.

import mesa

# Data visualization tools.
import seaborn as sns

# Has multi-dimensional arrays and matrices. Has a large collection of
# mathematical functions to operate on these arrays.
import numpy as np

# Data manipulation and analysis.
import pandas as pd

Create Agent#

First create the agent. As the tutorial progresses, more functionality will be added to the agent.

Background: Agents are the individual entities that act in the model. It is a good modeling practice to make certain each Agent can be uniquely identified.

Model-specific information: Agents are the individuals that exchange money, in this case the amount of money an individual agent has is represented as wealth. Additionally, agents each have a unique identifier.

Code implementation: This is done by creating a new class (or object) that extends mesa.Agent creating a subclass of the Agent class from mesa. The new class is named MoneyAgent. The technical details about the agent object can be found in the mesa repo.

The MoneyAgent class is created with the following code:

class MoneyAgent(mesa.Agent):
    """An agent with fixed initial wealth."""

    def __init__(self, unique_id, model):
        # Pass the parameters to the parent class.
        super().__init__(unique_id, model)

        # Create the agent's variable and set the initial values.
        self.wealth = 1

Create Model#

Next, create the model. Again, as the tutorial progresses, more functionality will be added to the model.

Background: The model can be visualized as a grid containing all the agents. The model creates, holds and manages all the agents on the grid. The model evolves in discrete time steps.

Model-specific information: When a model is created the number of agents within the model is specified. The model then creates the agents and places them on the grid. The model also contains a scheduler which controls the order in which agents are activated. The scheduler is also responsible for advancing the model by one step. The model also contains a data collector which collects data from the model. These topics will be covered in more detail later in the tutorial.

Code implementation: This is done by creating a new class (or object) that extends mesa.Model and calls super().__init__(), creating a subclass of the Model class from mesa. The new class is named MoneyModel. The technical details about the agent object can be found in the mesa repo.

The MoneyModel class is created with the following code:

class MoneyModel(mesa.Model):
    """A model with some number of agents."""

    def __init__(self, N):
        super().__init__()
        self.num_agents = N
        # Create agents
        for i in range(self.num_agents):
            a = MoneyAgent(i, self)

Adding the Scheduler#

Now the model will be modified to add a scheduler.

Background: The scheduler controls the order in which agents are activated, causing the agent to take their defined action. The scheduler is also responsible for advancing the model by one step. A step is the smallest unit of time in the model, and is often referred to as a tick. The scheduler can be configured to activate agents in different orders. This can be important as the order in which agents are activated can impact the results of the model [Comer2014]. At each step of the model, one or more of the agents – usually all of them – are activated and take their own step, changing internally and/or interacting with one another or the environment.

Model-specific information: A new class is named RandomActivationByAgent is created which extends mesa.time.RandomActivation creating a subclass of the RandomActivation class from Mesa. This class activates all the agents once per step, in random order. Every agent is expected to have a step method. The step method is the action the agent takes when it is activated by the model schedule. We add an agent to the schedule using the add method; when we call the schedule’s step method, the model shuffles the order of the agents, then activates and executes each agent’s step method. The scheduler is then added to the model.

Code implementation: The technical details about the timer object can be found in the mesa repo. Mesa offers a few different built-in scheduler classes, with a common interface. That makes it easy to change the activation regime a given model uses, and see whether it changes the model behavior. The details pertaining to the scheduler interface can be located the same mesa repo.

With that in mind, the MoneyAgent code is modified below to visually show when a new agent is created. The MoneyModel code is modified by adding the RandomActivation method to the model. with the scheduler added looks like this:

class MoneyAgent(mesa.Agent):
    """An agent with fixed initial wealth."""

    def __init__(self, unique_id, model):
        # Pass the parameters to the parent class.
        super().__init__(unique_id, model)

        # Create the agent's attribute and set the initial values.
        self.wealth = 1

    def step(self):
        # The agent's step will go here.
        # For demonstration purposes we will print the agent's unique_id
        print(f"Hi, I am an agent, you can call me {str(self.unique_id)}.")


class MoneyModel(mesa.Model):
    """A model with some number of agents."""

    def __init__(self, N):
        super().__init__()
        self.num_agents = N
        # Create scheduler and assign it to the model
        self.schedule = mesa.time.RandomActivation(self)

        # Create agents
        for i in range(self.num_agents):
            a = MoneyAgent(i, self)
            # Add the agent to the scheduler
            self.schedule.add(a)

    def step(self):
        """Advance the model by one step."""

        # The model's step will go here for now this will call the step method of each agent and print the agent's unique_id
        self.schedule.step()

Running the Model#

A basic model has now been created. The model can be run by creating a model object and calling the step method. The model will run for one step and print the unique_id of each agent. You may run the model for multiple steps by calling the step method multiple times.

Note: If you are using .py (script) files instead of .ipynb (Jupyter), the common convention is to have a run.py in the same directory as your model code. You then (1) import the MoneyModel class, (2) create a model object and (3) run it for a few steps. As shown below:

from money_model import MoneyModel

starter_model = MoneyModel(10)
starter_model.step()

Create the model object, and run it for one step:

starter_model = MoneyModel(10)
starter_model.step()
Hi, I am an agent, you can call me 9.
Hi, I am an agent, you can call me 6.
Hi, I am an agent, you can call me 4.
Hi, I am an agent, you can call me 5.
Hi, I am an agent, you can call me 3.
Hi, I am an agent, you can call me 8.
Hi, I am an agent, you can call me 1.
Hi, I am an agent, you can call me 7.
Hi, I am an agent, you can call me 2.
Hi, I am an agent, you can call me 0.
# Run this step overnight and see what happens! Notice the order of the agents changes each time.
starter_model.step()
Hi, I am an agent, you can call me 3.
Hi, I am an agent, you can call me 9.
Hi, I am an agent, you can call me 0.
Hi, I am an agent, you can call me 5.
Hi, I am an agent, you can call me 7.
Hi, I am an agent, you can call me 8.
Hi, I am an agent, you can call me 2.
Hi, I am an agent, you can call me 6.
Hi, I am an agent, you can call me 4.
Hi, I am an agent, you can call me 1.

Exercise#

Modifying the code below to have every agent print out its wealth when it is activated.

class MoneyAgent(mesa.Agent):
    """An agent with fixed initial wealth."""

    def __init__(self, unique_id, model):
        # Pass the parameters to the parent class.
        super().__init__(unique_id, model)

        # Create the agent's variable and set the initial values.
        self.wealth = 1

    def step(self):
        # The agent's step will go here.
        # FIXME: Need to print the agent's wealth
        print(f"Hi, I am an agent and I am broke!")

Create a model for 12 Agents, and run it for a few steps to see the output.

# Fixme: Create the model object, and run it

Agent Step#

Returning back to the MoneyAgent the actual step process is now going to be created.

Background: This is where the agent’s behavior as it relates to each step or tick of the model is defined.

Model-specific information: In this case, the agent will check its wealth, and if it has money, give one unit of it away to another random agent.

Code implementation: The agent’s step method is called by the scheduler during each step of the model. To allow the agent to choose another agent at random, we use the model.random random-number generator. This works just like Python’s random module, but with a fixed seed set when the model is instantiated, that can be used to replicate a specific model run later.

To pick an agent at random, we need a list of all agents. Notice that there isn’t such a list explicitly in the model. The scheduler, however, does have an internal list of all the agents it is scheduled to activate.

With that in mind, we rewrite the agent step method as shown below:

import copy


class MoneyAgent(mesa.Agent):
    """An agent with fixed initial wealth."""

    def __init__(self, unique_id, model):
        # Pass the parameters to the parent class.
        super().__init__(unique_id, model)

        # Create the agent's variable and set the initial values.
        self.wealth = 1

    def step(self):
        # Verify agent has some wealth
        if self.wealth > 0:
            other_agent = self.random.choice(self.model.schedule.agents)
            if other_agent is not None:
                other_agent.wealth += 1
                self.wealth -= 1

Running your first model#

With that last piece in hand, it’s time for the first rudimentary run of the model.

If you’ve written the code in its own script file (money_model.py or a different name) you can now modify your run.py or even launch a Jupyter Notebook. You then just follow the same three steps of (1) import your model class MoneyModel, (2) create the model object and (3) run it for a few steps. If you wrote the code in one Notebook then step 1, importing, is not necessary.

from money_model import MoneyModel

Now let’s create a model with 10 agents, and run it for 10 steps.

model = MoneyModel(10)
for i in range(10):
    model.step()

Next, we need to get some data out of the model. Specifically, we want to see the distribution of the agent’s wealth. We can get the wealth values with list comprehension, and then use seaborn (or another graphics library) to visualize the data in a histogram.

agent_wealth = [a.wealth for a in model.schedule.agents]
# Create a histogram with seaborn
g = sns.histplot(agent_wealth, discrete=True)
g.set(
    title="Wealth distribution", xlabel="Wealth", ylabel="Number of agents"
);  # The semicolon is just to avoid printing the object representation
../_images/a44b9e445ea1507255a227291946a74b81a46ef6b9fa9b405479692af16d0ca8.png

You’ll should see something like the distribution above. Yours will almost certainly look at least slightly different, since each run of the model is random, after all.

To get a better idea of how a model behaves, we can create multiple model runs and see the distribution that emerges from all of them. We can do this with a nested for loop:

all_wealth = []
# This runs the model 100 times, each model executing 10 steps.
for j in range(100):
    # Run the model
    model = MoneyModel(10)
    for i in range(10):
        model.step()

    # Store the results
    for agent in model.schedule.agents:
        all_wealth.append(agent.wealth)

# Use seaborn
g = sns.histplot(all_wealth, discrete=True)
g.set(title="Wealth distribution", xlabel="Wealth", ylabel="Number of agents");
../_images/1df7a6091510d9e4b2a6c52c89f5da26eccc71d414f4540dbd3486f8c2a59fbf.png

This runs 100 instantiations of the model, and runs each for 10 steps. (Notice that we set the histogram bins to be integers, since agents can only have whole numbers of wealth). This distribution looks a lot smoother. By running the model 100 times, we smooth out some of the ‘noise’ of randomness, and get to the model’s overall expected behavior.

This outcome might be surprising. Despite the fact that all agents, on average, give and receive one unit of money every step, the model converges to a state where most agents have a small amount of money and a small number have a lot of money.

Adding space#

Many ABMs have a spatial element, with agents moving around and interacting with nearby neighbors. Mesa currently supports two overall kinds of spaces: grid, and continuous. Grids are divided into cells, and agents can only be on a particular cell, like pieces on a chess board. Continuous space, in contrast, allows agents to have any arbitrary position. Both grids and continuous spaces are frequently toroidal, meaning that the edges wrap around, with cells on the right edge connected to those on the left edge, and the top to the bottom. This prevents some cells having fewer neighbors than others, or agents being able to go off the edge of the environment.

Let’s add a simple spatial element to our model by putting our agents on a grid and make them walk around at random. Instead of giving their unit of money to any random agent, they’ll give it to an agent on the same cell.

Mesa has two main types of grids: SingleGrid and MultiGrid*. SingleGrid enforces at most one agent per cell; MultiGrid allows multiple agents to be in the same cell. Since we want agents to be able to share a cell, we use MultiGrid.

*However there are more types of space to include HexGrid, NetworkGrid, and the previously mentioned ContinuousSpace. Similar to mesa.time context is retained with mesa.space.[enter class]. You can see the different classes as mesa.space

We instantiate a grid with width and height parameters, and a boolean as to whether the grid is toroidal. Let’s make width and height model parameters, in addition to the number of agents, and have the grid always be toroidal. We can place agents on a grid with the grid’s place_agent method, which takes an agent and an (x, y) tuple of the coordinates to place the agent.

class MoneyModel(mesa.Model):
    """A model with some number of agents."""

    def __init__(self, N, width, height):
        super().__init__()
        self.num_agents = N
        self.grid = mesa.space.MultiGrid(width, height, True)
        self.schedule = mesa.time.RandomActivation(self)

        # Create agents
        for i in range(self.num_agents):
            a = MoneyAgent(i, self)
            self.schedule.add(a)

            # Add the agent to a random grid cell
            x = self.random.randrange(self.grid.width)
            y = self.random.randrange(self.grid.height)
            self.grid.place_agent(a, (x, y))

Under the hood, each agent’s position is stored in two ways: the agent is contained in the grid in the cell it is currently in, and the agent has a pos variable with an (x, y) coordinate tuple. The place_agent method adds the coordinate to the agent automatically.

Now we need to add to the agents’ behaviors, letting them move around and only give money to other agents in the same cell.

First let’s handle movement, and have the agents move to a neighboring cell. The grid object provides a move_agent method, which like you’d imagine, moves an agent to a given cell. That still leaves us to get the possible neighboring cells to move to. There are a couple ways to do this. One is to use the current coordinates, and loop over all coordinates +/- 1 away from it. For example:

neighbors = []
x, y = self.pos
for dx in [-1, 0, 1]:
    for dy in [-1, 0, 1]:
        neighbors.append((x+dx, y+dy))

But there’s an even simpler way, using the grid’s built-in get_neighborhood method, which returns all the neighbors of a given cell. This method can get two types of cell neighborhoods: Moore (includes all 8 surrounding squares), and Von Neumann(only up/down/left/right). It also needs an argument as to whether to include the center cell itself as one of the neighbors.

With that in mind, the agent’s move method looks like this:

class MoneyAgent(mesa.Agent):
   #...
    def move(self):
        possible_steps = self.model.grid.get_neighborhood(
            self.pos,
            moore=True,
            include_center=False)
        new_position = self.random.choice(possible_steps)
        self.model.grid.move_agent(self, new_position)

Next, we need to get all the other agents present in a cell, and give one of them some money. We can get the contents of one or more cells using the grid’s get_cell_list_contents method, or by accessing a cell directly. The method accepts a list of cell coordinate tuples, or a single tuple if we only care about one cell.

class MoneyAgent(mesa.Agent):
    #...
    def give_money(self):
        cellmates = self.model.grid.get_cell_list_contents([self.pos])
        if len(cellmates) > 1:
            other = self.random.choice(cellmates)
            other.wealth += 1
            self.wealth -= 1

And with those two methods, the agent’s step method becomes:

class MoneyAgent(mesa.Agent):
    # ...
    def step(self):
        self.move()
        if self.wealth > 0:
            self.give_money()

Now, putting that all together should look like this:

class MoneyAgent(mesa.Agent):
    """An agent with fixed initial wealth."""

    def __init__(self, unique_id, model):
        super().__init__(unique_id, model)
        self.wealth = 1

    def move(self):
        possible_steps = self.model.grid.get_neighborhood(
            self.pos, moore=True, include_center=False
        )
        new_position = self.random.choice(possible_steps)
        self.model.grid.move_agent(self, new_position)

    def give_money(self):
        cellmates = self.model.grid.get_cell_list_contents([self.pos])
        if len(cellmates) > 1:
            other_agent = self.random.choice(cellmates)
            other_agent.wealth += 1
            self.wealth -= 1

    def step(self):
        self.move()
        if self.wealth > 0:
            self.give_money()


class MoneyModel(mesa.Model):
    """A model with some number of agents."""

    def __init__(self, N, width, height):
        super().__init__()
        self.num_agents = N
        self.grid = mesa.space.MultiGrid(width, height, True)
        self.schedule = mesa.time.RandomActivation(self)
        # Create agents
        for i in range(self.num_agents):
            a = MoneyAgent(i, self)
            self.schedule.add(a)
            # Add the agent to a random grid cell
            x = self.random.randrange(self.grid.width)
            y = self.random.randrange(self.grid.height)
            self.grid.place_agent(a, (x, y))

    def step(self):
        self.schedule.step()

Let’s create a model with 100 agents on a 10x10 grid, and run it for 20 steps.

model = MoneyModel(100, 10, 10)
for i in range(20):
    model.step()

Now let’s use seaborn and numpy to visualize the number of agents residing in each cell. To do that, we create a numpy array of the same size as the grid, filled with zeros. Then we use the grid object’s coord_iter() feature, which lets us loop over every cell in the grid, giving us each cell’s positions and contents in turn.

agent_counts = np.zeros((model.grid.width, model.grid.height))
for cell_content, (x, y) in model.grid.coord_iter():
    agent_count = len(cell_content)
    agent_counts[x][y] = agent_count
# Plot using seaborn, with a size of 5x5
g = sns.heatmap(agent_counts, cmap="viridis", annot=True, cbar=False, square=True)
g.figure.set_size_inches(4, 4)
g.set(title="Number of agents on each cell of the grid");
../_images/737a3adaaa632f41d5aba191d0942d5934ffe8f09c6bcdb097fb524915793797.png

Collecting Data#

So far, at the end of every model run, we’ve had to go and write our own code to get the data out of the model. This has two problems: it isn’t very efficient, and it only gives us end results. If we wanted to know the wealth of each agent at each step, we’d have to add that to the loop of executing steps, and figure out some way to store the data.

Since one of the main goals of agent-based modeling is generating data for analysis, Mesa provides a class which can handle data collection and storage for us and make it easier to analyze.

The data collector stores three categories of data: model-level variables, agent-level variables, and tables (which are a catch-all for everything else). Model- and agent-level variables are added to the data collector along with a function for collecting them. Model-level collection functions take a model object as an input, while agent-level collection functions take an agent object as an input. Both then return a value computed from the model or each agent at their current state. When the data collector’s collect method is called, with a model object as its argument, it applies each model-level collection function to the model, and stores the results in a dictionary, associating the current value with the current step of the model. Similarly, the method applies each agent-level collection function to each agent currently in the schedule, associating the resulting value with the step of the model, and the agent’s unique_id.

Let’s add a DataCollector to the model with mesa.DataCollector, and collect two variables. At the agent level, we want to collect every agent’s wealth at every step. At the model level, let’s measure the model’s Gini Coefficient, a measure of wealth inequality.

def compute_gini(model):
    agent_wealths = [agent.wealth for agent in model.schedule.agents]
    x = sorted(agent_wealths)
    N = model.num_agents
    B = sum(xi * (N - i) for i, xi in enumerate(x)) / (N * sum(x))
    return 1 + (1 / N) - 2 * B


class MoneyAgent(mesa.Agent):
    """An agent with fixed initial wealth."""

    def __init__(self, unique_id, model):
        super().__init__(unique_id, model)
        self.wealth = 1

    def move(self):
        possible_steps = self.model.grid.get_neighborhood(
            self.pos, moore=True, include_center=False
        )
        new_position = self.random.choice(possible_steps)
        self.model.grid.move_agent(self, new_position)

    def give_money(self):
        cellmates = self.model.grid.get_cell_list_contents([self.pos])
        cellmates.pop(
            cellmates.index(self)
        )  # Ensure agent is not giving money to itself
        if len(cellmates) > 1:
            other = self.random.choice(cellmates)
            other.wealth += 1
            self.wealth -= 1
            if other == self:
                print("I JUST GAVE MONEY TO MYSELF HEHEHE!")

    def step(self):
        self.move()
        if self.wealth > 0:
            self.give_money()


class MoneyModel(mesa.Model):
    """A model with some number of agents."""

    def __init__(self, N, width, height):
        super().__init__()
        self.num_agents = N
        self.grid = mesa.space.MultiGrid(width, height, True)
        self.schedule = mesa.time.RandomActivation(self)

        # Create agents
        for i in range(self.num_agents):
            a = MoneyAgent(i, self)
            self.schedule.add(a)
            # Add the agent to a random grid cell
            x = self.random.randrange(self.grid.width)
            y = self.random.randrange(self.grid.height)
            self.grid.place_agent(a, (x, y))

        self.datacollector = mesa.DataCollector(
            model_reporters={"Gini": compute_gini}, agent_reporters={"Wealth": "wealth"}
        )

    def step(self):
        self.datacollector.collect(self)
        self.schedule.step()

At every step of the model, the datacollector will collect and store the model-level current Gini coefficient, as well as each agent’s wealth, associating each with the current step.

We run the model just as we did above. Now is when an interactive session, especially via a Notebook, comes in handy: the DataCollector can export the data its collected as a pandas* DataFrame, for easy interactive analysis.

*If you are new to Python, please be aware that pandas is already installed as a dependency of Mesa and that pandas is a “fast, powerful, flexible and easy to use open source data analysis and manipulation tool”. pandas is great resource to help analyze the data collected in your models

model = MoneyModel(100, 10, 10)
for i in range(100):
    model.step()

To get the series of Gini coefficients as a pandas DataFrame:

gini = model.datacollector.get_model_vars_dataframe()
# Plot the Gini coefficient over time
g = sns.lineplot(data=gini)
g.set(title="Gini Coefficient over Time", ylabel="Gini Coefficient");
../_images/f56b79804a9d37f992c55d74dd73ca8ddc789da7f535c203c6b76f5024a21aa4.png

Similarly, we can get the agent-wealth data:

agent_wealth = model.datacollector.get_agent_vars_dataframe()
agent_wealth.head()
Wealth
Step AgentID
0 0 1
1 1
2 1
3 1
4 1

You’ll see that the DataFrame’s index is pairings of model step and agent ID. This is because the data collector stores the data in a dictionary, with the step number as the key, and a dictionary of agent ID and variable value pairs as the value. The data collector then converts this dictionary into a DataFrame, which is why the index is a pair of (model step, agent ID). You can analyze it the way you would any other DataFrame. For example, to get a histogram of agent wealth at the model’s end:

last_step = agent_wealth.index.get_level_values("Step").max()
end_wealth = agent_wealth.xs(last_step, level="Step")["Wealth"]
# Create a histogram of wealth at the last step
g = sns.histplot(end_wealth, discrete=True)
g.set(
    title="Distribution of wealth at the end of simulation",
    xlabel="Wealth",
    ylabel="Number of agents",
);
../_images/871025560ca48a30b8751dca36c6c4e552bbfc345779bb3bdd5a63c8348c889a.png

Or to plot the wealth of a given agent (in this example, agent 14):

# Get the wealth of agent 14 over time
one_agent_wealth = agent_wealth.xs(14, level="AgentID")

# Plot the wealth of agent 14 over time
g = sns.lineplot(data=one_agent_wealth, x="Step", y="Wealth")
g.set(title="Wealth of agent 14 over time");
../_images/e962c886cbac1cc36ecba289b2bb2a0fe131f786347ecd55a9231803942c281c.png

You can also plot a reporter of multiple agents over time.

agent_list = [3, 14, 25]

# Get the wealth of multiple agents over time
multiple_agents_wealth = agent_wealth[
    agent_wealth.index.get_level_values("AgentID").isin(agent_list)
]
# Plot the wealth of multiple agents over time
g = sns.lineplot(data=multiple_agents_wealth, x="Step", y="Wealth", hue="AgentID")
g.set(title="Wealth of agents 3, 14 and 25 over time");
../_images/e71517e5a321df6f9116607375d431459b2292b26f69b5128c53ef980a807db5.png

We can also plot the average of all agents, with a 95% confidence interval for that average.

# Transform the data to a long format
agent_wealth_long = agent_wealth.T.unstack().reset_index()
agent_wealth_long.columns = ["Step", "AgentID", "Variable", "Value"]
agent_wealth_long.head(3)

# Plot the average wealth over time
g = sns.lineplot(data=agent_wealth_long, x="Step", y="Value", errorbar=("ci", 95))
g.set(title="Average wealth over time")
[Text(0.5, 1.0, 'Average wealth over time')]
../_images/d22c5380f833e399407a5e83c25dc68200da9250e9956af4ad4ac478d1e64162.png

Which is exactly 1, as expected in this model, since each agent starts with one wealth unit, and each agent gives one wealth unit to another agent at each step.

You can also use pandas to export the data to a CSV (comma separated value), which can be opened by any common spreadsheet application or opened by pandas.

If you do not specify a file path, the file will be saved in the local directory. After you run the code below you will see two files appear (model_data.csv and agent_data.csv)

# save the model data (stored in the pandas gini object) to CSV
gini.to_csv("model_data.csv")

# save the agent data (stored in the pandas agent_wealth object) to CSV
agent_wealth.to_csv("agent_data.csv")

Batch Run#

Like we mentioned above, you usually won’t run a model only once, but multiple times, with fixed parameters to find the overall distributions the model generates, and with varying parameters to analyze how they drive the model’s outputs and behaviors. Instead of needing to write nested for-loops for each model, Mesa provides a batch_run function which automates it for you.

The batch runner also requires an additional variable self.running for the MoneyModel class. This variable enables conditional shut off of the model once a condition is met. In this example it will be set as True indefinitely.

Additional agent reporter#

To make the results a little bit more interesting, we will also calculate the number of consecutive time steps an agent hasn’t given any wealth as an agent reporter.

This way we can see how data is handled when multiple reporters are used.

def compute_gini(model):
    agent_wealths = [agent.wealth for agent in model.schedule.agents]
    x = sorted(agent_wealths)
    N = model.num_agents
    B = sum(xi * (N - i) for i, xi in enumerate(x)) / (N * sum(x))
    return 1 + (1 / N) - 2 * B


class MoneyModel(mesa.Model):
    """A model with some number of agents."""

    def __init__(self, N, width, height):
        super().__init__()
        self.num_agents = N
        self.grid = mesa.space.MultiGrid(width, height, True)
        self.schedule = mesa.time.RandomActivation(self)
        self.running = True

        # Create agents
        for i in range(self.num_agents):
            a = MoneyAgent(i, self)
            self.schedule.add(a)
            # Add the agent to a random grid cell
            x = self.random.randrange(self.grid.width)
            y = self.random.randrange(self.grid.height)
            self.grid.place_agent(a, (x, y))

        self.datacollector = mesa.DataCollector(
            model_reporters={"Gini": compute_gini},
            agent_reporters={"Wealth": "wealth", "Steps_not_given": "steps_not_given"},
        )

    def step(self):
        self.datacollector.collect(self)
        self.schedule.step()


class MoneyAgent(mesa.Agent):
    """An agent with fixed initial wealth."""

    def __init__(self, unique_id, model):
        super().__init__(unique_id, model)
        self.wealth = 1
        self.steps_not_given = 0

    def move(self):
        possible_steps = self.model.grid.get_neighborhood(
            self.pos, moore=True, include_center=False
        )
        new_position = self.random.choice(possible_steps)
        self.model.grid.move_agent(self, new_position)

    def give_money(self):
        cellmates = self.model.grid.get_cell_list_contents([self.pos])
        if len(cellmates) > 1:
            other = self.random.choice(cellmates)
            other.wealth += 1
            self.wealth -= 1
            self.steps_not_given = 0
        else:
            self.steps_not_given += 1

    def step(self):
        self.move()
        if self.wealth > 0:
            self.give_money()
        else:
            self.steps_not_given += 1

Batch run#

We call batch_run with the following arguments:

  • model_cls The model class that is used for the batch run.

  • parameters A dictionary containing all the parameters of the model class and desired values to use for the batch run as key-value pairs. Each value can either be fixed ( e.g. {"height": 10, "width": 10}) or an iterable (e.g. {"N": range(10, 500, 10)}). batch_run will then generate all possible parameter combinations based on this dictionary and run the model iterations times for each combination.

  • number_processes If not specified, defaults to 1. Set it to None to use all the available processors. Note: Multiprocessing does make debugging challenging. If your parameter sweeps are resulting in unexpected errors set number_processes=1.

  • iterations The number of iterations to run each parameter combination for. Optional. If not specified, defaults to 1.

  • data_collection_period The length of the period (number of steps) after which the model and agent reporters collect data. Optional. If not specified, defaults to -1, i.e. only at the end of each episode.

  • max_steps The maximum number of time steps after which the model halts. An episode does either end when self.running of the model class is set to False or when model.schedule.steps == max_steps is reached. Optional. If not specified, defaults to 1000.

  • display_progress Display the batch run progress. Optional. If not specified, defaults to True.

In the following example, we hold the height and width fixed, and vary the number of agents. We tell the batch runner to run 5 instantiations of the model with each number of agents, and to run each for 100 steps.

We want to keep track of

  1. the Gini coefficient value at each time step, and

  2. the individual agent’s wealth development and steps without giving money.

Important: Since for the latter changes at each time step might be interesting, we set data_collection_period=1. By default, it only collects data at the end of each episode.

Note: The total number of runs is 245 (= 49 different populations * 5 iterations per population). However, the resulting list of dictionaries will be of length 6186250 (= 250 average agents per population * 49 different populations * 5 iterations per population * 101 steps per iteration).

Note for Windows OS users: If you are running this tutorial in Jupyter, make sure that you set number_processes = 1 (single process). If number_processes is greater than 1, it is less straightforward to set up. You can read Mesa’s how-to guide, in ‘Using multi-process batch_run on Windows’ section for how to do it.

params = {"width": 10, "height": 10, "N": range(5, 100, 5)}

results = mesa.batch_run(
    MoneyModel,
    parameters=params,
    iterations=7,
    max_steps=100,
    number_processes=1,
    data_collection_period=1,
    display_progress=True,
)

To further analyze the return of the batch_run function, we convert the list of dictionaries to a Pandas DataFrame and print its keys.

results_df = pd.DataFrame(results)
print(results_df.keys())
Index(['RunId', 'iteration', 'Step', 'width', 'height', 'N', 'Gini', 'AgentID',
       'Wealth', 'Steps_not_given'],
      dtype='object')

First, we want to take a closer look at how the Gini coefficient at the end of each episode changes as we increase the size of the population. For this, we filter our results to only contain the data of one agent (the Gini coefficient will be the same for the entire population at any time) at the 100th step of each episode and then scatter-plot the values for the Gini coefficient over the the number of agents. Notice there are five values for each population size since we set iterations=5 when calling the batch run.

# Filter the results to only contain the data of one agent (the Gini coefficient will be the same for the entire population at any time) at the 100th step of each episode
results_filtered = results_df[(results_df.AgentID == 0) & (results_df.Step == 100)]
results_filtered[["iteration", "N", "Gini"]].reset_index(
    drop=True
).head()  # Create a scatter plot
g = sns.scatterplot(data=results_filtered, x="N", y="Gini")
g.set(
    xlabel="Number of agents",
    ylabel="Gini coefficient",
    title="Gini coefficient vs. number of agents",
);
../_images/e57cfa48136b1c65e15465a036edbbfcba960e5806edc7cfe02332a6fe1a53c1.png

We can create different kinds of plot from this filtered DataFrame. For example, a point plot with error bars.

# Create a point plot with error bars
g = sns.pointplot(data=results_filtered, x="N", y="Gini", linestyle='none')
g.figure.set_size_inches(8, 4)
g.set(
    xlabel="Number of agents",
    ylabel="Gini coefficient",
    title="Gini coefficient vs. number of agents",
);
../_images/401ed8fc709f63e3c610a10b34e38c6e270e672e4f8eab809dcc68490f6c180a.png

Second, we want to display the agent’s wealth at each time step of one specific episode. To do this, we again filter our large data frame, this time with a fixed number of agents and only for a specific iteration of that population. To print the results, we convert the filtered data frame to a string specifying the desired columns to print.

Pandas has built-in functions to convert to a lot of different data formats. For example, to display as a table in a Jupyter Notebook, we can use the to_html() function which takes the same arguments as to_string() (see commented lines).

# First, we filter the results
one_episode_wealth = results_df[(results_df.N == 10) & (results_df.iteration == 2)]
# Then, print the columns of interest of the filtered data frame
print(
    one_episode_wealth.to_string(
        index=False, columns=["Step", "AgentID", "Wealth"], max_rows=10
    )
)
# For a prettier display we can also convert the data frame to html, uncomment to test in a Jupyter Notebook
# from IPython.display import display, HTML
# display(HTML(one_episode_wealth.to_html(index=False, columns=['Step', 'AgentID', 'Wealth'], max_rows=25)))
 Step  AgentID  Wealth
    0        0       1
    0        1       1
    0        2       1
    0        3       1
    0        4       1
  ...      ...     ...
  100        6       2
  100        4       0
  100        8       1
  100        3       0
  100        9       0

Lastly, we want to take a look at the development of the Gini coefficient over the course of one iteration. Filtering and printing looks almost the same as above, only this time we choose a different episode.

results_one_episode = results_df[
    (results_df.N == 10) & (results_df.iteration == 1) & (results_df.AgentID == 0)
]
print(results_one_episode.to_string(index=False, columns=["Step", "Gini"], max_rows=10))
 Step  Gini
    0  0.00
    1  0.00
    2  0.18
    3  0.18
    4  0.18
  ...   ...
   96  0.72
   97  0.72
   98  0.72
   99  0.72
  100  0.72

Analyzing model reporters: Comparing 5 scenarios#

Other insight might be gathered when we compare the Gini coefficient of different scenarios. For example, we can compare the Gini coefficient of a population with 25 agents to the Gini coefficient of a population with 400 agents. While doing this, we increase the number of iterations to 25 to get a better estimate of the Gini coefficient for each population size and get usable error estimations.

params = {"width": 10, "height": 10, "N": [5, 10, 20, 40, 80]}

results_5s = mesa.batch_run(
    MoneyModel,
    parameters=params,
    iterations=100,
    max_steps=120,
    number_processes=1,
    data_collection_period=1,  # Important, otherwise the datacollector will only collect data of the last time step
    display_progress=True,
)

results_5s_df = pd.DataFrame(results_5s)
# Again filter the results to only contain the data of one agent (the Gini coefficient will be the same for the entire population at any time)
results_5s_df_filtered = results_5s_df[(results_5s_df.AgentID == 0)]
results_5s_df_filtered.head(3)
RunId iteration Step width height N Gini AgentID Wealth Steps_not_given
0 0 0 0 10 10 5 0.00 0 1 0
8 0 0 1 10 10 5 0.00 0 1 1
11 0 0 2 10 10 5 0.32 0 0 0
# Create a lineplot with error bars
g = sns.lineplot(
    data=results_5s_df,
    x="Step",
    y="Gini",
    hue="N",
    errorbar=("ci", 95),
    palette="tab10",
)
g.figure.set_size_inches(8, 4)
plot_title = "Gini coefficient for different population sizes\n(mean over 100 runs, with 95% confidence interval)"
g.set(title=plot_title, ylabel="Gini coefficient");
../_images/49d47f9e197252993178c28692fd6077ee93a3c2951b8a46d50bd21daa52cda5.png

In this case it looks like the Gini coefficient increases slower for smaller populations. This can be because of different things, either because the Gini coefficient is a measure of inequality and the smaller the population, the more likely it is that the agents are all in the same wealth class, or because there are less interactions between agents in smaller populations, which means that the wealth of an agent is less likely to change.

Analyzing agent reporters#

From the agents we collected the wealth and the number of consecutive rounds without a transaction. We can compare the 5 different population sizes by plotting the average number of consecutive rounds without a transaction for each population size.

Note that we’re aggregating multiple times here: First we take the average of all agents for each single replication. Then we plot the averages for all replications, with the error band showing the 95% confidence interval of that first average (over all agents). So this error band is representing the uncertainty of the mean value of the number of consecutive rounds without a transaction for each population size.

# Calculate the mean of the wealth and the number of consecutive rounds for all agents in each episode
agg_results_df = (
    results_5s_df.groupby(["iteration", "N", "Step"])
    .agg({"Wealth": "mean", "Steps_not_given": "mean"})
    .reset_index()
)
agg_results_df.head(3)
iteration N Step Wealth Steps_not_given
0 0 5 0 1.0 0.0
1 0 5 1 1.0 0.8
2 0 5 2 1.0 1.4
# Create a line plot with error bars
g = sns.lineplot(
    data=agg_results_df, x="Step", y="Steps_not_given", hue="N", palette="tab10"
)
g.figure.set_size_inches(8, 4)
g.set(
    title="Average number of consecutive rounds without a transaction for different population sizes\n(mean with 95% confidence interval)",
    ylabel="Consecutive rounds without a transaction",
);
../_images/7b94a2bf3b17d97b99041beec9650fdd117adc0430ea81c077d742bc3b8b4055.png

It can be clearly seen that the lower the number of agents, the higher the number of consecutive rounds without a transaction. This is because the agents have fewer interactions with each other and therefore the wealth of an agent is less likely to change.

General steps for analyzing results#

Many other analysis are possible based on the policies, scenarios and uncertainties that you might be interested in. In general, you can follow these steps to do your own analysis:

  1. Determine which metrics you want to analyse. Add these as model and agent reporters to the datacollector of your model.

  2. Determine the input parameters you want to vary. Add these as parameters to the batch_run function, using ranges or lists to test different values.

  3. Determine the hyperparameters of the batch_run function. Define the number of iterations, the number of processes, the number of steps, the data collection period, etc.

  4. Run the batch_run function and save the results.

  5. Transform, filter and aggregate the results to get the data you want to analyze. Make sure it’s in long format, so that each row represents a single value.

  6. Choose a plot type, what to plot on the x and y axis, which columns to use for the hue. Seaborn also has an amazing Example Gallery.

  7. Plot the data and analyze the results.

Happy Modeling!#

This document is a work in progress. If you see any errors, exclusions or have any problems please contact us.

[Comer2014] Comer, Kenneth W. “Who Goes First? An Examination of the Impact of Activation on Outcome Behavior in AgentBased Models.” George Mason University, 2014. http://mars.gmu.edu/bitstream/handle/1920/9070/Comer_gmu_0883E_10539.pdf

[Dragulescu2002] Drăgulescu, Adrian A., and Victor M. Yakovenko. “Statistical Mechanics of Money, Income, and Wealth: A Short Survey.” arXiv Preprint Cond-mat/0211175, 2002. http://arxiv.org/abs/cond-mat/0211175.