Fixing a Gravity Model That Sends Everyone to the Biggest Store
Problem statement
You fit a Huff model to a set of stores and map each neighbourhood's most likely store. One store — the largest — wins almost everywhere, including in towns with their own supermarket a few minutes away. Or its predicted market share looks far too large. Nothing raised an error.
Measured with 39 supermarkets around Chittenden County, Vermont, building footprints as attractiveness and network drive times from 2,241 populated census blocks, the largest store is genuinely the nearest store for 9.3% of residents. Five common mistakes each inflate it:
setting homes whose likeliest store is the largest
power decay, β = 2, minutes (reference) 12.7%
no distance decay, β = 0 100.0%
β = 0.5 27.4%
exponential λ = 0.2 per minute, times in hours 59.1%
power decay on (hours + 1) 38.7%
And two inflate its market share instead: raising attractiveness to the power 2 took its share from 8.5% to 13.2%, and filling missing footprints with 1 m² took it to 11.2%.
Quick answer
Check the parameters against the units and the data, in this order:
import numpy as np
A = footprints.fillna(footprints.median()).to_numpy() # a stated fill rule, not an accident
T = np.maximum(T_minutes, 0.5) # no zero travel times
P = A[None, :] ** 1.0 * T ** -2.0 # power decay: unit-free
P /= P.sum(axis=1, keepdims=True)
largest = A.argmax()
print("likeliest is largest:", pop[P.argmax(axis=1) == largest].sum() / pop.sum())
print("nearest is largest: ", pop[T.argmin(axis=1) == largest].sum() / pop.sum())If the first figure is far above the second, distance decay is too weak or in the wrong units. If the largest store's share is far above its share of floor area, α is too high or missing attractiveness is filled too low.
Step-by-step solution
1. Measure the symptom against the nearest-store baseline
"Everyone goes to the biggest store" needs a number. Compute the share of population whose likeliest store is the largest, and compare it with the share whose nearest store is the largest — 9.3% here. A well-specified model with reasonable decay lands somewhat above the nearest-store share, because a large store does pull custom past smaller ones: 12.7% at β = 2.
2. Check that distance decay is switched on
With β = 0, travel time drops out: every home's probabilities are proportional to floor area, so the largest store is everyone's likeliest store — 100% — and its share equals its share of floor area, 6.7%. This happens when β is read from a configuration as zero, or when a variable meant to hold travel time holds something constant.
3. Check that decay is strong enough
Small exponents leave distance too weak to overcome size. At β = 0.5 the largest store was the likeliest for 27.4% of residents. Unless you have calibrated β, test the range you consider plausible and show how the map changes across it:
β = 0.5 27.4%
β = 1.0 14.4%
β = 2.0 12.7%
β = 5.0 11.3%
4. Check the units of an exponential or offset decay
The power form, T^−β, gives the same probabilities in minutes, seconds or hours. The exponential form does not: λ is per unit of time.
exponential λ = 0.2 per minute likeliest is largest 18.4%
same λ applied to times in seconds 9.3% (decay 60× too strong)
same λ applied to times in hours 59.1% (decay 60× too weak)
An offset such as (T + 1)^−β breaks unit invariance too: in minutes it gave 12.9%, in hours 38.7%, because one hour is not one minute. Convert times to the unit the parameter was fitted in, or use the pure power form.
5. Check the attractiveness exponent
α = 1 treats a store twice the size as twice as attractive. Larger α lets large stores capture custom far beyond their surroundings:
α = 0.5 largest store's share 6.2% top three 18.7%
α = 1.0 8.5% 21.7%
α = 2.0 13.2% 32.4%
Values above 1 need evidence. Without calibration, α = 1 is the neutral choice.
6. Check how missing attractiveness was filled
Fifteen of the 39 supermarkets had no building footprint in OpenStreetMap. The fill rule changes who wins:
missing footprints filled with the median largest store's share 8.5%
missing footprints filled with the maximum 5.9%
missing footprints filled with 1 m² 11.2%
Filling with a tiny value — or with zero, which removes the stores entirely — hands their customers to the stores that do have footprints, including the largest. Use a real measure for every store if at all possible.
7. Floor zero travel times instead of dropping rows
Seven origin–store pairs had a travel time of zero. Without a floor, those rows became NaN in the power form; if they are later dropped or filled with zeros, the result is biased in whichever direction the rows happened to be. Floor times at a small positive value in the same unit.
8. Recalibrate if you have observed data
If loyalty-card, survey or turnover data exist, fit β and α to them. Two uncalibrated parameters can reproduce almost any market story; calibrated ones are a model.
Code examples
Example 1 — a diagnostic for the largest-store effect
import numpy as np
def largest_store_diagnostic(A, T, population, alpha=1.0, beta=2.0, floor=0.5):
"""Compare the model's pull of the largest store with the nearest-store baseline."""
A = np.asarray(A, float)
Tf = np.maximum(np.asarray(T, float), floor)
U = A[None, :] ** alpha * Tf ** -beta
P = U / U.sum(axis=1, keepdims=True)
largest = int(np.argmax(A))
w = population / population.sum()
report = {
"floor_area_share": A[largest] / A.sum(),
"market_share": float((P[:, largest] * w).sum()),
"likeliest_is_largest": float(w[P.argmax(axis=1) == largest].sum()),
"nearest_is_largest": float(w[Tf.argmin(axis=1) == largest].sum()),
}
print(" | ".join(f"{k} {v:.1%}" for k, v in report.items()))
return reportfloor_area_share 6.7% | market_share 8.5% | likeliest_is_largest 12.7% | nearest_is_largest 9.3%
These four numbers are the sanity check. A likeliest-is-largest share several times the nearest-is-largest share, or a market share several times the floor-area share, points to one of the mistakes above.
Example 2 — sweeping the decay parameter
def decay_sweep(A, T, population, betas=(0.0, 0.5, 1.0, 2.0, 5.0)):
for beta in betas:
r = largest_store_diagnostic(A, T, population, beta=beta)floor_area_share 6.7% | market_share 6.7% | likeliest_is_largest 100.0% | nearest_is_largest 9.3%
floor_area_share 6.7% | market_share 7.8% | likeliest_is_largest 27.4% | nearest_is_largest 9.3%
floor_area_share 6.7% | market_share 8.4% | likeliest_is_largest 14.4% | nearest_is_largest 9.3%
floor_area_share 6.7% | market_share 8.5% | likeliest_is_largest 12.7% | nearest_is_largest 9.3%
floor_area_share 6.7% | market_share 9.3% | likeliest_is_largest 11.3% | nearest_is_largest 9.3%
The rows are β = 0, 0.5, 1, 2 and 5. A curve that flattens between 1 and 2, as here, says that the exact value in that range matters little for the household map; a curve that is still falling steeply says it matters a lot.
Example 3 — refusing unit-dependent decay without a unit
def distance_decay(T, form="power", beta=None, lam=None, unit=None):
"""Distance decay that refuses to guess units for the forms that depend on them."""
T = np.asarray(T, float)
if form == "power":
return T ** -beta
if form == "exponential":
if unit not in ("s", "min", "h"):
raise ValueError("exponential decay needs the unit lambda was fitted in: 's', 'min' or 'h'")
return np.exp(-lam * T)
raise ValueError(f"unknown form {form!r}")Making the unit a required argument does not convert anything; it makes whoever calls the function state what λ means, which is where the 59.1% error came from.
Explanation
Why weak decay favours size
In the Huff model a store's pull is its attractiveness divided by a power of distance. When the power is small, distance barely discounts a store, so the store with the largest numerator wins every household's comparison. As β grows, nearby stores' denominators become much smaller than distant stores', and local stores win their own neighbourhoods.
Why exponential decay is fragile
exp(−λT) compares λ × T with one. Change the unit of T by a factor of 60 and the product changes by 60, so the same λ describes decay that is 60 times stronger or weaker. The power form T^−β rescales every store's utility by the same factor, which cancels when probabilities are normalised. That is the whole reason the power form survived unit mistakes here and the exponential form did not.
Why missing values move shares
Attractiveness enters every household's normalisation. Giving 15 stores a near-zero attractiveness removes them from the competition for everyone, and their customers spread to the remaining stores in proportion to those stores' pull. The largest store takes the largest part of that windfall.
Why the share and the map can disagree
A store's share sums probabilities; the likeliest-store map takes each household's maximum. A store can be every household's second choice and still earn a large share, or be the likeliest for many households with modest probabilities each. Checking only one of the two lets the other go wrong unseen.
Edge cases or notes
- Zero attractiveness removes a store. A store with
A = 0has zero probability for every household. - Distances in degrees behave like a unit error. Project coordinates or use network times.
- A single dominant store can be real. A regional hypermarket may win widely; compare with observed data before "fixing" it.
- Stores outside the study area must be included or in-area stores absorb their customers.
- Constrained models behave differently. If store capacity is enforced, the largest store cannot absorb unlimited demand.
- Aggregated origins hide the effect. Tract-level origins average away the household map; use small units.
- Report α, β, the decay form, the time unit and the fill rule with any published result.
Internal links
- Gravity and Huff models explained: predicting where people go — what the parameters mean
- How to estimate market share with a Huff model in Python — the complete workflow
- Accessibility measures explained: nearest, cumulative and gravity — the same decay in access scores
- How to measure distance to the nearest facility for every home — the nearest-store baseline
- Fixing accessibility scores that are wrong near the study area edge — missing stores beyond the boundary
- Catchment areas explained: buffers, isochrones and Voronoi compared — deterministic catchments for comparison
- How to handle missing and null values in spatial datasets — filling attractiveness honestly
- Location analysis explained: catchments, accessibility and site selection — where market share fits
FAQ
Why does my Huff model assign everyone to the largest store?
Distance decay is too weak or not working. With β = 0 the largest store was every household's likeliest store; at β = 0.5 it was for 27.4%; at β = 2 for 12.7%, close to the 9.3% for whom it is actually nearest.
Does the unit of travel time matter?
For power decay, no: minutes and seconds gave identical probabilities. For exponential decay, yes: a per-minute λ applied to hours made the largest store the likeliest for 59.1% of residents.
What attractiveness exponent should I use?
α = 1 unless data support something else. At α = 2 the largest store's share rose from 8.5% to 13.2%, and the top three stores took a third of the market.
How should I handle stores with no size data?
Fill them with a stated, defensible value, or find a real measure. Filling 15 missing footprints with 1 m² raised the largest store's share from 8.5% to 11.2%.
Why do some rows come out as NaN?
A travel time of zero raised to a negative power is infinite. Seven origin–store pairs had zero time; floor times at a small positive value.
How do I know the model is reasonable?
Compare the share of people whose likeliest store is the largest with the share for whom it is the nearest, and the largest store's market share with its share of floor area. Better still, calibrate against observed shopping data.