Skip To Content ArcGIS for Developers Sign In Dashboard

ArcGIS API for Python

Guide to Network Analysis (Part 3 - Generate Service Area)

Introduction

Now we have learned about Network Datasets and Network Analysis Services, and how to find routes from one point to another, and among multiple points, let's move onto the third topic - how to generate Service Area. Please refer to the road map below if you want to revisit the previous topics or jump to the next topic -

  • Network Dataset and Network Analysis services (Part 1)
  • Find Routes (Part 2)
  • Generate Service Area (You are here!)
  • Find Closest Facility (Part 4)
  • Generate Origin Destination Cost Matrix (Part 5)
  • Solve Location Allocation (Part 6)
  • Vehicle Routing Problem Service (Part 7)

This part of guide to Networking Analysis will help you find out what's occurring within a set distance of a feature, or what's within traveling range (measured using distance, time, or cost), and bring you these benefits [1]:

  • allow you to identify the area - and the features inside the area - affected by an event or activity. For instance, a city planner needs to notify residents within 500 feet of a proposed liquor store of its forthcoming.
  • let you monitor activities in the area. For instance, staff working at state forestry department needs to monitor logging that should not happen within 100 meters buffer along streams.
  • help define the area served by a facility. For example, the chief of fire department wants to know which streets are within 3-minute drive of a fire station, and/or a GIS analyst in search of a good retail store location wants to find out how many people live within a 15-minute drive of a proposed store location.
  • help delineate areas that are suitable for/capable of supporting a specific use. Take a wildlife biologist as an example, he/she might need to map the area within a half-mile of streams combined with vegetation type, slope, and other factors that could be used to identify prime deer habitats.

Fig 1. The department of fire and rescue in Prince William County , VA, mapped areas within 5, 10, and 15 minutes of drive to fire stations, in order to help decide where to build new stations (source: [1]).

To find what's nearby, the first step is to decide how to measure the "nearness". You can measure by straight-line distance, distance or cost over a network, or cost over a surface [1].

  • When what's nearby is defined by straight-line distance, you need to specify the source feature and the distance, then the ArcGIS API for Python finds the area or the surrounding features within the distance. This approach is good for creating a boundary or selecting features at a set distance around a source (as shown in Fig 2(L)).
  • When what's nearby is defined by distance or cost over a network, you need to specify the source locations and a distance (or travel cost) along each linear feature, then the API finds which segments of the network are within the distance/cost. This approach is good for finding what's within a travel distance or cost of a location, over a fixed network (as shown in Fig 2(M)).
  • When what's nearby is defined by cost over a surface, you need to specify the location of the source features and a travel cost, then the API creates a new layer showing the travel cost from each source feature. This approach is good for calculating overland travel cost (as shown in Fig 2(R)).

Fig 2. (Left) Creates a buffer by straight-line distance measurement; (Middle) Calculates distance or cost over a network; (Right) Uses the travel cost map based on slopes from source to destinations (Source: [1]).

Service areas are commonly used to visualize and measure accessibility. For example, a three-minute, drive-time polygon around a grocery store can determine which residents are able to reach the store within three minutes and are thus more likely to shop there [2]. Service areas also help evaluate accessibility. Concentric service areas show how accessibility varies. Once service areas are created, you can use them to identify how much land, how many people, or how much of anything else is within the neighborhood or region. Service area solver provides functionality for finding out how far a vehicle could go within a specified time or distance limit.

Now that we have learned why generating service area is needed, and the basis of measuring the nearness, it is time to play with some code!

Data

The first step to everything is always importing all required modules and establishing a GIS connection to ArcGIS Online organization or your enterprise for ArcGIS.

If you have already set up a profile to connect to your ArcGIS Online organization, execute the cell below to load the profile and create the GIS class object. If not, use a traditional username/password log-in e.g. my_gis = GIS('https://www.arcgis.com', 'arcgis_python', 'P@ssword123', verify_cert=False, set_active=True)

