Fixing Accessibility Scores That Are Wrong Near the Study Area Edge

Problem statement

Accessibility scores computed from the road network, facilities and population inside a study area are wrong near its edge. The nearest facility may be just over the line; the fastest road to an inside facility may leave the area; and a facility's catchment contains people who were never loaded. The symptom is a ring of apparently poorly served places around the boundary — and sometimes a ring of apparently excellent service, if the inputs were cut inconsistently.

Measured with supermarket accessibility for the 168,323 residents of Chittenden County, Vermont, against a calculation with the network, the 41 supermarkets and the population for 10 km around the county:

• Within 2 km of the county's land boundary, the county-only nearest-store time was 29.0% too long: 11.36 minutes against 8.81.
• The county-only two-step floating catchment (2SFCA) score there was 58.6% too low, and 14.9% of all residents had a score off by 10% or more.
• Adding the outside supermarkets but not the outside population made it 326.8% too high near the edge.
• A 2 km margin fixed nearest-store times to within 2.9% but left the 2SFCA score 19.8% low.
• Along the county's 53 km of boundary in Lake Champlain there was almost nothing to fix.

Quick answer

Extend all three inputs — network, facilities and population — beyond the study area by at least the reach of the largest catchment, compute, and keep only the study area's results:

truncated = with_margin(graph, node_xy, stores, blocks, county, residents, margin_m=0)
buffered = with_margin(graph, node_xy, stores, blocks, county, residents, margin_m=10_000)

with_margin (Example 1) keeps the road nodes, stores and population blocks within the margin, computes nearest-store time and 2SFCA, and returns scores for the study area's residents only. Compare the two before trusting either.

Confirm the edge pattern, then extend everything together

1. Confirm that the errors follow the boundary

Compute the scores twice — county-only, and with a generous margin — and compare them by distance to the boundary. For residents nearer the county's land boundary than its lake shore:

distance to land edge   nearest store: county-only vs 10 km margin   2SFCA: county-only vs 10 km margin
0-2 km                   11.36 vs 8.81 min  (+29.0%)                  0.584 vs 1.410  (-58.6%)
2-5 km                    7.75 vs 7.26 min                            0.702 vs 1.227
5-10 km                   4.93 vs 4.87 min                            1.004 vs 1.193

The nearest-store error faded within about 5 km. The 2SFCA error reached further, because a 15-minute catchment reaches further than the nearest store. County-wide, the population-weighted mean 2SFCA score was 1.782 county-only against 1.852 with the margin (Example 3 prints all of this).

2. Separate land edges from water edges

The county boundary is 201 km long: 148 km on land and 53 km in Lake Champlain. Across the lake there are no roads or shops within reach, so the margin changes nothing there: for residents 2–10 km from the shore, nearest-store times were identical with and without the margin, and 2SFCA scores differed by under 5%. Split the boundary before deciding where a margin is needed (Example 3).

3. Extend the network, the facilities and the population together

With a 10 km margin the road network grew from 6,498 to 9,818 nodes, the supermarkets from 30 to 41, and the population blocks from 2,241 to 3,452. Each matters for a different reason: the network for routes that leave the area, the facilities for nearer options across the line, and the population for the competition in each facility's catchment.

4. Do not extend only some of the inputs

Each partial fix was measured against the full 10 km margin, within 2 km of the land edge (Example 4):

inputs extended by 10 km            nearest store   2SFCA near edge   2SFCA county mean
none (county only)                     +29.0%          -58.6%             -3.8%
network only                           +27.6%          -54.9%             -3.8%
network + population                   +27.6%          -59.5%             -5.0%
network + facilities                    0.0%          +326.8%            +15.5%
network + facilities + population       0.0%             0.0%              0.0%

Adding the outside facilities fixed nearest-store times but, without the people who also use them, made those facilities look almost empty — so everyone near the boundary appeared to have them nearly to themselves.

5. Choose the margin from the catchment

A margin narrower than the catchment's reach can make things worse. With a 2 km margin the nearest-store error near the edge fell to 2.9% but the 2SFCA error stayed at −19.9%; with 5 km it was −29.1%, because the wider ring added population without adding any more supermarkets. Across all residents, 9.5% had a 2SFCA score off by 10% or more with a 2 km margin and 4.8% with 5 km. For a 15-minute drive catchment, 10 km was the widest margin the network data here allowed, so treat it as the reference rather than proof that nothing changes further out.

