SimPy Simulation — Dynamic Resources

Introduction

Do you like complex network problems? DO YOU LOVE THINKING ABOUT SUPPLY AND DEMAND? Do you HATE LINES? If you answered yes to any of those questions, you have come to the right place.

In my last two articles we went through what a simulation is and some key features around the SimPy simulation we will be leveraging to experiment in our Marketplace.That’s really the end product here, our desire is to test:

• The best routing logic.
• How new supply OR demand impacts our wait times.
• The best way to experiment — what kind of switchback windows do we need and how long before we see significant results.
• Really any change that will impact the network — the simulation will help us unpack how a change will impact our core metrics.

Today, we’re venturing deeper into the rabbit hole, the goal of this article is to give you the power to spin up your own simulation OR pitch it to your boss (so I can build it).

Why SimPy for Healthcare Simulation?

SimPy, a Python-based discrete-event simulation framework, is our chosen instrument for this complex task. Why SimPy? Because it’s adept at mimicking the unpredictable nature of healthcare systems. Its ability to model and manage time-dependent behaviors makes it an excellent tool for simulations where timing and resource allocation are crucial.

Overall Flow

1. Data Preparation and Ingestion: Data from CSV files to set up our inputs and generate probability distributions.
2. Clinician and Consult Classes — The Building Blocks: Clinician and Consult classes reflect real-world entities with attributes.
3. Core Simulation Processes: Functions like generate arrivals, process_consult, sync_online, manage consult creation, completion and clinician availability within the simulation.
4. Routing Logic: Depending on your situation this could be like a FIFO, or as unique as routing to the tallest clinician.
5. Running the Simulation: The run method activates the simulation, establishing its timeline and process flow.
6. Measuring Outcomes — Metrics and Outputs: Key performance metrics are captured to assess and refine the simulation’s effectiveness.

Data Preparation and Ingestion

We aren’t going to spend much time here — as everyone’s data is different. But just so you can wrap your head around our data structures, here’s the highlights:

• Clinicians — they can service many states and many clients. They have attributes like [TX, CA, GA], [Client A, Client B, Client C], [11/17/2023 7AM, 11/17/2023 10PM],
• We translate every clinician shift to simulation time — so if our simulation starts at 11/17/2023 7AM, that translates to 0th second in the simulation. The shift end time is the 54000th second, [0.54000].
• Consults — these represent patients needing care. They have attributes like State, Client, Arrival Time, Modality, SLA

Wait Times

We have a ton of data around each action a clinician takes to complete a consult. We could spend an entire article talking about probability distributions but we will just put it in a black box here. We calculate a discrete probability distribution for each category of wait times (assign to accept, accept to start, start to complete) for each client-modality pair. This distribution is then used to generate random wait times for our simulation, aiming to mimic real-world variability. We also remove some wild outliers to keep it from going bonkers.

Clinician and Consult Classes — The Building Blocks

In our simulation, the Clinician and Consult classes are not just passive data structures; they are active participants in the SimPy simulation environment. SimPy excels in handling complex interactions and events, and these classes are designed to plug directly into this environment, allowing us to model real-world scenarios with impressive fidelity.

Clinician Class:

Each clinician has attributes like the ones listed above, dictating their specialization and availability. But crucially, they also have SimPy resources which are key to managing their availability in the simulation. We leverage the class to store what happened in the simulation, so any consult taken, how long each one took, etc. This allows us to validate the simulation in the end and to verify if our experiment had the correct or ill effects.

Notice in our class pass a bunch of information to it. A key piece of information we must have is self, env that we pass to the class. We will do this for each of our classes.

self: A reference to the specific instance of the class itself, used to access variables and methods of that instance.
env: In SimPy, it’s the simulation environment where the processes and events of the simulation are managed and executed.

