Get Vaccinated With Python

Are you wondering how long it will take you to go and get that dosage of the COVID-19 vaccine at your nearest health facility in your community? Or maybe you want to help your community get vaccinated by providing them with the much needed resources that are essential during this time just like the COVID-19 Hub for Zimbabwe which has been helping Zimbabweans better understand the pandemic in their local context.

Well, in this article (in tutorial form) we are going to be using some very useful tools and show how Python can just help save the world somewhere out there in the world. All we want is to get back to the way things were before Coronavirus showed up.

We are going to be leveraging on three things:

The first thing I like about Python is how much you can get to manipulate data and get to do some really amazing and crazy stuff with it. This time we are not going to be getting any crazy ideas other than just trying to save the world starting with our very own communities.

Just a lil prep

First of all before we start getting deep into coding, make sure this checklist is marked complete otherwise you are going to face a little challenges along the way: Make sure you

  1. have Python installed
  2. have Jupyter Notebooks up and running
  3. have an ArcGIS Developer Account subscription (it is totally free do not worry)
  4. like saving the world (you really need to love this)

Let’s Start

Now that your checklist is marked at 100% let us begin.

Using Jupyter Notebook you need to copy and paste the code snippets that I am going to be placing along the way. I will also do a lil explanation below each and every code snippet that is involved in this process along the way.

from arcgis.gis import GIS
import arcgis.network as network
from arcgis.features import FeatureLayer, Feature, FeatureSet, use_proximity, FeatureCollection
import pandas as pd
import time
import datetime as dt

# The environment settings library
from decouple import config

Here we summoned the ‘soldiers’ that we need in order to prepare for our journey. From the GIS library we have imported:

  • GIS – for authentication
  • Network – for peforming network analysis functions
  • …. much of the modules are self explanatory here
  • Pandas – we will use this for displaying our raw data and a little manipulation

and Decouple is going to be my little secret. I am going to be placing my ArcGIS Developer account credentials in environment variables so that no-one else can have access to my account other than me.

To know how to hide your credentials too I recommend reading the official project documentation and guide.

# Connect to the GIS
gis = GIS(url='https://arcgis.com', username=config('GIS_USERNAME'), password=config('GIS_PASSWORD'))

We have been authenticated and authorized to access and use the API with our account after having provided username and password above. Some prefer using the Token instead, you can always get that from the Developer Dashboard but I am a bit lazy navigating around so I just use the easy method.

vaccination_resources = gis.content.get("f464c7f46aea483694cf2ed13e6471f4")
display(vaccination_resources)

The Coronavirus Hub also contains locations of all vaccination centers in Zimbabwe. We just called the Hosted Feature Layer using its ID.

The COVID-19 vaccination centers feature layer

You should be seeing the same display as above if you are using Jupyter Notebooks.

vac_centers = vaccination_resources.layers[0]
vac_centers

Let us get the services URL of the Feature layer in case anyone else might just need it. The URL:

<FeatureLayer url:"https://services8.arcgis.com/oTalEaSXAuyNT7xf/arcgis/rest/services/zwe_vaccination_centers_/FeatureServer/0">
from arcgis.geocoding import geocode
my_location = input("What is your current location?")

We imported the geocoding module which helps us turn an address into an X, Y coordinate location. So in this case, the user is going to be asked to type their location and the code basically translates their input into X, Y since we cannot directly ask anyone to input their coordinates.

results = geocode(my_location)
# query the first matched result
coordinates = results[0]['location']
coordinates

Now we have called the geocode function which takes one parameter which is the address. We are only going to be taking the first result in this case hence we used the index of zero [0]. By default, the best result gets the top priority index value so that’s what we basically need for now. Just to display the coordinates call the variable:

{'x': 27.750630000000058, 'y': -19.76296999999994}
# --> For Tsholotsho (in Zimbabwe)
x_coordinate = coordinates['x']
y_coordinate = coordinates['y']

Now let’s place the coordinates into their respective variables individually so we can use them later. We are extracting the coordinate value of both X and Y from the JSON dictionary provided earlier on.

incidents_json = {
                    "features": [{"attributes": {"CurbApproach": 0,
                                                 "ID": "COVID001",
                                                 "Name": "User Location"},
                                  "geometry": {"x": x_coordinate, "y": y_coordinate}}
                                 ],
                    "spatialReference": {"wkid": 4326, "latestWkid": 4326},
                    "geometryType": "esriGeometryPoint",
                    "fields" : [
                        {"name" : "ID", "type" : "esriFieldTypeString", "alias" : "ID", "length" : "50"},
                        {"name" : "Name", "type" : "esriFieldTypeString", "alias" : "Name", "length" : "50"},
                        {"name" : "CurbApproach", "type" : "esriFieldTypeInteger", "alias" : "CurbApproach"}
                    ]
                }