6. Report only the study area, and say what margin was used

Scores for places inside the margin are themselves edge-affected. Keep them for the calculation and drop them from the results, and state the margin with the map.

Scores with and without a margin

The data used on this page

Chittenden County, Vermont, from the 2020 Census: the TIGER/Line 2020 blocks for Vermont, which carry each block's POP20, dissolved to give the county (land and lake), with populated blocks reduced to a representative point; and the county's TIGER area-water layer for Lake Champlain. Supermarkets are OpenStreetMap shop=supermarket within 10 km of the county, and the road network is OpenStreetMap's drive network for the county plus 10.5 km, with OSMnx's imputed free-flow travel times, held as a sparse matrix of seconds. Only Vermont blocks are loaded: the New York shore has no road link inside the margin. OpenStreetMap comes from the maps.mail.ru Overpass mirror, as saved on 9 October 2026.

import io
import os
import requests
import urllib.request
import zipfile

import geopandas as gpd
import numpy as np
import osmnx as ox
import scipy.sparse as sp
import shapely

# OpenStreetMap as Overpass returned it on 9 October 2026, saved as osmnx's own cache so the counts
# below do not drift as OSM is edited; delete the cache folder to query live instead
if not os.path.exists("cache"):
    reply = requests.get("https://www.spatialworkflow.io/data/osm-cache-accessibility-edge-effects-2026-10-09.zip", timeout=300)
    reply.raise_for_status()
    with zipfile.ZipFile(io.BytesIO(reply.content)) as archive:
        archive.extractall("cache")

ox.settings.overpass_url = "https://maps.mail.ru/osm/tools/overpass/api"
ox.settings.overpass_rate_limit = False
ox.settings.http_user_agent = "spatialworkflow.io guide examples (https://www.spatialworkflow.io/)"
CRS = 32145                                               # NAD83 / Vermont, metres

TIGER = "https://www2.census.gov/geo/tiger/TIGER2020/"
for name, url in [("tl_2020_50_tabblock20.zip", TIGER + "TABBLOCK20/tl_2020_50_tabblock20.zip"),
                  ("tl_2020_50007_areawater.zip", TIGER + "AREAWATER/tl_2020_50007_areawater.zip")]:
    if not os.path.exists(name):
        urllib.request.urlretrieve(url, name)

vt_blocks = gpd.read_file("tl_2020_50_tabblock20.zip").to_crs(CRS)
county = vt_blocks[vt_blocks.COUNTYFP20 == "007"].union_all()           # Chittenden, land and lake
populated = vt_blocks[vt_blocks.POP20 > 0].copy()
populated["geometry"] = populated.geometry.representative_point()
residents = populated[populated.COUNTYFP20 == "007"].reset_index(drop=True)
blocks = populated[shapely.distance(populated.geometry.values, county) <= 10_000].reset_index(drop=True)
water_polygons = gpd.read_file("tl_2020_50007_areawater.zip").to_crs(CRS)

reach = gpd.GeoSeries([county.buffer(10_500)], crs=CRS).to_crs(4326).iloc[0]
stores = ox.features_from_polygon(reach, tags={"shop": "supermarket"}).to_crs(CRS)
stores["geometry"] = stores.geometry.representative_point()
stores = stores[shapely.distance(stores.geometry.values, county) <= 10_000].reset_index()

G = ox.add_edge_travel_times(ox.add_edge_speeds(ox.graph_from_polygon(reach, network_type="drive")))
G = ox.project_graph(G, to_crs=CRS)
index = {n: i for i, n in enumerate(G.nodes)}
node_xy = np.array([(d["x"], d["y"]) for _, d in G.nodes(data=True)])
seconds = {}
for u, v, d in G.edges(data=True):                      # shortest of any parallel edges
    seconds[(index[u], index[v])] = min(seconds.get((index[u], index[v]), np.inf), d["travel_time"])
rows, cols = zip(*seconds)
graph = sp.csr_array((list(seconds.values()), (rows, cols)), shape=(len(index), len(index)))
print(f"{residents.POP20.sum():,} residents in {len(residents):,} blocks; {len(blocks):,} populated blocks "
      f"within 10 km; {len(stores)} supermarkets; {len(index):,} road nodes")