Class Clinician:
    def __init__(self, env,queue,clin_id, name, licenses, sync_schedule,clients):
        self.env = env
        self.id = clin_id
        self.name = name
        self.licenses = [] if licenses is None else licenses
        self.schedule = [] if sync_schedule is None else sync_schedule
        self.clients = [] if clients is None else clients
        self.sync_cmpl_time = -1
        self.busy = None
        self.sync_busy = None
        self.sync_resource = simpy.Resource(env, capacity=1)
        self.consults_completed = 0
        self.expected_lengths = 0
        self.queue = []
        self.bookings =[]
        self.consults = []
        self.time_in_consults = 0
        self.earnings = []
        self.time_to_accept = []
        self.time_to_complete = []

Consult Class

Each consult carries attributes like state, client, modality but also contains dynamic attributes like priority, that change over time depending on the routing logic you want to use. This is also where we pull in the discrete probability distribution talked about above — we generate a number (representing seconds) based on the distribution of the client/modality. Just like in the clinician class we store when a consult is created, started, completed, etc.

 def __init__(self, state, client,sla,length,clinician,dispatch_time,p_id,created_time):
        self.state = state
        self.client = client
        self.id = p_id
        self.sla = sla
        self.modality = modality
        self.length = length
        self.clinician = clinician
        self.dispatch_time = dispatch_time
        self.wait_accept_duration = 0.0
        self.wait_start_duration = 0.0
        self.wait_complete_duration = 0.0
        self.created_time = created_time
        self.priority = 0.0
        self.booking_time = 0.0
        self.assign_time = 0.0
        self.start_time = 0.0
        self.complete_time = 0.0
        self.queue_time = 0.0
        self.created_start_time = 0.0
        self.create_assign_time = 0.0
        self.assign_accept_time = 0.0
        self.accept_start_time = 0.0
        self.start_end_time = 0.0
        self.time_in_system = 0.0
        self.breached = 0
        self.ignored = 0

Wheel Model Class

In SimPy, the environment (env) is a core component that orchestrates the timing and interaction of events and processes. In the Wheel model class here’s why

1. Central Simulation Coordinator: The SimPy environment created within the wheel_model class (self.env = simpy.Environment()) acts as the central coordinator for all simulation activities. It keeps track of the simulation clock and controls the execution of events over time.
2. Initialization (__init__ method): Sets up the simulation environment (env), initializes various counters (like consult_counter, sync_clin_counter), creates a store (staging_store), and sets up lists and distributions related to consults (queue, consults_list, group_bin_probabilities).

class wheel_model:
    #Set the environment here, Also initialize the store from the simpy Store class
    def __init__(self,run_number):
        self.env = simpy.Environment()
        self.consult_counter = 0
        self.run_number = run_number
        self.sync_clin_counter = 0
        self.staging_store = simpy.Store(self.env)
        self.queue = []
        self.consults_list = []
        self.group_bin_probabilities = wait_time_distribution.group_bin_probabilities

Interplay Between Classes and SimPy:

The true magic happens in how these classes interact within the SimPy environment. For instance, when a Consult instance is created, it’s not just a static record. It’s an event in the simulation, triggering processes like consult creations, clinician assignment, routing calculation, and more. Similarly, each Clinician isn’t just a set of attributes; they’re agents in the simulation, with their availability and actions governed by SimPy’s event-driven architecture.

The classes define the data structures, creating ways for us to have dynamic resources with very specific characteristics ALL orchestrated by SimPy’s powerful simulation capabilities.

Section 4: Core Simulation Processes

Several processes keep our simulation ticking:

generate_arrivals: Here, we simulate the arrival of consultations, adding them to our system with a splash of randomness.

