Diagnosed by Claude (Sonnet 5) while debugging _geographic-key-concepts.qmd in hsma_lambda_geographic on 2026-08-05.
Summary
solution.get_hotspots(what="travel_equity") can raise an opaque esda/scipy error when the travel matrix passed to add_travel_matrix() doesn't cover every location in demand_data / region_geometry_layer. The real problem — mismatched ID coverage, usually from a stale or wrong-vintage travel matrix — is invisible until it surfaces three call frames deep as a scipy.sparse matmul shape error.
Repro
problem = SiteProblem()
problem.add_demand("demand.csv", demand_col="...", location_id_col="LSOA 2021 Name")
problem.add_sites(sites_gdf, candidate_id_col="Facility_Name")
problem.add_travel_matrix("travel_matrix.csv", source_col="from_id", unit="minutes")
problem.add_equity_data("equity.csv", equity_col="...", common_col="LSOA name (2021)", ...)
problem.add_region_geometry_layer(region_gdf, common_col="LSOA 2021 Name")
solution = problem.solve(p=15)
solution.get_hotspots(what="travel_equity")
ValueError: matmul: dimension mismatch with signature (n,k=729),(k=683,1?)->(n,1?)
Root cause
In our case: demand_data, equity_data, and region_geometry_layer all agreed on 729 LSOAs. The travel matrix only had 707 from_id values, and 24 of those used a different sub-area LSOA vintage (e.g. Torbay 013F, Plymouth 029A–029F) that doesn't exist in the 2021 boundary file we were using elsewhere — a boundary-vintage/provenance mismatch, not a whitespace/casing issue. Overlap between the travel matrix and the official geometry: 683 LSOAs.
Two independent code paths in lokigi/mixins/site_eda.py then disagree on how many rows survive:
_prepare_analysis_dataframe(what="travel_equity") (~lines 610–681) builds its result via demand_data ⋈ travel_df ⋈ equity_data, both how="inner" — silently drops the 46 uncovered LSOAs, giving a 683-row result.
get_hotspots() (~lines 292–376) builds the spatial weights matrix w by filtering region_geometry_layer only against demand_data ids (never against travel/equity coverage) — stays at 729 rows.
esda.moran.Moran_Local(result[df_col], w) then fails on the 683-vs-729 mismatch.
add_travel_matrix() itself does no cross-validation against IDs already registered via add_demand()/add_region_geometry_layer() — it only checks for source_col presence and NaN values (_reject_missing_travel_values), not for which IDs are present.
Suggested fix
Primarily: warn, don't silently drop. When a travel matrix (or equity dataset) is added whose ID coverage doesn't fully match IDs already loaded via add_demand() / add_region_geometry_layer() (or vice versa, if added in the other order), emit a UserWarning naming:
- how many IDs are missing/uncovered,
- a few example IDs,
- a nudge that this commonly means a stale/wrong-vintage travel matrix or boundary file (e.g. mismatched LSOA boundary years).
This should fire at add_travel_matrix() / add_equity_data() / add_region_geometry_layer() time — whichever call completes the mismatch — rather than only surfacing when a downstream method happens to need full coverage.
Secondary/defensive: get_hotspots()'s geometry filter (site_eda.py:361–366) should align to the same ID set actually used to build result in _prepare_analysis_dataframe, not just demand_data. That way, even with the warning suppressed or missed, get_hotspots() degrades to the
Diagnosed by Claude (Sonnet 5) while debugging
_geographic-key-concepts.qmdinhsma_lambda_geographicon 2026-08-05.Summary
solution.get_hotspots(what="travel_equity")can raise an opaqueesda/scipyerror when the travel matrix passed toadd_travel_matrix()doesn't cover every location indemand_data/region_geometry_layer. The real problem — mismatched ID coverage, usually from a stale or wrong-vintage travel matrix — is invisible until it surfaces three call frames deep as ascipy.sparsematmul shape error.Repro
Root cause
In our case:
demand_data,equity_data, andregion_geometry_layerall agreed on 729 LSOAs. The travel matrix only had 707from_idvalues, and 24 of those used a different sub-area LSOA vintage (e.g.Torbay 013F,Plymouth 029A–029F) that doesn't exist in the 2021 boundary file we were using elsewhere — a boundary-vintage/provenance mismatch, not a whitespace/casing issue. Overlap between the travel matrix and the official geometry: 683 LSOAs.Two independent code paths in
lokigi/mixins/site_eda.pythen disagree on how many rows survive:_prepare_analysis_dataframe(what="travel_equity")(~lines 610–681) builds its result viademand_data ⋈ travel_df ⋈ equity_data, bothhow="inner"— silently drops the 46 uncovered LSOAs, giving a 683-rowresult.get_hotspots()(~lines 292–376) builds the spatial weights matrixwby filteringregion_geometry_layeronly againstdemand_dataids (never against travel/equity coverage) — stays at 729 rows.esda.moran.Moran_Local(result[df_col], w)then fails on the 683-vs-729 mismatch.add_travel_matrix()itself does no cross-validation against IDs already registered viaadd_demand()/add_region_geometry_layer()— it only checks forsource_colpresence and NaN values (_reject_missing_travel_values), not for which IDs are present.Suggested fix
Primarily: warn, don't silently drop. When a travel matrix (or equity dataset) is added whose ID coverage doesn't fully match IDs already loaded via
add_demand()/add_region_geometry_layer()(or vice versa, if added in the other order), emit aUserWarningnaming:This should fire at
add_travel_matrix()/add_equity_data()/add_region_geometry_layer()time — whichever call completes the mismatch — rather than only surfacing when a downstream method happens to need full coverage.Secondary/defensive:
get_hotspots()'s geometry filter (site_eda.py:361–366) should align to the same ID set actually used to buildresultin_prepare_analysis_dataframe, not justdemand_data. That way, even with the warning suppressed or missed,get_hotspots()degrades to the