168,323 residents in 2,241 blocks; 3,452 populated blocks within 10 km; 41 supermarkets; 10,055 road nodes

Example 1 — scores with a margin

import numpy as np
import shapely
from scipy.sparse.csgraph import connected_components, dijkstra
from scipy.spatial import cKDTree


def accessibility(graph, node_xy, stores, demand, residents, catchment=15):
    """Nearest-store minutes and 2SFCA (stores per 10,000) for residents; graph weights in seconds."""
    tree = cKDTree(node_xy)
    s = tree.query(np.column_stack([stores.geometry.x, stores.geometry.y]))[1]
    d = tree.query(np.column_stack([demand.geometry.x, demand.geometry.y]))[1]
    r = tree.query(np.column_stack([residents.geometry.x, residents.geometry.y]))[1]
    unique, inverse = np.unique(s, return_inverse=True)
    to_store = dijkstra(graph.T.tocsr(), directed=True, indices=unique)[inverse] / 60   # stores x nodes
    load = (to_store[:, d] <= catchment) @ demand.POP20.to_numpy(float)
    ratio = np.divide(1.0, load, out=np.zeros(len(load)), where=load > 0)
    reach = to_store[:, r].T
    return reach.min(axis=1), ((reach <= catchment) * ratio).sum(axis=1) * 10_000


def with_margin(graph, node_xy, stores, blocks, study_area, residents, margin_m, *,
                stores_m=None, blocks_m=None, verbose=False, **kwargs):
    """Keep network nodes, stores and population blocks within a margin of the study area.

    stores_m and blocks_m default to margin_m; set them apart to extend only some inputs."""
    stores_m = margin_m if stores_m is None else stores_m
    blocks_m = margin_m if blocks_m is None else blocks_m
    near = shapely.distance(shapely.points(node_xy), study_area) <= margin_m
    idx = np.flatnonzero(near)
    _, label = connected_components(graph[idx][:, idx], directed=True, connection="strong")
    keep = idx[label == np.bincount(label).argmax()]

    def inside(gdf, m):
        return gdf[shapely.distance(gdf.geometry.values, study_area) <= m]

    s, b = inside(stores, stores_m), inside(blocks, blocks_m)
    if verbose:
        print(f"  margin {margin_m / 1000:g} km: {len(keep):,} road nodes, {len(s)} supermarkets, "
              f"{len(b):,} populated blocks")
    return accessibility(graph[keep][:, keep], node_xy[keep], s, b, residents, **kwargs)

graph is a sparse matrix of road travel times in seconds, node_xy the node coordinates, stores and blocks point layers covering the study area and the widest margin, and residents the study area's own blocks. Keeping the largest strongly connected component stops a trimmed network from stranding homes on fragments.

Example 2 — split the boundary into land and water

def split_boundary(study_area, water, tolerance=300):
    """Boundary pieces away from water (land edge) and along it (shore)."""
    wet = water.buffer(tolerance)
    return study_area.boundary.difference(wet), study_area.boundary.intersection(wet)


lake = water_polygons[water_polygons.area > 50e6].union_all()
land_edge, shore = split_boundary(county, lake)
print(f"boundary {county.boundary.length / 1000:.0f} km: land {land_edge.length / 1000:.0f} km, "
      f"shore {shore.length / 1000:.0f} km")

boundary 201 km: land 148 km, shore 53 km

The 300 m tolerance absorbs the small gap between a boundary drawn in the lake and the mapped water polygon; the 50 km² filter keeps only Lake Champlain.

Example 3 — errors by distance to the edge

def edge_bias(residents, edge, other_edge, truncated, buffered, bands=(0, 2, 5, 10)):
    km = shapely.distance(residents.geometry.values, edge) / 1000
    nearer = km <= shapely.distance(residents.geometry.values, other_edge) / 1000
    w = residents.POP20.to_numpy(float)
    for lo, hi in zip(bands[:-1], bands[1:]):
        m = nearer & (km >= lo) & (km < hi)
        a, b = np.average(truncated[m], weights=w[m]), np.average(buffered[m], weights=w[m])
        print(f"{lo}-{hi} km: {m.sum()} blocks, truncated {a:.3f} vs buffered {b:.3f} ({a / b - 1:+.1%})")