while self.env.now < run_time:
  if self.consult_counter < len(df_consults):
      for index, row in df_consults.iterrows():
          # Schedule the consult to arrive at the specified time
          self.consult_counter += 1 
          consult_init = Consult(row['State'],row['Name'],row['sla'],row['Length'],row['dispatch_time'],self.consult_counter,self.env.now)
          print(f'Consult {consult_init.id}, {consult_init.modality} has been created at {self.env.now}')
          self.consults_list.append(consult_init)

          #Kicks of finding a clinician that isn't busy & didn't complete a consult in this second
          clinician = None
          # Sort the clinicians by priority
          sorted_clinicians = sorted(clinicians, key=lambda x: (random.randint(1, 10000)))
          # Check if there are any available clinicians, make sure the clinician doesn't have a queue (modality filtered specific to clinician license and clients)
          for c in sorted_clinicians:
              if c.available(consult_init,self.env.now) and len([con for con in self.queue if and con.state in c.licenses and con.client in c.clients]) == 0:
                  clinician = c
                  break
          # If there are available clinicians, kick off the new consult to clinician process
          if clinician is not None:
              self.new_consult_to_clinician(consult_init,clinicians,clinician)
          # If there are no available clinicians, add the consult to the queue
          else:
              #Append the consult to a queue to be processed when a clinician becomes available
              self.queue.append(consult_init)
              yield self.env.timeout(row['duration'])
  else:
      #we need this to keep the simulation running even after the consults finish arriving
      yield self.env.timeout(1)

sync_online: These processes manage when clinicians clock in for their shift, ready to take on consults.

def sync_online(self,clinicians,run_time):
  while self.env.now < run_time:
      if self.sync_clin_counter < len(df_sync_clins):
          for index, row in df_sync_clins.iterrows():
              self.sync_clin_counter += 1 
              #the input data is clinician IDs and times they came online. We need to find the clinician object that matches the ID.
              for clinician in clinicians:
                  if clinician.id == row['ID']:
                      clinician = clinician
                      break
              #wait a period of time before the next clinician arrives for async.
              print(f'{clinician.name} is now online for sync {self.env.now}')
              #Go check to see if a consult is in the queue
              self.queued_consult(clinician,clinicians,"sync")
              yield self.env.timeout(row['duration'])
      else:
          #This statement will keep the while loop alive and the clock ticking.
          yield self.env.timeout(1)

process_consult: This is where the magic happens, each consult is processed, mirroring the real-world.

 consult.value.wait_accept_duration = self.generate_wait_time(self.group_bin_probabilities, consult.value.client, consult.value.modality,'assign_to_accept_duration')
        consult.value.wait_start_duration = self.generate_wait_time(self.group_bin_probabilities, consult.value.client, consult.value.modality,'accept_to_start_duration')
        consult.value.wait_complete_duration= self.generate_wait_time(self.group_bin_probabilities, consult.value.client, consult.value.modality,'start_to_complete_duration')

        
        if consult.value.modality == 'sync':
            resource = clinician.sync_resource
            busy_attr = "sync_busy"
            rplc_modality ='sync'

            # Request the clinician resource
        with resource.request() as req:
            # Yield statements are like return statements but more like a pause to do something.
            yield req
            # Set sync busy not busy states.
            setattr(clinician, busy_attr, consult)
            # Trigger clinician to wait
            consult.value.assign_time = self.env.now
            yield self.env.timeout( consult.value.wait_accept_duration)
            #Record Start Time
            consult.value.start_time = self.env.now
            yield self.env.timeout( consult.value.wait_complete_duration)
            # After clinician waits the timeout period, return here and free the clinician up
            setattr(clinician, busy_attr, None)
            consult.value.complete_time = self.env.now
            consult.value.create_assign_time = consult.value.assign_time - consult.value.created_time
            consult.value.queue_time = consult.value.start_time - consult.value.created_time
            consult.value.time_in_system = consult.value.complete_time - consult.value.created_time
            consult.value.start_end_time = consult.value.complete_time - consult.value.start_time
            consult.value.created_start_time = consult.value.start_time - consult.value.created_time
            time_on_clinician = consult.value.complete_time - consult.value.assign_time
            consult.value.clinician = clinician.name
            clinician.consults.append(consult.value.state+consult.value.modality+str(consult.value.id))
            clinician.time_in_consults += time_on_clinician
            clinician.consults_completed += 1
            clinician.last_consult_time = self.env.now
            clinician.expected_lengths += consult.value.length
           
            elif consult.value.created_start_time > consult.value.sla:
                consult.value.breached = 1

            print(f'Consult {consult.value.id} for {consult.value.state}, {consult.value.modality}, {consult.value.priority}, {consult.value.modality} completed at {self.env.now:.1f} with {clinician.name}{new_line}')
            #Checks the consult queue to see if a clinician is available
            self.queued_consult(clinician,clinicians,rplc_modality)

