How to Join a Census Table to Its Boundaries Without Losing Rows
Problem statement
Attaching a census table to its polygons is one line of pandas. The difficulty is everything that line does not tell you.
Measured with the 2023 ACS tract population table (85,381 rows) and the 2023 cartographic tract boundaries (85,186 polygons), a join that is entirely correct still leaves:
- 141 polygons with no data โ 126 in island areas the ACS does not survey, and 15 in New York.
- 336 table rows with no polygon โ water-only and empty tracts the boundary file omits.
An inner join hides both groups; a left join from the boundaries hides the second. Neither says why. And the same line, run on inputs with a fault, reports nothing either:
- Identifiers that lost their leading zeros: 15,972 table rows unmatched and 64,740,957 people โ 19.3% of the population โ missing from the map.
- Boundaries from 2020 instead of 2023: 1,214 table rows unmatched, 884 of them in Connecticut, 3,598,348 people unmapped.
The goal of this guide is a join in which every unmatched row has a known reason.
Quick answer
Merge from the boundaries, keep both sides, and check the key is unique:
import pandas as pd
import geopandas as gpd
tracts = gpd.read_file("cb_2023_us_tract_500k.zip")
table = pd.read_csv("acsdt5y2023-b01003.dat", sep="|", dtype=str)
table = table[table.GEO_ID.str.startswith("1400000US")].copy()
table["GEOID"] = table.GEO_ID.str[9:]
table["population"] = pd.to_numeric(table.B01003_E001)
joined = tracts.merge(table[["GEOID", "population"]], on="GEOID",
how="outer", indicator=True, validate="one_to_one")
print(joined["_merge"].value_counts().to_dict())
{'both': 85045, 'right_only': 336, 'left_only': 141}
how="outer" shows the unmatched rows from both files, indicator=True labels them, and validate raises if either side has a duplicate key. Explain the left_only and right_only rows, then do the left join you actually map.
Step-by-step solution
1. Read both identifiers as text
Pass dtype=str to every CSV or summary-file read. The boundary shapefile stores GEOID as text and GeoPandas keeps it that way, so the table side is where damage enters. An integer key does at least fail loudly:
TypeError: table.GEOID is int64; read ids as text
That message comes from the guard in Example 2 below. Converting the integers back to strings without padding silences it and loses 18.5% of the boundaries instead.
2. Derive the key the right way
A table's GEO_ID is the boundary's GEOID with a nine-character prefix. Three ways to strip it all matched 85,045 tracts; an off-by-two slice matched none:
table["GEOID"] = table.GEO_ID.str[9:] # 85,045 matches
table["GEOID"] = table.GEO_ID.str.replace("1400000US", "") # 85,045 matches
table["GEOID"] = table.GEO_ID.str[7:] # 0 matches
The 2023 boundary files also carry GEOIDFQ, identical to GEO_ID, so no key needs deriving at all:
direct = tracts.merge(table, left_on="GEOIDFQ", right_on="GEO_ID", how="inner")
print(len(direct))
85045
3. Merge from the GeoDataFrame
The left operand decides the result type:
print(type(tracts.merge(table, on="GEOID")).__name__,
type(table.merge(tracts, on="GEOID")).__name__)
GeoDataFrame DataFrame
A merge written from the table side returns a plain DataFrame, and .plot() and .to_file() stop working. If you must merge that way, rebuild it: gpd.GeoDataFrame(back, geometry="geometry", crs=tracts.crs).
4. Account for every unmatched row
left = joined[joined["_merge"] == "left_only"]
right = joined[joined["_merge"] == "right_only"]
print(left.STUSPS.value_counts().to_dict())
print(len(right), right.population.sum(), right.GEOID.str[5:7].eq("99").sum())
{'GU': 56, 'VI': 29, 'MP': 23, 'AS': 18, 'NY': 15}
336 0.0 335
The boundary file covers Guam, the US Virgin Islands, the Northern Mariana Islands and American Samoa; the ACS does not. The 336 table rows without polygons hold nobody, and 335 carry 99xxxx water-tract codes. That leaves 15 New York polygons โ 14 in Suffolk County and one in Ulster County โ that are genuinely unexplained by design, and worth a note in any analysis of those counties.
5. Strip whitespace and check for duplicates
A single trailing space on every key:
padded = table[["GEOID", "population"]].assign(GEOID=table.GEOID + " ")
print("plain merge on padded ids:", len(tracts.merge(padded, on="GEOID")))
plain merge on padded ids: 0
Duplicates are the opposite failure: three repeated table rows added three extra polygons to the result with no error, while validate="one_to_one" refused:
MergeError: Merge keys are not unique in right dataset; not a one-to-one merge
6. Check the boundary vintage
tracts2020 = gpd.read_file("cb_2020_us_tract_500k.zip", columns=["GEOID", "STUSPS"],
ignore_geometry=True)
v = tracts2020.merge(table[["GEOID", "population"]], on="GEOID", how="outer", indicator=True)
print(v["_merge"].value_counts().to_dict())
{'both': 84167, 'right_only': 1214, 'left_only': 1020}
Of the 1,214 unmatched table rows, 884 are in Connecticut, whose county-level codes changed in 2022 when planning regions replaced its counties in census geography; 879 Connecticut polygons are unmatched from the other side. Everything else about the join looks normal, which is why a vintage mismatch survives into published maps.
7. Do the final join as a left join, and record the numbers
The map needs every polygon, so the delivered join is a left join from the boundaries. Store the outer-merge counts beside it; Example 3 prints them in one call.
Code examples
Example 1 โ classify the unmatched rows
ISLAND_AREAS = {"60", "66", "69", "78"} # AS, GU, MP, VI
def classify_unmatched(merged, id_col="GEOID", width=11, value_col="population"):
"""Give every unmatched row of an outer merge a probable reason."""
out = merged[merged["_merge"] != "both"].copy()
ids = out[id_col].astype("string")
reason = pd.Series("check the vintage", index=out.index)
reason[ids.str.len() != width] = "wrong width: leading zeros or whitespace"
reason[(out["_merge"] == "left_only") & ids.str[:2].isin(ISLAND_AREAS)] = "island area: not in the ACS"
reason[(out["_merge"] == "right_only") & ids.str[5:7].eq("99")] = "water tract: no polygon"
reason[(out["_merge"] == "right_only") & ~ids.str[5:7].eq("99")
& out[value_col].eq(0)] = "empty tract: no polygon"
out["reason"] = reason
summary = out.groupby(["_merge", "reason"], observed=True).size()
print(summary.to_string())
return out
_merge reason
left_only check the vintage 15
island area: not in the ACS 126
right_only empty tract: no polygon 1
water tract: no polygon 335
The categories are ordered so the most specific reason wins. Anything still labelled "check the vintage" is the residue that needs a person.
Example 2 โ a join that refuses to return a broken result
def attach_table(boundaries, table, key="GEOID", max_unmatched_share=0.01):
"""Left-join a table to boundaries, refusing to return a quietly broken result."""
for frame, name in ((boundaries, "boundaries"), (table, "table")):
if frame[key].duplicated().any():
raise ValueError(f"duplicate {key} values in {name}")
if not pd.api.types.is_string_dtype(frame[key]):
raise TypeError(f"{name}.{key} is {frame[key].dtype}; read ids as text")
table = table.assign(**{key: table[key].str.strip()})
merged = boundaries.merge(table, on=key, how="left", indicator=True, validate="one_to_one")
unmatched = (merged["_merge"] == "left_only").mean()
print(f"{len(boundaries):,} boundaries, {len(table):,} table rows, "
f"{unmatched:.2%} of boundaries without data")
if unmatched > max_unmatched_share:
raise ValueError(f"{unmatched:.1%} of boundaries unmatched (limit {max_unmatched_share:.1%})")
return merged.drop(columns="_merge")
On the correct table, on ids that had lost their zeros and been turned back into strings, and on ids padded with a trailing space:
85,186 boundaries, 85,381 table rows, 0.17% of boundaries without data
85,186 boundaries, 85,381 table rows, 18.52% of boundaries without data
ValueError: 18.5% of boundaries unmatched (limit 1.0%)
85,186 boundaries, 85,381 table rows, 0.17% of boundaries without data
Set the tolerance from a known-good run of the same geography. For US tracts, 0.17% of boundaries is the floor that island areas and the New York tracts set.
Example 3 โ the report to keep beside any census join
def join_report(boundaries, table, key="GEOID", value_col="population"):
"""The numbers to record next to any census join."""
m = boundaries.merge(table, on=key, how="outer", indicator=True)
counts = m["_merge"].value_counts()
lost = m.loc[m["_merge"] == "right_only", value_col].sum()
total = table[value_col].sum()
report = {
"boundaries": len(boundaries), "table rows": len(table),
"matched": int(counts.get("both", 0)),
"boundaries without data": int(counts.get("left_only", 0)),
"table rows without boundary": int(counts.get("right_only", 0)),
f"{value_col} not mapped": f"{lost:,.0f} of {total:,.0f} ({lost / total:.3%})",
}
for k, val in report.items():
print(f"{k:28} {val}")
return report
boundaries 85186
table rows 85381
matched 85045
boundaries without data 141
table rows without boundary 336
population not mapped 0 of 335,559,225 (0.000%)
The whole report took 0.09 s for 85,000 tracts. The last line is the one that matters: the 336 unmatched rows hold nobody, so the correct join maps every person in the table. The same report on damaged ids ends with 64,740,957 of 335,559,225 (19.293%).
Explanation
Why the outer merge is the diagnostic join
An inner join returns only the rows both files agree on, so it cannot show you what was lost. A left join shows lost boundaries as empty values but drops table rows without trace. Only an outer merge returns both kinds of disagreement, and indicator=True labels them for free.
It costs nothing to run, and it turns "the join worked" from an assumption into a count.
Why a correct join still has unmatched rows
A table and a boundary file are made for different purposes. The table lists every tabulated area, including water-only tracts nobody lives in; the cartographic boundary file drops shapes with nothing to draw once clipped to the shoreline, and it includes territories surveyed by a different programme. Measured, those two facts explain 461 of the 477 unmatched rows, and none of them carry population.
Why the share of rows understates the damage
A tolerance on the share of unmatched rows is necessary but not sufficient. With a vintage mismatch only 1.4% of table rows failed to join, yet they held 3.6 million people and every tract in one state. With damaged identifiers 18.5% of boundaries failed and 19.3% of people went unmapped, because the affected states include California.
Failures concentrate in whole states, because identifiers and boundary changes are organised by state. That is why the report includes the population not mapped, and why Example 1 is worth running even when the share looks small.
Why the left operand matters
DataFrame.merge builds its result from the class of the object it was called on. Called on a GeoDataFrame, the result keeps the geometry column's special status and the CRS; called on a DataFrame, the geometry column comes back as an ordinary column of shapes. Nothing fails until the first spatial operation, often in a different notebook cell.
Edge cases or notes
- Unmatched is not the same as zero. An empty tract has population 0 in the table; a tract that failed to join has NaN. Do not fill one with the other.
- Never spatially join a census table. Tables have no geometry; the identifier is the join, and a spatial join of centroids to polygons only adds error.
GEOIDFQjoins without string work in the 2023 files; older boundary files may not carry it.- Island area codes are 60, 66, 69 and 78. Expect them unmatched against any ACS table.
- Water-only tracts are coded from 990000 upwards and are expected to lack polygons in cartographic files.
- Overlapping column names get suffixes. If both files carry
NAMEorALAND, pandas returns_xand_ycolumns; drop what you do not need before merging. - The same rules hold for UK tables. A Nomis LSOA table joins to ONS LSOA boundaries on
LSOA21CD; a 2011 lookup against 2021 boundaries fails in the way a 2020 US boundary file does. - Record the join report with the output. The unmatched counts are the provenance of the map.
Internal links
- Census identifiers explained: GEOIDs, codes and the leading zero problem โ the key this join depends on
- Fixing census joins that fail because leading zeros were dropped โ the 19.3% failure, fixed
- Fixing census data that does not match the boundary file year โ the Connecticut failure, in depth
- How to download census boundaries in Python โ getting the right boundary file
- How to download census tables from an API in Python โ getting the table
- How to join attribute data to a GeoDataFrame in Python โ the general attribute join
- GeoPandas merge returns NaN or no matches โ the general form of an empty join
- Attribute join or spatial join? โ why an identifier join is the right tool here
FAQ
Why do some census tracts have no data after a join?
In a correct 2023 join, 126 of the 141 polygons without data are in island areas the ACS does not survey. If far more are empty, suspect identifiers that lost their leading zeros or a boundary file from another year.
Should I use an inner join or a left join?
Run an outer merge first to see unmatched rows on both sides, then deliver a left join from the boundaries so every polygon stays on the map.
Why does my merged result no longer plot?
It was merged from the table side, which returns a plain DataFrame. Call merge on the GeoDataFrame, or rebuild the result with gpd.GeoDataFrame and the original CRS.
How do I join ACS GEO_ID to a shapefile GEOID?
Strip the nine-character prefix with GEO_ID.str[9:], or join GEO_ID to the GEOIDFQ column that the 2023 boundary files carry. Both matched 85,045 tracts.
How many unmatched rows are normal for US tracts?
With 2023 tables and 2023 cartographic boundaries, 141 boundaries and 336 table rows, holding no population. Treat anything above that as a fault to explain.
Why did Connecticut disappear from my map?
The boundary file is probably from 2020 or earlier. Connecticut's county-level codes changed in 2022, and 884 of its tract rows in the 2023 table did not match 2020 boundaries.