Spatial Epi Public Health GIS

Building an Urban Heat Island Index by Census Block

This guide, part of Remote Sensing Indices for Health Studies, builds the composite index that a city heat programme actually uses: one number per census block, ranked, with enough documentation that the ranking survives a council meeting.

Problem Context & Constraints

A heat vulnerability index is a policy instrument wearing a statistical costume. It combines several layers — surface temperature, tree canopy, impervious surface, sometimes air conditioning prevalence and age structure — into a score that decides where shade trees and cooling centres go. Each combination step embeds a value judgement about how much one factor matters relative to another, and those judgements are usually buried in an unstated weighting.

Three specific problems follow. Scale mismatch: surface temperature at 30 m, canopy at 1 m, population at block level, and air conditioning at survey-area level do not share a geography, and the index inherits the coarsest one whatever resolution it is rendered at. Standardisation choice: converting each layer to a comparable scale by z-score, by percentile rank or by min–max produces different rankings from identical inputs, and the difference is largest in the tails, which is where the programme acts. Reference ambiguity: “urban heat island” is a difference from a rural baseline, and the baseline — which rural area, measured when — is frequently left undefined, which makes the magnitude uninterpretable even when the ranking is useful.

The constraint that keeps this honest is separating the measurement from the value judgement. The layers and their standardisation are technical; the weights are policy. Both must be recorded, but only one of them is the analyst’s to choose.

Where the Technical Work Ends and the Policy Begins Four input layers with their native resolutions listed: surface temperature at thirty metres, tree canopy at one metre, impervious surface at thirty metres, and population and age structure at block level. Each is summarised to census blocks and standardised, which is marked as the technical stage. The standardised layers then pass through a weighting step which is marked as a policy decision requiring sign-off, not an analytic default. The weighted combination produces the block index and its rank, which feeds the cooling programme. A note observes that changing the weights changes which blocks receive resources. The weights are not an analytic choice — they decide who gets the trees Surface temperature 30 m, night composite Tree canopy 1 m classification Impervious surface 30 m, percent cover Population and age block level, decennial Technical: summarise each layer to blocks, standardise on one declared scale the index inherits the coarsest input geography regardless of render resolution Policy: the weights, signed off and published not a default, not a principal component, not the analyst's alone Block index and rank

Prerequisites

  • python 3.11
  • numpy 1.26.4, pandas 2.2.2, geopandas 0.14.4, rasterio 1.3.10
  • A night-time surface temperature composite derived as in the sibling guide
  • Tree canopy and impervious surface layers with documented vintages
  • Census block geometry and population, matched to the same vintage
  • A declared rural reference area for the heat island differential

Step-by-Step Solution

Summarise each layer to blocks, standardise on one declared scale, apply published weights, and test how much the ranking moves under reasonable alternatives.

# Block-level heat vulnerability index with explicit weights and a rank stability check.
# Pinned: numpy==1.26.4, pandas==2.2.2, geopandas==0.14.4
import logging

import numpy as np
import pandas as pd
import geopandas as gpd

logging.basicConfig(level=logging.INFO, format="%(asctime)s %(levelname)s %(message)s")
log = logging.getLogger("exposure.heat")

# Published weights. These come from the programme's governance body, not from
# this script, and they are part of the released documentation.
WEIGHTS = {"lst_night": 0.40, "canopy_deficit": 0.25,
           "impervious": 0.20, "age65_share": 0.15}


def standardise(s: pd.Series, method: str = "percentile") -> pd.Series:
    """Put a layer on a common 0-1 scale. The METHOD changes the ranking tails."""
    x = s.astype(float)
    if method == "percentile":
        return x.rank(pct=True, na_option="keep")
    if method == "zscore":
        return (x - x.mean()) / x.std(ddof=0)
    if method == "minmax":
        rng = x.max() - x.min()
        return (x - x.min()) / rng if rng > 0 else x * 0.0
    raise ValueError(f"unknown standardisation: {method}")


def heat_index(blocks: pd.DataFrame, weights: dict = WEIGHTS,
               method: str = "percentile", min_pop: int = 20) -> pd.DataFrame:
    """Weighted composite index per block, with sparse blocks excluded."""
    missing = [c for c in weights if c not in blocks.columns]
    if missing:
        raise KeyError(f"index inputs missing: {missing}")
    if abs(sum(weights.values()) - 1.0) > 1e-9:
        raise ValueError("weights must sum to one so the score is interpretable")

    df = blocks.copy()
    keep = df["population"] >= min_pop
    log.info("blocks: %d total, %d below population floor excluded",
             len(df), int((~keep).sum()))
    df = df[keep].copy()

    z = pd.DataFrame({c: standardise(df[c], method) for c in weights}, index=df.index)
    df["heat_index"] = sum(z[c] * w for c, w in weights.items())
    df["heat_rank"] = df["heat_index"].rank(ascending=False, method="min").astype(int)
    df["heat_decile"] = pd.qcut(df["heat_index"], 10, labels=False, duplicates="drop") + 1
    log.info("index computed for %d blocks, method=%s, weights=%s",
             len(df), method, weights)
    return df