Section 5: Assigning a Clinician

We designed this simulation to be able to assign a clinician based on their state licenses, clients, schedule and some other factors. We first do a check to make sure that the clinician meets all these conditions and if they do, then we choose a clinician based on specific routing logic. This could range from choosing the lowest utilized clinician a high performer or even just at random.

This code does all the checks to make sure we can assign a consult based on the clinicians shift, licenses, clients. This takes into account the end of shift as well.

def available_schedule(self, consult,now, check_type):

    if consult.modality == 'sync':
        #make sure they aren't busy & that they haven't completed another sync consult in the same instance of another consult being created
        if self.schedule != 0:
            for time in self.schedule:
                start = time[0]
                end = time[1]
                if now >= start and now <= end:
                    return True
        return False
    def available(self, consult,now):
        # Check if the clinician has the required licenses
        if consult.state not in self.licenses:
            return False 
        # Check if the clinician is able to serve the client
        if consult.client not in self.clients:
            return False

Section 6: Routing

Creating custom routing logic in SimPy was a bit of a hassle. There are ways to create prioritized consults but I couldn’t find a way to reorder the queue in a specific way. FIFO is the standard but one of my key requirements was to test routing logics. The way we do it in SimPy is to create a staging queue that we order in the exact way we want to. Then we pop the top consult out of the staging queue, and place it into the SimPy queue for the system to process.

 def queued_consult(self,clinician,clinicians,last_con_mod):
    consult = None
    
    #Filter the queue for the specific clinician (this is a temp queue, it is reset everytime) + Filter by modality
    self.filtered_consults = [c for c in self.queue if c.modality == last_con_mod and c.state in clinician.licenses and c.client in clinician.clients]

    #Check to see if a clinician has a queue at all, filter first.
    if self.filtered_consults:
        # Filter the consults for the clinician
        for con in self.filtered_consults:
            #THIS IS WHERE YOU WOULD CREATE CUSTOM LOGIC
            #You could assign a score based on when it was created
            #You could assign a score based on what state it is in
            #I removed my logic (secret sauce)
            con.priority = random.randint(1, 10000)
        #Sort the consults, descending order
        self.filtered_consults = sorted(self.filtered_consults, key=lambda x: x.priority, reverse=True)
        
        #This is round about, but big goal is to get the system to recognize the next consult.
        # 1. pop out the highest priority consult | 2. remove the consult from the big queue | 
        # 3. Put it in the Simpy Store, (Stores are predefined simpy classes) | 4. Get the consult from the store
        if clinician.available_schedule(self.filtered_consults[0],self.env.now,"clinician"):
        
            next_consult = self.filtered_consults.pop(0)
            self.queue.remove(next_consult)
            self.staging_store.put(next_consult)
            consult = self.staging_store.get()
            self.filtered_consults = []
            
            # Remember that function we made above process_consult. This is what we call to get us there
            self.env.process(self.process_consult(clinician,consult,clinicians))

Section 7: Running the Simulation

The run method is where our simulation comes to life. Here, we set the timeline, kick-starting the processes and recording key moments to measure. This is where we set the # of runs and call our data processing step to set up all of our inputs. Here’s how we set everything in motion.

for run in range(1,30):
    #print (f"Run {run+1} of {2}")
    my_wheel_model = wheel_model(run)
    my_wheel_model.run()

Section 8: Measuring Outcomes — Metrics and Outputs

What’s a simulation without metrics? We track various outcomes like consultation completion times, clinician utilization rates, and more. These metrics are not just numbers; they’re insights, helping us understand the efficiency and effectiveness of our simulated environment. We export these into CSV files for a deep dive into the analytics, we can do more of this in Python or hand off to another team to unpack.

