Building Weighted Site Suitability Scores

Turning a shortlist of candidate storefronts into a defensible ranking hinges on how rigorously you assemble the criteria matrix, choose the weights, and normalize attributes measured on wildly different scales before collapsing them into a single suitability score.

Retail expansion committees do not fund maps; they fund rankings they can defend. Multi-criteria decision analysis (MCDA) is the formal machinery that converts reachable population, demographic fit, competitive saturation, and occupancy cost into one comparable number per site. This stage of the Suitability Scoring & Site Ranking Models workflow consumes the catchment geometry and demographic attributes produced upstream and emits a scored GeoDataFrame that the ranking and shortlisting stage sorts, filters, and stress-tests. Get the weighting or normalization wrong here and every downstream shortlist inherits a bias that no amount of later sensitivity analysis can fully unwind.

Concept: the weighted-sum scoring model #

MCDA formalizes a judgment planners already make informally: some criteria matter more than others, and a site strong on the important ones should outrank a site strong only on the marginal ones. The weighted linear combination (WLC) model expresses this directly. For candidate site $i$ scored across criteria $j$, the composite suitability score is

$$S_i = \sum_{j} w_j \, x_{ij}, \qquad \sum_{j} w_j = 1, \quad w_j \ge 0$$

where $x_{ij}$ is the normalized value of criterion $j$ at site $i$ and $w_j$ is that criterion’s weight. The constraint $\sum_j w_j = 1$ is what makes scores comparable across sites and across scoring runs — without it, adding a criterion silently inflates every score.

The model only behaves if every $x_{ij}$ lives on a common, dimensionless scale. Raw inputs never do: reachable population is a count in the tens of thousands, rent is dollars per square foot, competitor distance is kilometers, and median income is dollars. Feeding those directly into the sum lets the variable with the largest raw magnitude dominate regardless of its real importance. Normalization to a shared $[0, 1]$ range — the subject of normalizing mixed-scale site attributes for scoring — is therefore a precondition of the sum, not a cosmetic step.

Directionality is the other subtlety. Criteria split into benefit criteria (more is better: reachable population, income fit) and cost criteria (more is worse: rent, competitive saturation). Cost criteria must be inverted during normalization so that a high $x_{ij}$ always means “more suitable.” Skip the inversion and your model will reward the most expensive, most saturated sites.

Assembling the scoring matrix #

The matrix $X$ has one row per candidate site and one column per criterion, with the geometry carried alongside so the result stays a spatial object. Each column originates in a different upstream process: reachable population comes from intersecting drive-time catchments (see Isochrone Generation & Network Analysis) with block-group demographics; demographic fit comes from weighting demographic variables for target audiences; competitive saturation comes from the competitor mapping and cannibalization analysis stage; and occupancy cost comes from the deal pipeline or a rent surface.

Because these columns arrive from separate systems, the assembly step is where join keys, missing values, and unit mismatches surface. Persist the assembled matrix — geometry included — to PostGIS or GeoParquet keyed by site_id so a scoring run is reproducible from a fixed snapshot rather than a live join that shifts between runs.

Weighting schemes #

The weights $w_j$ encode strategy, and three defensible methods dominate practice. Each produces a vector that must satisfy $\sum_j w_j = 1$.

• Expert / direct weights. A committee assigns weights by consensus. Fast and transparent, but prone to anchoring and hard to audit.
• Analytic Hierarchy Process (AHP). Experts perform pairwise comparisons of criteria on Saaty’s 1–9 scale; the priority vector is the principal eigenvector of the reciprocal comparison matrix. AHP surfaces inconsistency numerically via a consistency ratio (CR), and a CR above 0.10 signals contradictory judgments that must be revised.
• Regression-derived weights. When historical store performance exists, regress a performance metric (sales per square foot, first-year revenue) on the normalized criteria; the standardized coefficients become empirically grounded weights. This ties the model to observed outcomes rather than opinion, at the cost of needing a clean, sufficiently large training set.

For AHP, given a pairwise comparison matrix $A$ with $a_{jk}$ expressing how much more important criterion $j$ is than $k$, the weight vector $w$ satisfies $A w = \lambda_{\max} w$, and the consistency ratio is

$$CR = \frac{CI}{RI}, \qquad CI = \frac{\lambda_{\max} - n}{n - 1}$$

where $n$ is the number of criteria and $RI$ is the random index for that $n$. In practice teams often blend methods: AHP for a defensible starting vector, then a regression check to confirm the empirically strongest criteria are not badly under-weighted.

Configuration parameters #

The table below defines the default criteria set, their direction, plausible value ranges, and a starting weight vector for general grocery-anchored retail. Treat the weights as a documented starting point to be revised per format, not a universal truth — the sensitivity analysis for site ranking weights stage exists precisely to test how fragile a ranking is to these numbers.