In [1]:
from arcgis.gis import GIS
import arcgis.network as network
from arcgis.features import FeatureLayer, Feature, FeatureSet, use_proximity
import pandas as pd
import datetime as dt
import time
In [2]:
my_gis = GIS(profile="your_online_profile")

Define Facilities Layer

If you do not have access to the hospital layer provided by Esri_livingatlas as referenced in the cell below, an alternative approach is for you to download the hospital listing of San Bernadino County from this source and publish the csv to the organization before proceeding forward.

In [3]:
""" This try-except block will help you download the CSV and publish to current GIS object, if 
    "You do not have permissions to access this resource or perform this operation."
"""
try:
    hospital_item = my_gis.content.get("a2817bf9632a43f5ad1c6b0c153b0fab")
except RuntimeError as ne:
    try:
        print("Trying from an alternative source...")
        hospital_item = my_gis.content.get("50fb63f303304835a048d16bd86c3024")
    except RuntimeError as ne:
        print("Trying to publish from csv...")
        import requests
        import csv
        import os

        out_file_name = 'hospitals_SB_county.csv'
        url = "https://data.chhs.ca.gov/datastore/dump/641c5557-7d65-4379-8fea-6b7dedbda40b?q=&sort=_id+asc&fields=OSHPD_ID%2CFACILITY_NAME%2CLICENSE_NUM%2CFACILITY_LEVEL_DESC%2CDBA_ADDRESS1%2CDBA_CITY%2CDBA_ZIP_CODE%2CCOUNTY_CODE%2CCOUNTY_NAME%2CER_SERVICE_LEVEL_DESC%2CTOTAL_NUMBER_BEDS%2CFACILITY_STATUS_DESC%2CFACILITY_STATUS_DATE%2CLICENSE_TYPE_DESC%2CLICENSE_CATEGORY_DESC%2CLATITUDE%2CLONGITUDE&filters=%7B%22COUNTY_CODE%22%3A+%5B%2236%22%5D%7D&format=csv"
        download = requests.get(url)

        with open(out_file_name, 'w') as out_file:
            out_file.writelines(download.text)
            print(out_file_name)
        csv_item = my_gis.content.add({}, out_file_name)
        hospital_item = csv_item.publish()
display(hospital_item)
You do not have permissions to access this resource or perform this operation.
Trying from an alternative source...
hospitals_SB_county
Feature Layer Collection by arcgis_python
Last Modified: September 27, 2019
0 comments, 10 views

From the Feature Layer item, we would derive the url of its first Feature Service, and create into a Feature Layer class object, which will then be used to make queries.

In [4]:
hospital_fl = FeatureLayer(hospital_item.url + "/0")
hospital_fl
Out[4]:
<FeatureLayer url:"https://services7.arcgis.com/JEwYeAy2cc8qOe3o/arcgis/rest/services/hospitals_SB_county/FeatureServer/0">
In [5]:
""" If you are using the exisiting layer from Esri_LivngAtlas, there is a "County" column in the dataset;
    or else, the feature layer collection published from the downloaded CSV file is already targetted at SB County.
"""
try:
    facilities = hospital_fl.query(where="County='SAN BERNARDINO' AND State='CA'", as_df=False)
except RuntimeError as re:
    """ when seeing 'Invalid field: County' parameter is invalid
    """
    print("Trying from an alternative approach...")
    facilities = hospital_fl.query(where="Dba_city='REDLANDS' or Dba_city='LOMA LINDA'", as_df=False)
display(facilities)
'Invalid field: County' parameter is invalid
Trying from an alternative approach...
<FeatureSet> 33 features

If the procedures all go well, you will be seeing a facilities object of a FeatureSet that contains 33 features. With such, we can draw the output FeatureSet with customized symbology in the map widget below.

In [25]:
map1 = my_gis.map('SAN BERNARDINO, CA', zoomlevel=12)
map1