def consults_to_csv(self,consults, filename):
  with open(filename, mode='w', newline='') as csv_file:
      fieldnames = ['state', 'client', 'id', 'sla', 'modality', 'length', 'clinician', 'dispatch_time', 'p_id', 'created_time', 'priority', 'booking_time','assign_time', 'start_time', 'complete_time', 'queue_time', 'create_assign_time', 'start_end_time', 'time_in_system', 'breached',]
      writer = csv.DictWriter(csv_file, fieldnames=fieldnames)

      writer.writeheader()
      for consult in consults:
          writer.writerow({
              'state': consult.state,
              'client': consult.client,
              'id': consult.id,
              'sla': consult.sla,
              'modality': consult.modality,
              'length': consult.length,
              'clinician': consult.clinician,
              'dispatch_time': consult.dispatch_time,
              'p_id': consult.id,
              'created_time': consult.created_time,
              'priority': consult.priority,
              'assign_time': consult.assign_time,
              'start_time': consult.start_time,
              'complete_time': consult.complete_time,
              'queue_time': consult.queue_time,
              'create_assign_time': consult.create_assign_time,
              'accept_start_time': consult.accept_start_time,
              'start_end_time': consult.start_end_time,
              'time_in_system': consult.time_in_system,
              'breached': consult.breached,
          })

Section 9: Putting it all together.

I think it is helpful to see the full code structure to help visualize how everything flows. I am not going to put each sections full code, but you will see the function headers and be able to tie everything back.

import simpy
import csv
import random
import pandas as pd
import math
from operator import attrgetter
from input_data import arrival_times, clinicians_table
import numpy as np
import wait_time_distribution
import gspread
from oauth2client.service_account import ServiceAccountCredentials


#Data INPUTS

new_line = '\n'

#Define Clinicians
class Clinician:
    def __init__(self, env,queue,clin_id, name, licenses, sync_schedule,clients):
    ....
    .... etc.

    def available_schedule(self, consult,now, check_type):
    ....
    .... etc.   

    #Available function defines if a clinician can take a consult.
    def available(self, consult,now):
 

#Define Consults

class Consult:
    def __init__(self, state, client,sla,modality,length,clinician,dispatch_time,p_id,created_time):
    ....
    .... etc.    
            
class wheel_model:
    #Set the environment here, Also initialize the store from the simpy Store class
    def __init__(self,run_number):
        .... etc        

    def generate_wait_time(self,group_bin_probabilities, client,duration_label):                   
       ...etc

    def process_consult(self, clinician, consult,clinicians):
       ....etc 
         
    #SYNC This triggers when a new clinician comes online. It goes and checks a queue first.
    def sync_online(self,clinicians,run_time):
        .... etc


    def generate_wl_arrivals(self,clinicians,run_time):
        .....etc
   
#The Queued Consult Process. If a clinician finishes a consult, go to this process to check if they have a queue.
    def queued_consult(self,clinician,clinicians,last_con_mod):
        ....etc

    def consults_to_csv(self,consults, filename):
        ....etc
    
    
    def run(self):
        run_time = 1208400
        clinicians = []
        for index, row in df_clinician_table.iterrows():
            clinician = Clinician(self.env, self.queue, row['ID'], row['name'],
                            row['licenses'], row['sync_schedule'],
                            row['clients'])
            clinicians.append(clinician)

        self.env.process(self.sync_online(clinicians,run_time))
        self.env.process(self.toggle_online(clinicians,run_time))
        #kick of the first process. This initiates the sim
        self.env.process(self.generate_wl_arrivals(clinicians,run_time))
        #how long should the sim run
        self.env.run(until=run_time)
        self.consults_to_csv(self.consults_list, 'consult_results.csv')
             
for run in range(1,30):
    my_wheel_model = wheel_model(run)
    my_wheel_model.run()
    #print ()

Closing Remarks

This should give you a good idea of how to spin up a simulation in SimPy! Refer back to my last article where I link some helpful resources in getting started in a simulation. Coding along side someone is probably one of the most helpful exercises you can do.

My recommendation as you go on this journey is to start small. I built the simulation first with a super simple version, then I kept adding in complexities.

What’s next??!?

We are going to test some switchbacks! I think it would be brilliant to use a simulation to help us determine time to significance, switchback windows and all that jazz.