Two configuration rules prevent the most common silent errors. First, the Direction column drives inversion — a cost criterion normalized as if it were a benefit criterion inverts the model’s intent without raising any error. Second, the weight vector is validated as data, not trusted as a constant: it is asserted to be non-negative and to sum to 1 before any score is computed.

Step-by-step Python implementation #

The pipeline below assembles a scored GeoDataFrame. It normalizes each criterion according to its direction, applies the validated weight vector, and computes $S_i$. Geometry is carried throughout, and because saturation is expressed per capita over a catchment area, any area computation is done in an equal-area projection (EPSG:5070, NAD83 / Conus Albers) with an explicit CRS assertion — never in geographic degrees.

import numpy as np
import pandas as pd
import geopandas as gpd
from pyproj import CRS

EQUAL_AREA_CRS = CRS.from_epsg(5070)   # NAD83 / Conus Albers, for any area math

# Criterion config: (column, direction) — 'benefit' higher-is-better, 'cost' inverted.
CRITERIA = {
    "reach_pop":  "benefit",
    "demo_fit":   "benefit",
    "comp_sat":   "cost",
    "med_income": "benefit",
    "rent_sqft":  "cost",
}

WEIGHTS = {
    "reach_pop": 0.30, "demo_fit": 0.25, "comp_sat": 0.20,
    "med_income": 0.10, "rent_sqft": 0.15,
}


def validate_weights(weights: dict, criteria: dict) -> np.ndarray:
    """Fail loudly before scoring: weights must be non-negative and sum to 1."""
    if set(weights) != set(criteria):
        raise ValueError(f"weight keys {set(weights)} != criteria {set(criteria)}")
    w = np.array([weights[c] for c in criteria])
    if (w < 0).any():
        raise ValueError("negative weight detected")
    if not np.isclose(w.sum(), 1.0, atol=1e-9):
        raise ValueError(f"weights sum to {w.sum():.6f}, expected 1.0")
    return w


def _minmax(col: pd.Series, invert: bool) -> pd.Series:
    lo, hi = col.min(), col.max()
    if np.isclose(hi, lo):                 # constant column -> neutral 0.5
        return pd.Series(0.5, index=col.index)
    scaled = (col - lo) / (hi - lo)
    return 1.0 - scaled if invert else scaled


def score_sites(sites: gpd.GeoDataFrame) -> gpd.GeoDataFrame:
    """Return `sites` with normalized criteria columns and a composite `score`."""
    if sites.crs is None:
        raise ValueError("candidate sites have no CRS; refusing to score")

    # Saturation depends on catchment area -> compute it in an equal-area CRS.
    projected = sites.to_crs(EQUAL_AREA_CRS)
    assert CRS.from_user_input(projected.crs) == EQUAL_AREA_CRS
    sites = sites.copy()
    sites["catchment_km2"] = projected.geometry.area / 1_000_000.0

    missing = [c for c in CRITERIA if c not in sites.columns]
    if missing:
        raise KeyError(f"scoring matrix missing criteria: {missing}")
    if sites[list(CRITERIA)].isna().any().any():
        raise ValueError("NaN in criteria matrix; impute or drop before scoring")

    w = validate_weights(WEIGHTS, CRITERIA)

    norm_cols = []
    for col, direction in CRITERIA.items():
        norm = _minmax(sites[col], invert=(direction == "cost"))
        sites[f"n_{col}"] = norm.clip(0.0, 1.0)   # guard against float drift
        norm_cols.append(f"n_{col}")

    matrix = sites[norm_cols].to_numpy()          # shape (n_sites, n_criteria)
    sites["score"] = matrix @ w                    # S_i = sum_j w_j x_ij
    return sites

The @ matrix-vector product is the literal implementation of $S_i = \sum_j w_j x_{ij}$: an (n_sites × n_criteria) normalized matrix multiplied by the (n_criteria,) weight vector yields the (n_sites,) score vector in one vectorized operation. Every score lands in $[0, 1]$ because each normalized column is clipped to $[0, 1]$ and the weights are a convex combination.

Integrating the result into a ranking is then a sort with a stable tie-break:

def rank_sites(scored: gpd.GeoDataFrame) -> gpd.GeoDataFrame:
    ranked = scored.sort_values(
        ["score", "reach_pop"], ascending=[False, False]
    ).reset_index(drop=True)
    ranked["rank"] = ranked.index + 1
    return ranked

Sorting on reach_pop as the secondary key makes ties deterministic, so re-running the pipeline never reshuffles equally-scored sites — a property the Huff-model catchment work relies on when it consumes ranked output.

Edge cases and failure modes #