near_0, tsfca_0 = with_margin(graph, node_xy, stores, blocks, county, residents, 0, verbose=True)
near_10, tsfca_10 = with_margin(graph, node_xy, stores, blocks, county, residents, 10_000, verbose=True)
print("land edge — nearest-store minutes, then 2SFCA")
edge_bias(residents, land_edge, shore, near_0, near_10)
edge_bias(residents, land_edge, shore, tsfca_0, tsfca_10)
print("lake shore — nearest-store minutes, then 2SFCA")
edge_bias(residents, shore, land_edge, near_0, near_10)
edge_bias(residents, shore, land_edge, tsfca_0, tsfca_10)
w = residents.POP20.to_numpy(float)
print(f"county mean 2SFCA: truncated {np.average(tsfca_0, weights=w):.3f} "
      f"vs buffered {np.average(tsfca_10, weights=w):.3f}")

  margin 0 km: 6,498 road nodes, 30 supermarkets, 2,241 populated blocks
  margin 10 km: 9,818 road nodes, 41 supermarkets, 3,452 populated blocks
land edge — nearest-store minutes, then 2SFCA
0-2 km: 161 blocks, truncated 11.358 vs buffered 8.826 (+28.7%)
2-5 km: 218 blocks, truncated 7.746 vs buffered 7.254 (+6.8%)
5-10 km: 328 blocks, truncated 4.929 vs buffered 4.873 (+1.2%)
0-2 km: 161 blocks, truncated 0.584 vs buffered 1.411 (-58.6%)
2-5 km: 218 blocks, truncated 0.702 vs buffered 1.228 (-42.8%)
5-10 km: 328 blocks, truncated 1.004 vs buffered 1.193 (-15.8%)
lake shore — nearest-store minutes, then 2SFCA
0-2 km: 5 blocks, truncated 10.557 vs buffered 10.557 (+0.0%)
2-5 km: 124 blocks, truncated 3.408 vs buffered 3.408 (+0.0%)
5-10 km: 700 blocks, truncated 2.432 vs buffered 2.432 (+0.0%)
0-2 km: 5 blocks, truncated 0.574 vs buffered 0.612 (-6.2%)
2-5 km: 124 blocks, truncated 1.118 vs buffered 1.161 (-3.7%)
5-10 km: 700 blocks, truncated 2.066 vs buffered 2.076 (-0.5%)
county mean 2SFCA: truncated 1.782 vs buffered 1.852

Under each heading the first three lines are nearest-store minutes and the next three 2SFCA scores; land_edge and shore come from Example 2. Along the lake shore the nearest-store times did not change at all.

Example 4 — partial fixes and narrower margins

km = shapely.distance(residents.geometry.values, land_edge) / 1000
edge = (km < 2) & (km <= shapely.distance(residents.geometry.values, shore) / 1000)
everyone = np.ones(len(residents), bool)


def err(a, b, m):
    return np.average(a[m], weights=w[m]) / np.average(b[m], weights=w[m]) - 1


scenarios = {
    "none (county only)": dict(margin_m=0),
    "network only": dict(margin_m=10_000, stores_m=0, blocks_m=0),
    "network + population": dict(margin_m=10_000, stores_m=0),
    "network + facilities": dict(margin_m=10_000, blocks_m=0),
    "all three, 2 km": dict(margin_m=2_000),
    "all three, 5 km": dict(margin_m=5_000),
    "all three, 10 km": dict(margin_m=10_000),
}
print(f"{'inputs extended':22} {'nearest store':>13} {'2SFCA near edge':>16} {'2SFCA mean':>11} {'off by 10%+':>12}")
for label, kw in scenarios.items():
    near, tsfca = with_margin(graph, node_xy, stores, blocks, county, residents, **kw)
    with np.errstate(divide="ignore", invalid="ignore"):
        off = np.average(np.abs(tsfca / tsfca_10 - 1) >= 0.1, weights=w)
    print(f"{label:22} {err(near, near_10, edge):+13.1%} {err(tsfca, tsfca_10, edge):+16.1%} "
          f"{err(tsfca, tsfca_10, everyone):+11.1%} {off:12.1%}")