def rank_stability(blocks: pd.DataFrame, weights: dict = WEIGHTS,
                   methods=("percentile", "zscore", "minmax"),
                   jitter: float = 0.10, n_draws: int = 200,
                   seed: int = 20260918) -> pd.DataFrame:
    """How often does a block stay in the top decile under reasonable alternatives?"""
    rng = np.random.default_rng(seed)
    base = heat_index(blocks, weights, methods[0])
    top = set(base.loc[base["heat_decile"] == 10, "geoid"])
    hits = {g: 0 for g in base["geoid"]}
    trials = 0
    for m in methods:
        for _ in range(n_draws):
            w = {k: max(0.0, v + rng.uniform(-jitter, jitter)) for k, v in weights.items()}
            total = sum(w.values())
            w = {k: v / total for k, v in w.items()}
            d = heat_index(blocks, w, m)
            for g in d.loc[d["heat_decile"] == 10, "geoid"]:
                hits[g] = hits.get(g, 0) + 1
            trials += 1
    out = pd.DataFrame({"geoid": list(hits), "top_decile_frequency":
                        [hits[g] / trials for g in hits]})
    out["in_base_top_decile"] = out["geoid"].isin(top)
    stable = out[(out["in_base_top_decile"]) & (out["top_decile_frequency"] > 0.8)]
    log.info("rank stability over %d variants: %d/%d base top-decile blocks stable",
             trials, len(stable), len(top))
    return out

The rank stability check is what turns this from an index into a defensible allocation. It perturbs the weights within a stated tolerance and swaps the standardisation method, then reports how often each block stays in the top decile. Blocks that remain there across every variant are robust priorities; blocks that appear in half the variants are on a boundary, and a programme that treats them as identical to the robust ones is claiming a precision the data do not support.

The population floor exists because a block with four residents can top the index on a hot parking lot. Excluding them is a policy choice too, and the floor belongs in the documentation alongside the weights.

Which Priority Blocks Survive a Change of Method A distribution of blocks by how often they land in the top decile across six hundred index variants with perturbed weights and different standardisation methods. A cluster at the right stays in the top decile in almost every variant and is labelled robust priority. A broad middle group appears between thirty and seventy per cent of the time and is labelled boundary: their inclusion depends on the weights chosen. A cluster at the left almost never appears. A note records that a programme should fund the robust group first and treat the boundary group as requiring a separate decision. Not every top-decile block is equally a priority 80% threshold boundary blocks inclusion depends on the weights robust 0% 50% 100% share of index variants placing the block in the top decile Fund the robust group first; the boundary group needs a decision, not an algorithm

Validation & Edge Cases

Blocks with no valid temperature cells. Small downtown blocks after masking can end up empty. Report them rather than filling them, and consider aggregating to block groups where the emptiness is widespread.

A canopy layer from a different year. High-resolution canopy products are updated irregularly. A 2018 canopy layer combined with a 2024 temperature composite will misrepresent any block where planting or removal happened, which is exactly where a programme wants to measure change.

Correlated inputs. Impervious surface and surface temperature correlate strongly, so weighting both effectively double-weights the same construct. Report the correlation matrix alongside the index; where two inputs correlate above roughly 0.8, the programme should decide whether it means to weight that construct twice.

Standardisation on a clipped extent. Percentile ranks are relative to whatever area is included. Adding the suburbs to a city-only index moves every rank. Fix the extent before computing and state it.

The rural reference. Where a heat island differential is reported rather than a rank, the rural reference area and its measurement dates must be stated. A differential against a nearby irrigated agricultural area and one against dry scrubland differ by several degrees.

INFO exposure.heat blocks: 14208 total, 2166 below population floor excluded
INFO exposure.heat index computed for 12042 blocks, method=percentile, weights={'lst_night': 0.4, 'canopy_deficit': 0.25, 'impervious': 0.2, 'age65_share': 0.15}
INFO exposure.heat rank stability over 600 variants: 812/1204 base top-decile blocks stable

The standardisation method matters most in the tails, which is exactly where a cooling programme allocates its resources:

Two Standardisations, Two Different Top Tens Ten blocks ranked by the composite heat index under two standardisation methods. The middle of the ranking is almost identical. At the top, the z-score method gives extra weight to a block with one extreme value on a single input, pushing it into first place, while the percentile method caps that extreme and ranks a block that is consistently high on every input first. The difference in the top positions is where resources are allocated, which is why the method must be chosen deliberately and tested for stability. The methods agree in the middle and disagree where the money goes percentile ranks z-scores 1 · block 214 (high on all) 2 · block 087 3 · block 331 4 · block 102 5 · block 059 1 · block 418 (one extreme input) 2 · block 214 3 · block 087 4 · block 102 5 · block 059 Positions 4 and 5 agree; positions 1 to 3 do not choose the method deliberately and report the stability check same inputs, same weights, only the standardisation changed

Compliance Notes

Publish the index with its recipe. The record needs every input layer with vintage and native resolution, the summarisation method to blocks, the standardisation method, the weights with the body that approved them and the date, the population floor, the excluded-block count, the correlation matrix among inputs, and the rank stability results.

Heat indices drive resource allocation, so they attract scrutiny from the communities they rank. Two practices reduce avoidable conflict. Publish the block-level inputs alongside the score, so a neighbourhood association can see why their block scored as it did rather than only the score. And publish the stability analysis, because a block that sits on the boundary deserves to be described as such rather than presented as a precise fifty-first place.