• Constant criterion column. If every candidate has the same value for a criterion (e.g., a regional rent survey returns one number), min-max division by zero yields NaN. The _minmax guard returns a neutral 0.5 so the criterion neither helps nor hurts any site rather than poisoning the whole matrix.
• Outlier-driven compression. A single site with an extreme rent or population value stretches the min-max range so the other sites cluster near one end, flattening real differences. Winsorizing or a robust scaler mitigates this — the trade-offs are detailed in normalizing mixed-scale site attributes for scoring.
• Direction inversion errors. Labeling comp_sat a benefit criterion silently ranks the most saturated markets highest. There is no runtime error; only a monotonicity test on a synthetic site (below) catches it.
• Weight drift. Editing weights by hand until they sum to 0.98 produces scores that are quietly no longer comparable to a prior run. The validate_weights gate rejects this before any score is written.
• CRS-less geometry. A GeoDataFrame loaded from raw GeoJSON without a declared CRS will compute a nonsensical catchment area. The pipeline refuses to proceed until a CRS is asserted.

Performance and scaling #

The scoring math is trivially cheap — a single matrix multiply over a few thousand rows is microseconds. The cost is entirely in assembling the matrix: the upstream spatial joins, the catchment intersections, and the reprojection to EPSG:5070. For portfolios in the tens of thousands of candidate parcels, three tactics keep runs fast:

1. Reproject once. Convert the full candidate layer to the equal-area CRS a single time and reuse it, rather than reprojecting inside a per-site loop.
2. Vectorize normalization. Use pandas column operations, never a Python for loop over rows; the _minmax function operates on the whole Series at once.
3. Snapshot the matrix. Cache the assembled, pre-normalization matrix to GeoParquet keyed by site_id and a data_version tag. Re-scoring under new weights then reads the snapshot instead of re-running every spatial join — which is exactly what a sensitivity sweep over hundreds of weight vectors needs.

Validation and QA gates #

No score reaches the ranking stage until it clears these hard gates. Run them as assertions in the pipeline, not as advisory logs.

1. Weights sum to one. assert np.isclose(sum(weights.values()), 1.0) and all weights non-negative. This is the single most important gate; it guarantees cross-run comparability.
2. No NaN in the matrix. Every criterion cell is populated; missing values are imputed or their sites dropped before scoring, never silently coerced to zero (which would read as “worst possible” for benefit criteria).
3. Normalized range. Every n_* column lies within $[0, 1]$ after clipping; a value outside signals a normalization bug.
4. Monotonicity. Inject a synthetic site that strictly dominates another on every criterion and assert it scores higher. This one test catches direction-inversion and weight-sign errors that no range check will.

def assert_monotonic(score_fn, template: gpd.GeoDataFrame):
    """A site better on every criterion must score higher — catches sign errors."""
    base = template.iloc[[0]].copy()
    better = base.copy()
    for col, direction in CRITERIA.items():
        bump = base[col].iloc[0] * 0.10 + 1e-6
        better[col] = base[col].iloc[0] + (bump if direction == "benefit" else -bump)
    pair = gpd.GeoDataFrame(pd.concat([base, better]), crs=template.crs)
    scored = score_fn(pair)
    assert scored["score"].iloc[1] > scored["score"].iloc[0], "monotonicity violated"

Integration notes #

The scored GeoDataFrame is the handoff to the rest of the ranking stack. Persist it to PostGIS keyed by site_id, scoring_run_id, and the weight vector used, so a committee can trace any ranking back to the exact criteria and weights that produced it. The ranking and shortlisting candidate sites stage sorts and thresholds this output; the automated reporting and GIS export stage renders it into stakeholder PDFs and BI extracts.

Two upstream inputs deserve extra care because they dominate most scores. Reachable population is only as good as the catchment geometry it aggregates, and its decay-weighted refinement is the subject of calibrating distance-decay functions for trade areas. Demographic fit depends on how the underlying variables were combined, which follows the normalization discipline in normalizing mixed-scale site attributes for scoring and the target-audience weighting in the demographic section. To explore deeper, both child topics extend this page: distance decay makes the population term honest, and mixed-scale normalization makes every term commensurable.

Frequently Asked Questions #

Related #

• Suitability Scoring & Site Ranking Models — the scoring and ranking stage this work belongs to.
• Calibrating Distance-Decay Functions for Trade Areas — make the reachable-population term reflect real customer friction.
• Normalizing Mixed-Scale Site Attributes for Scoring — put every criterion on one commensurable scale.
• Ranking and Shortlisting Candidate Sites — turn scores into a defensible shortlist.
• Weighting Demographic Variables for Target Audiences — how the demographic-fit criterion is built.

← Back to Suitability Scoring & Site Ranking Models