Site Suitability Explained: Constraints, Factors and Weights
Problem statement
"Where is the best place for it?" is the oldest question in applied GIS, and a suitability analysis is how a map answers it. Every version has the same three parts:
- Constraints โ yes/no rules that make a place ineligible: open water, protected land, slopes too steep to build on.
- Factors โ graded qualities where more or less is better: gentler slope, a shorter distance to a main road, more people nearby.
- Weights โ how much each factor counts relative to the others.
The method is simple arithmetic on a grid. Most of its mistakes come from treating those three parts as interchangeable. That was measured on a real 30 m grid of Chittenden County, Vermont: 1,704,541 land cells, two constraints (water, slope over 15ยฐ), three factors (slope, distance to a primary or secondary road, population density), and one question โ which places come out on top.
what went wrong effect on the answer
constraint left to the factors (water scored, not masked) 8.4% of the top 5% is open water
factors summed in raw units (degrees + metres + people/kmยฒ) 7 of the 10 best sites change
0.1 of weight moved between any pair of factors 0 or 1 of the 10 best sites change
The weights โ the part people argue about in meetings โ mattered least. Treating a constraint as a factor and forgetting to standardise mattered most.
Quick answer
Mask with the constraints, rescale each factor to 0โ1 with thresholds you can defend, then take a weighted sum:
import numpy as np
def decreasing(x, good, bad):
"""1 at or below `good`, 0 at or beyond `bad`, linear between."""
return np.clip((bad - x) / (bad - good), 0, 1)
def increasing(x, bad, good):
return np.clip((x - bad) / (good - bad), 0, 1)
factors = {
"slope": decreasing(slope_deg, good=0, bad=15),
"road": decreasing(dist_road_m, good=0, bad=3000),
"demand": increasing(people_per_km2, bad=0, good=500),
}
weights = {"slope": 1/3, "road": 1/3, "demand": 1/3}
allowed = in_county & ~water & (slope_deg <= 15) # constraints: a mask
score = sum(weights[k] * factors[k] for k in factors) # factors: a weighted sum
score = np.where(allowed, score, np.nan)The two constraints removed 14.9% of the county's land; the weighted sum then ranks the remaining 1,450,723 cells. The arithmetic is cheap โ the judgement in the thresholds and weights is not.
Step-by-step solution
1. Sort every criterion into a constraint or a factor
Ask of each criterion: is there a value beyond which the place is simply not a candidate? If so, it is a constraint. If a worse value can be made up for by a better value elsewhere, it is a factor.
Some criteria are both. Slope was used twice here: slopes over 15ยฐ were excluded outright, and below 15ยฐ gentler was better.
2. Express constraints as a mask, never as a low score
A constraint that only lowers a score can still be outscored. Leaving Lake Champlain and the county's other water bodies unmasked, and relying on the factors to rank them low, measured:
no mask: top 5% = 89,209 cells; in water blocks 7,482 (8.4%)
mean score: water cells 0.518 land cells 0.395
Water scored higher than land. On the slope factor it is perfect by definition (0.944 against 0.528 for land), and the 4 km around it is more densely settled than the county average (0.270 against 0.166). On road access the two were about level (0.339 against 0.364). Two of three factors said the lake was an excellent site.
3. Rescale each factor to a common 0โ1 range
Slope is measured in degrees, distance in metres and density in people per kmยฒ. Adding them as they are is adding numbers in different units, and whichever has the widest spread decides the answer (see Explanation).
Use thresholds that mean something: 0ยฐ is ideal and 15ยฐ is as bad as it gets; 0 m from a main road is ideal and 3 km is as bad as it gets. Thresholds make a factor score readable โ 0.5 means "halfway between ideal and unacceptable" โ and they do not move when the study area is extended.
4. Choose weights, and write down why
Direct weights are fine for three or four factors. For more, a pairwise comparison (the Analytic Hierarchy Process) forces one judgement at a time: is slope more important than road access, and by how much?
The judgements used here โ slope moderately more important than roads (3), slightly more than demand (2), demand slightly more than roads (2) โ give:
AHP weights slope/road/demand = [0.54 0.163 0.297] lambda_max 3.0092 CR 0.0079
A consistency ratio (CR) below 0.1 means the judgements agree with each other. A deliberately circular set โ slope beats roads, roads beat demand, demand beats slope โ gave CR 1.149, which is the method telling you the weights are incoherent.
5. Combine, then apply the mask
The weighted linear combination is score = ฮฃ weight ร factor, and the mask sets excluded cells to NoData. Keep NoData as NaN in a float array so it cannot be mistaken for a score of zero.
6. Select sites at the scale of the decision
A single 30 m cell is not a site. Summarise the score over the footprint you would actually build on โ here, a moving 1 km window with at least 90% eligible land โ and pick the best non-overlapping windows.
That also makes the result more stable than a cell-level top 5%, because one unusual cell cannot promote a whole site.
7. Test how much the answer depends on the weights
Move 0.1 of weight between each pair of factors, starting from equal weights, and count how many of the ten best sites survive:
+0.1 slope -0.1 road : 9/10 kept
+0.1 slope -0.1 demand: 9/10 kept
+0.1 road -0.1 slope : 10/10 kept
+0.1 road -0.1 demand: 9/10 kept
+0.1 demand -0.1 slope : 10/10 kept
+0.1 demand -0.1 road : 10/10 kept
AHP weights: 8/10 kept
raw units: 3/10 kept
If a site only appears under one weighting, report it as sensitive rather than as a recommendation.
Code examples
Example 1 โ factors with explicit thresholds and a constraint mask
import numpy as np
from scipy import ndimage
def distance_to(mask, cell_size):
"""Euclidean distance in map units from every cell to the nearest True cell."""
return ndimage.distance_transform_edt(~mask) * cell_size
def standardise(x, bad, good):
"""Linear rescale to 0..1; works for increasing (good > bad) and decreasing factors."""
out = (np.asarray(x, dtype="float32") - bad) / (good - bad)
return np.clip(out, 0, 1)
def weighted_overlay(factors, weights, allowed):
"""Weighted linear combination; excluded cells become NaN, not 0."""
total = sum(weights.values())
if not np.isclose(total, 1.0):
raise ValueError(f"weights sum to {total:.3f}, not 1")
missing = set(weights) ^ set(factors)
if missing:
raise KeyError(f"factor/weight names do not match: {sorted(missing)}")
score = sum(weights[name] * factors[name] for name in weights)
return np.where(allowed, score, np.nan).astype("float32")factors = {
"slope": standardise(slope_deg, bad=15, good=0),
"road": standardise(distance_to(on_major_road, 30), bad=3000, good=0),
"demand": standardise(people_per_km2, bad=0, good=500),
}
allowed = in_county & ~water & (slope_deg <= 15)
score = weighted_overlay(factors, {"slope": 1/3, "road": 1/3, "demand": 1/3}, allowed)On the 2,220 ร 1,611 grid โ 3.6 million cells โ the three factors, the distance transform and the overlay took 0.13โ0.33 s across three runs on a shared machine. The two checks catch the silent errors that otherwise produce a plausible-looking map:
ValueError: weights sum to 1.100, not 1
Example 2 โ AHP weights and the consistency ratio
import numpy as np
RANDOM_INDEX = {3: 0.58, 4: 0.90, 5: 1.12, 6: 1.24, 7: 1.32, 8: 1.41}
def ahp_weights(matrix, names):
"""Weights from a reciprocal pairwise comparison matrix (principal eigenvector)."""
A = np.asarray(matrix, dtype=float)
n = A.shape[0]
if not np.allclose(A * A.T, 1.0):
raise ValueError("matrix is not reciprocal: A[j, i] must equal 1 / A[i, j]")
values, vectors = np.linalg.eig(A)
i = int(np.argmax(values.real))
w = np.abs(vectors[:, i].real)
w /= w.sum()
lam = values[i].real
cr = ((lam - n) / (n - 1)) / RANDOM_INDEX[n]
print(f"lambda_max {lam:.4f} CR {cr:.4f} "
f"{'consistent' if cr < 0.1 else 'INCONSISTENT - revisit the judgements'}")
return {name: round(float(v), 3) for name, v in zip(names, w)}>>> ahp_weights([[1, 3, 2], [1/3, 1, 1/2], [1/2, 2, 1]], ["slope", "road", "demand"])
lambda_max 3.0092 CR 0.0079 consistent
{'slope': 0.54, 'road': 0.163, 'demand': 0.297}
>>> ahp_weights([[1, 3, 1/3], [1/3, 1, 3], [3, 1/3, 1]], ["slope", "road", "demand"])
lambda_max 4.3333 CR 1.1494 INCONSISTENT - revisit the judgements
{'slope': 0.333, 'road': 0.333, 'demand': 0.333}
The circular matrix returns perfectly equal weights. Without the CR check, incoherent judgements look like a deliberate choice of equal weighting.
Example 3 โ how many of the best sites survive a change of weights
import numpy as np
from scipy import ndimage
def best_sites(score, allowed, k=10, window=33, min_eligible=0.9):
"""The k best non-overlapping windows (window x window cells) by mean score."""
filled = np.nan_to_num(score, nan=0.0)
share = ndimage.uniform_filter(allowed.astype("float32"), size=window)
mean = ndimage.uniform_filter(filled, size=window) / np.maximum(share, 1e-6)
mean = np.where(share >= min_eligible, mean, -1.0)
picks = []
for _ in range(k):
r, c = np.unravel_index(np.argmax(mean), mean.shape)
picks.append((int(r), int(c)))
mean[max(0, r - window):r + window, max(0, c - window):c + window] = -1.0
return picks
def sites_kept(base, other, tolerance=17):
return sum(any(abs(a[0] - b[0]) <= tolerance and abs(a[1] - b[1]) <= tolerance
for b in other) for a in base)
def weight_sensitivity(factors, weights, allowed, step=0.1):
base = best_sites(weighted_overlay(factors, weights, allowed), allowed)
for up in weights:
for down in weights:
if up == down or weights[down] < step:
continue
w = dict(weights, **{up: weights[up] + step, down: weights[down] - step})
kept = sites_kept(base, best_sites(weighted_overlay(factors, w, allowed), allowed))
print(f"+{step} {up:7s} -{step} {down:7s}: {kept}/{len(base)} kept")On a 30 m grid a 33-cell window is about 1 km, and a 17-cell tolerance treats two windows whose centres are within about 500 m as the same site. The baseline and six shifted runs took 0.5โ1.3 s together across two runs, so there is no reason to skip the test.
Explanation
Why a Boolean overlay and a weighted overlay disagree
A Boolean overlay turns every criterion into a threshold and keeps the cells that pass all of them. With slope โค 8ยฐ, within 1 km of a main road and at least 100 people per kmยฒ, it kept 13.6% of the eligible land (177.2 kmยฒ) โ and it treats a cell that barely passes all three exactly like one that passes them easily.
A weighted overlay lets a strength in one factor compensate for a weakness in another. Of the weighted top 5% (72,537 cells), 5.4% would have failed a Boolean threshold, nearly all of them (5.2%) because they were more than 1 km from a main road while being flat and densely settled.
Neither is wrong. Compensation is the point of a weighted overlay; if a shortfall should never be compensated, it belongs in the constraints.
Why unstandardised factors are decided by their units
Adding degrees, metres and people per kmยฒ gives a number dominated by whichever input has the largest spread:
factor range sd corr with score: all cells top 5%
slope 0 .. 15.0 deg 4.3 0.087 -0.005
road 0 .. 11,384.6 m 2,189.9 0.995 -0.094
demand 0 .. 2,354.3/kmยฒ 235.5 0.442 0.939
Across the whole county the raw score is almost exactly distance to the road, because metres vary by thousands. Among the cells that reach the top 5%, which all lie near a road, distance varies with a standard deviation of only 197 m โ and density, varying by hundreds, takes over.
So the raw-unit map is not merely weighted wrongly: which factor decides the answer changes depending on where you look, and slope decides nothing at all. The top 5% overlapped the standardised top 5% by a Jaccard index of only 0.539, and 3 of the 10 best sites matched.
Why the weights mattered less than expected
Shifting 0.1 of weight between any two factors kept 9 or 10 of the 10 best sites. Moving from equal weights to the AHP weights โ slope from 0.33 to 0.54 โ kept 8, and the two top-5% sets overlapped by a Jaccard index of 0.755.
The best cells were good on everything. Their mean factor scores were 0.85 for slope, 0.87 for roads and 0.93 for demand, against 0.62, 0.38 and 0.18 across all eligible land. When a site wins on every factor, reweighting reorders the winners without replacing them.
Part of that is geography: road access and demand correlate at 0.540 here, because people live along main roads. Slope was nearly independent of both (0.073 and 0.136). Where factors pull in opposite directions โ cheap land far from customers โ the same test shows much more movement, which is why it is worth running rather than assuming.
Why water scored better than land
A constraint encodes knowledge that none of the factors contain. Nothing in slope, road distance or population density says "you cannot build on a lake", and on two of the three the lake is attractive: flat by definition, and bordered by the county's densest settlement.
That is the general failure. Any exclusion that coincides with good factor values โ floodplains in flat valleys, protected land beside towns, airfields next to highways โ rises to the top unless it is masked.
Edge cases or notes
- Minโmax rescaling depends on the extent. Stretching each factor between its own minimum and maximum changes every score when the study area grows; fixed thresholds do not.
- NoData must survive the sum. A NaN factor makes that cell NaN, which is correct; filling NaN with 0 first silently makes missing data look unsuitable.
- A surface model is not bare ground. The slope here came from the Copernicus GLO-30 DSM, which includes tree canopy and buildings; check which elevation model a slope factor was derived from.
- Distance transforms assume square cells in a projected CRS. In degrees they measure nothing in particular.
- "Top 5%" is relative. It always returns 5% of the eligible cells, however poor they are; add an absolute floor if a site must meet a standard.
- Weights are value judgements. AHP makes them consistent, not correct; record who made them.
- Cell-level rankings are noisy. Aggregate to the footprint of the decision before naming sites.
- Suitability is not demand. A suitable site says nothing about whether customers will come; that is a gravity model's question.
Internal links
- How to run a weighted site suitability analysis in Python โ the full workflow on rasters
- Fixing a suitability map that comes out all one value or all NoData โ when the overlay collapses
- Location analysis explained: catchments, accessibility and site selection โ where suitability fits among the other questions
- Slope and aspect explained โ how the slope factor is derived
- How to rasterize a vector layer in Python โ turning roads and water into grids
- The raster data model explained โ NoData, dtype and the transform
- Location-allocation explained: p-median, coverage and what each optimises โ choosing several sites together
- Gravity and Huff models explained โ whether people will actually come
FAQ
What is the difference between a constraint and a factor?
A constraint excludes a place outright and is applied as a mask. A factor grades places that remain eligible and can be traded off against other factors through the weights.
Do I have to standardise factors before combining them?
Yes. Summed in raw units on the Chittenden grid, the county-wide score correlated 0.995 with distance to road alone, and only 3 of the 10 best sites matched the standardised result.
How should I choose the weights?
State them directly for a handful of factors, or use pairwise comparisons (AHP) and check the consistency ratio is below 0.1. Then test how many of the best sites survive a shift of 0.1 between factors.
Should I use a Boolean overlay or a weighted overlay?
Use Boolean rules for true exclusions and a weighted overlay for everything that can be traded off. A pure Boolean overlay kept 13.6% of the eligible land and could not rank any of it.
Why did my best sites end up in a lake?
Because an exclusion was left to the factors instead of being masked. Unmasked, water filled 8.4% of the top 5% here, since flat and central is what the factors reward.
What cell size should a suitability analysis use?
Fine enough to resolve the narrowest factor that matters, usually slope, and then aggregate the score to the footprint of the real decision before ranking sites.