incidents = FeatureSet.from_dict(incidents_json)

We have created an incidence layer. This layer is going to be based on the user’s input which has been geocoded and is now ready to be involved in some awesome network stuff. We do this by creating a JSON object which contains the coordinate values in the geometry parameter.

We also need to define a few items in this JSON file which will help for the actions we are going to be preforming next and that is the:

  • spatial reference
  • geometry type

In the end (last line of code), we converted that dictionary we just created into a layer which is going to be our incidence layer for the Closest Facility Network analysis we are going to be performing which requires an incidence layer and a facility layer.

map1 = gis.map('Zimbabwe', zoomlevel=6)
map1

user_location = {"type":"esriPMS",
                           "url":"https://image.flaticon.com/icons/png/512/10/10966.png", #User symbol from flaticon
                           "contentType": "image/png", "width":20, "height":20}

map1.draw(incidents, symbol=user_location)

Let’s just see what we have for now on a map display. You should be seeing a map displayed marking the user-location based on the location that you have entered in the input box.

%%time
current_time = dt.datetime.now()  
result1 = network.analysis.find_closest_facilities(incidents=incidents, facilities=vac_centers, 
                                                   cutoff=25, time_of_day=current_time, 
                                                   number_of_facilities_to_find=1,
                                                   save_output_network_analysis_layer=True,
                                                   gis=gis)

Running a Closest Facility Network Analysis with the snippet above. Remember our vac_centers layer that we had earlier on? We have used that to be our facilities layer just like we mentioned before and the user location as the incidence layer. We did a couple of things here. We,

  • called the network analysis method which then called the closest_facilities function which takes several parameters
  • in these parameters we have tasked the function to find just 1 facility for us (you can obviously change that)
  • we have also asked the function to save this network analysis layer and reminded it that we are the previously authenticated user through the gis method.
print("Is the tool executed successfully?", result1.solve_succeeded)

Let’s ask Python if the tool was executed successfully. If YES, then a TRUE result should be displayed below.

result1

The code snippet above is optional. Just for those who want to explore more about the generated network analysis result. It will provide you with information on the result generated and the URL’s to the layer.

result1.output_network_analysis_layer.url

Again, just maybe you might need this layer, the snippet above will help you grab it.

df4 = result1.output_routes.sdf
start_times = pd.to_datetime(df4["StartTime"], unit="ms")
end_times = pd.to_datetime(df4["EndTime"], unit="ms")
df4["StartTime"] = start_times.apply(lambda x: x.strftime("%H:%M:%S"))
df4["EndTime"] = end_times.apply(lambda x: x.strftime("%H:%M:%S"))
# df4[["Name", "StartTime", "EndTime", "IncidentID", "Total_Kilometers", "Total_Minutes"]]
df4.head()

Remember Pandas as pd? We will use this to extract the information we need from the result and be able to create human-readable format of this data.

The most important thing the users might need to know here is the distance (in kilometers) to the nearest facility and the time they will take in order for them to get there.

If you explore the dataFrame a little further, you will notice a column named SHAPE. This column contains the path and routes to the nearest facility. We will be visualizing that soon towards the end which we are close to by the way.

map2 = gis.map(my_location, zoomlevel=10)
map2

Let’s create a second map and this time focus the zoom extent to the user-location area.

# draw the facilities and incidents on this map
map2.draw(incidents, symbol=user_location, attributes={"title":"Incident Location"})

# draw the routes on the map
map2.add_layer(result1.output_routes)

We have drawn the routes onto the map now in addition to what we previously had before but just on a new map.

start_time = df4['StartTime'][0]
end_time = df4['EndTime'][0]
distance = df4['Total_Kilometers'][0]
total_time = round(df4['Total_Minutes'][0], 2)

We’ve then created some variables by extracting column fields from out DataFrame in a mission to provide human-readable format information to the user’s.

So we would need something like;

print(f"Based on your location it will take you {total_time} minutes to travel to the nearest Vaccination Center which is {distance} kilometers. You estimated arrival time is {end_time}")

Which will then display a message like

Based on your location it will take you 3.78 minutes to travel to the nearest Vaccination Center which is 1.6682732175993054 kilometers. You estimated arrival time is 09:29:47

…to our user’s.

The End

Wow. I hope you really had fun in using the COVID-19 Hub resources together with the ArcGIS API for Python to help people or someone close to you know how far they need to travel to get to the closest vaccination centers in order to stay safe against COVID-19.

Let’s hope you get to turn this into something real big and useful for the community and not just for COVID-19 but for anything which could be really fun and helpful too.