inputs extended        nearest store  2SFCA near edge  2SFCA mean  off by 10%+
none (county only)            +28.7%           -58.6%       -3.8%        14.9%
network only                  +27.3%           -55.0%       -3.8%        12.6%
network + population          +27.3%           -59.5%       -5.0%        10.5%
network + facilities           +0.0%          +326.5%      +15.5%        12.6%
all three, 2 km                +2.9%           -19.9%       +0.7%         9.5%
all three, 5 km                +2.9%           -29.1%       -1.2%         4.8%
all three, 10 km               +0.0%            +0.0%       +0.0%         0.0%

Every row is measured against the full 10 km margin, for the residents within 2 km of the land edge and, in the last two columns, for everyone. The share "off by 10%+" is population-weighted.

Why errors stay near the boundary, and catchments reach further

Why nearest-facility errors stay near the boundary

A home's nearest facility can only be outside the study area if the boundary is closer than every inside facility. A few kilometres in, there is almost always a nearer facility on the inside, so the error stops. Here it was 29% within 2 km of the land boundary and about 1% at 5–10 km.

Why catchment measures are affected further in

2SFCA divides each facility's supply by all the people within its catchment. A facility 10 km inside the boundary still has part of its 15-minute catchment across the line, and the people there were never counted — and a home's own catchment can include facilities across the line. Both errors extend as far as the catchment reaches.

Why partial fixes can be worse than none

Supply and demand have to be cut consistently. Facilities without their outside users look underused, so their ratios, and the scores of every home that reaches them, rise. People without their outside facilities look underserved. A margin narrower than the catchment includes some of each in proportions that depend on where the facilities happen to be.

When truncation is right

If a service is only available to residents — a county clinic, a school district — people across the line cannot use it and outside facilities are not options for residents. Then the right inputs are residents and eligible facilities, with the full road network for routing. Decide which case applies before adding a margin.

Mismatched sources, water and trimmed networks

• Borders with different data sources need the same population and facility definitions on both sides.
• Water is not always a barrier. Ferries and bridges connect shores; check the network.
• International and state borders can restrict use even where roads connect.
• Gravity measures are affected like 2SFCA, since distant facilities still contribute.
• Trimmed networks create stranded nodes. Keep the largest strongly connected component after trimming.
• Maps of margins need care. Do not publish the margin's own scores as if they were reliable.
• The margin is a parameter. Report it, and show that a wider margin does not change the result where data allow.

Internal links

• How to measure access to services with the two-step floating catchment method — the catchment measure tested here
• Accessibility measures explained: nearest, cumulative and gravity — which measures reach how far
• How to count the population within reach of each facility — demand beyond the boundary
• How to measure distance to the nearest facility for every home — nearest-facility times
• How to find the areas nobody can reach within a travel time — gaps that may be edge artefacts
• NetworkX Says No Path Exists Between Two Nodes
• How to build Voronoi service areas around facilities in Python — the same boundary problem for assignment
• How to measure network distance for many origin–destination pairs — the travel-time matrix

FAQ

What is an edge effect in accessibility analysis?

Error near a study area's boundary caused by leaving out the roads, facilities and people beyond it. Within 2 km of Chittenden County's land boundary, county-only nearest-supermarket times were 28.7% too long.

How wide should the buffer around my study area be?

At least as wide as the largest catchment or travel threshold reaches. A 2 km margin fixed nearest-store times but left 2SFCA scores near the edge 19.9% too low; the 10 km reference margin was needed for 15-minute catchments.

Do I need to include population outside the study area?

Yes, for any measure that shares facilities between people. Adding outside facilities without outside population made 2SFCA scores near the edge 326.5% too high.

Why are my 2SFCA scores wrong further from the edge than nearest-facility times?

Because a catchment reaches further than the nearest facility. County-only 2SFCA scores were still 16% low 5 to 10 km from the land boundary, where nearest-store times were almost right.

Does a lake or coastline cause edge effects?

Much less, when nothing is reachable across it. Along the county's lake shore, nearest-store times were unchanged and 2SFCA scores moved by at most 6%.

Should I always add a buffer?

Not if only residents may use the service. Then outside facilities are not real options, and a margin would overstate access.