🧭 Concept Air quality Census & demographics Air Quality Data · step 6 Notebook By · Updated

Exposure Assessment Explained

On this page

Problem statement

Epidemiology, health impact assessment and environmental justice all start from the same quantity: how much pollution each person breathes. Nobody measures that. What exists are concentration fields — station networks, models, satellite maps — and population data saying where people live. Exposure assessment is the set of assumptions that joins them, and each assumption moves the answer.

Three of them dominate. Where is a person assumed to be (at home, all the time)? At what spatial resolution is the concentration known (a 10 km cell, a street)? And at what scale is the result reported (a grid cell, a district, a city)? Change any one and the same data gives a different exposure — often in a direction that matters for who is counted as highly exposed.

Quick answer

The standard, practical definition is residential exposure: the outdoor concentration at each person's home, usually as an annual mean, summarised as a population-weighted mean and as the share of people above a threshold.

import numpy as np, rasterio

with rasterio.open("berlin_no2_on_pop_grid.tif") as src:     # from the population-weighted guide
    pop, no2 = src.read(1), src.read(2)
lived = pop > 0
print(f"{pop[lived].sum():,.0f} residents in {lived.sum():,} cells of ~1 km")
print(f"area mean of the NO2 field over populated cells: {no2[lived].mean():.1f} µg/m³")
print(f"population-weighted mean: {np.average(no2[lived], weights=pop[lived]):.1f} µg/m³")
for threshold in (10, 12, 13):
    share = pop[lived & (no2 > threshold)].sum() / pop[lived].sum()
    print(f"share of residents above {threshold} µg/m³: {share:.0%}")
3,583,694 residents in 1,609 cells of ~1 km
area mean of the NO2 field over populated cells: 11.5 µg/m³
population-weighted mean: 12.3 µg/m³
share of residents above 10 µg/m³: 86%
share of residents above 12 µg/m³: 63%
share of residents above 13 µg/m³: 36%

Every number here is conditional: on CAMS at about 10 km, on residents counted where they sleep, and on 1 km cells. The threshold shares in particular swing with each of those choices.

Define exposure, match resolution, weight by population

1. State the exposure definition

Residential annual mean, personal exposure, time-activity weighted, or peak exposure. Most studies use the first because it can be computed for everyone; it is a proxy, and saying so is part of the result.

2. Match the concentration field's resolution to the pollutant

NO2 varies within tens of metres of roads, so a 10 km model underestimates the spread of exposure. PM2.5 is more regional; coarse fields do better. Choose the field accordingly and report its resolution.

3. Weight by population, not by area

Area means count empty land. Population weighting counts people, and it is almost always higher in cities.

4. Report the distribution, not just the mean

The share of people above a guideline, and how that share changes with the threshold, is what health and equity questions need.

5. Beware of reporting scale

Aggregating to districts before computing threshold shares hides everyone whose cell is above the threshold in a district whose average is below it — the modifiable areal unit problem in its exposure form.

From a concentration field and a population grid through a location assumption and a reporting scale to an exposure figure. four assumptions, one exposure figureconcentrationfieldresolution?populationgridresidents, yearlocationat home, all dayreportingscalecell or district?exposuremean + distributionResidential exposure: measurable for everyone, a proxy for everyone spatialworkflow.io
Every link in the chain is an assumption that moves the answer.

Who counts as exposed in Berlin

The data used on this page

berlin_no2_on_pop_grid.tif is written by How to calculate population-weighted exposure in Python: GHSL 2020 residents on its 30 arc-second grid clipped to Berlin (band 1) and the 2023 annual mean NO2 from the CAMS European ensemble, via Open-Meteo, at the nearest 0.1° lattice point (band 2). If the file is not in the working directory, the setup rebuilds it the same way. District boundaries and roads come from OpenStreetMap through an Overpass mirror, read as it answered on 9 October 2026.

import io
import os
import zipfile

import numpy as np
import osmnx as ox
import rasterio
import rasterio.mask
import requests

# OpenStreetMap as Overpass returned it on 9 October 2026, saved as osmnx's own cache so the counts
# below do not drift as OSM is edited; delete the cache folder to query live instead
if not os.path.exists("cache"):
    reply = requests.get("https://data.spatialworkflow.io/osm-cache-exposure-assessment-explained-2026-10-09.zip", timeout=300)
    reply.raise_for_status()
    with zipfile.ZipFile(io.BytesIO(reply.content)) as archive:
        archive.extractall("cache")

ox.settings.overpass_url = "https://maps.mail.ru/osm/tools/overpass/api"
ox.settings.requests_timeout = 600
ox.settings.max_query_area_size = 50e9

if not os.path.exists("berlin_no2_on_pop_grid.tif"):
    berlin = ox.geocode_to_gdf("Berlin, Germany")
    tile = ("https://jeodpp.jrc.ec.europa.eu/ftp/jrc-opendata/GHSL/GHS_POP_GLOBE_R2023A/"
            "GHS_POP_E2020_GLOBE_R2023A_4326_30ss/V1-0/tiles/GHS_POP_E2020_GLOBE_R2023A_4326_30ss_V1_0_R4_C20.zip")
    zipfile.ZipFile(io.BytesIO(requests.get(tile, timeout=300).content)).extractall("ghsl")
    with rasterio.open("ghsl/GHS_POP_E2020_GLOBE_R2023A_4326_30ss_V1_0_R4_C20.tif") as src:
        grid, grid_tf = rasterio.mask.mask(src, berlin.geometry, crop=True, nodata=0)
    grid = np.where(grid[0] > 0, grid[0], 0.0)
    r, c = np.nonzero(grid)
    lon, lat = rasterio.transform.xy(grid_tf, r, c)
    lons, lats = np.round(np.arange(13.05, 13.80, 0.1), 2), np.round(np.arange(52.35, 52.70, 0.1), 2)
    pts = [(la, lo) for la in lats for lo in lons]
    annual = {}
    for i in range(0, len(pts), 20):
        batch = pts[i:i + 20]
        j = requests.get("https://air-quality-api.open-meteo.com/v1/air-quality", params={
            "latitude": ",".join(str(p[0]) for p in batch), "longitude": ",".join(str(p[1]) for p in batch),
            "hourly": "nitrogen_dioxide", "start_date": "2023-01-01", "end_date": "2023-12-31",
            "domains": "cams_europe"}, timeout=300).json()
        for p, reply in zip(batch, j):
            annual[p] = np.nanmean(np.array(reply["hourly"]["nitrogen_dioxide"], dtype=float))
    la_i = np.clip(np.round((np.array(lat) - lats[0]) / 0.1).astype(int), 0, len(lats) - 1)
    lo_i = np.clip(np.round((np.array(lon) - lons[0]) / 0.1).astype(int), 0, len(lons) - 1)
    field = np.zeros_like(grid, dtype="float32")
    field[r, c] = [annual[(lats[a], lons[b])] for a, b in zip(la_i, lo_i)]
    with rasterio.open("berlin_no2_on_pop_grid.tif", "w", driver="GTiff", height=grid.shape[0],
                       width=grid.shape[1], count=2, dtype="float32", crs="EPSG:4326",
                       transform=grid_tf, nodata=0) as dst:
        dst.write(grid.astype("float32"), 1)
        dst.write(field, 2)
print("berlin_no2_on_pop_grid.tif ready")
berlin_no2_on_pop_grid.tif ready

Example 1 — the reporting scale changes who is "exposed"

import geopandas as gpd, osmnx as ox, pandas as pd
from rasterio.features import rasterize

ox.settings.http_user_agent = "spatialworkflow-examples/1.0"
d = ox.features_from_place("Berlin, Germany", {"boundary": "administrative", "admin_level": "9"})
d = d[(d.admin_level == "9") & d.geom_type.isin(["Polygon", "MultiPolygon"])].reset_index()
with rasterio.open("berlin_no2_on_pop_grid.tif") as src:
    zone = rasterize(((g, i + 1) for i, g in enumerate(d.geometry)), out_shape=pop.shape,
                     transform=src.transform, fill=0)
df = pd.DataFrame({"zone": zone[lived], "pop": pop[lived], "no2": no2[lived]})
df = df[df.zone > 0]
district_mean = df.groupby("zone").apply(lambda g: np.average(g.no2, weights=g["pop"]),
                                         include_groups=False)
df["district_no2"] = df.zone.map(district_mean)
for threshold in (12, 12.5, 13):
    cell = df.loc[df.no2 > threshold, "pop"].sum() / df["pop"].sum()
    dist = df.loc[df.district_no2 > threshold, "pop"].sum() / df["pop"].sum()
    print(f"above {threshold:>4} µg/m³: {cell:4.0%} of residents by 1 km cell, "
          f"{dist:4.0%} by district average")
above   12 µg/m³:  63% of residents by 1 km cell,  67% by district average
above 12.5 µg/m³:  52% of residents by 1 km cell,  34% by district average
above   13 µg/m³:  36% of residents by 1 km cell,  23% by district average

The same people, the same field. Reported by district, a threshold either catches whole districts or none of their residents; reported by cell, it catches the people who actually live above it. The difference can run in either direction, and it is largest near the threshold.

Example 2 — how many people live next to a busy road

from scipy.signal import fftconvolve

roads = ox.features_from_bbox((13.08, 52.33, 13.77, 52.68),
                              {"highway": ["motorway", "trunk", "primary", "secondary"]})
roads = roads[roads.geom_type.isin(["LineString", "MultiLineString"])].to_crs(25833)
from shapely.geometry import box
import rasterio.warp, rasterio.enums

with rasterio.open("berlin_no2_on_pop_grid.tif") as src:
    bounds = gpd.GeoSeries([box(*src.bounds)], crs=4326).to_crs(25833).total_bounds
res = 10
x0, y0, x1, y1 = bounds
w, h = int((x1 - x0) / res), int((y1 - y0) / res)
tf10 = rasterio.transform.from_origin(x0, y1, res, res)
near = rasterize(((g.buffer(50), 1) for g in roads.geometry), out_shape=(h, w),
                 transform=tf10, fill=0).astype("float32")
# fraction of each ~1 km population cell that lies within 50 m of such a road
with rasterio.open("berlin_no2_on_pop_grid.tif") as src:
    frac = np.zeros(pop.shape, dtype="float32")
    rasterio.warp.reproject(near, frac, src_transform=tf10, src_crs="EPSG:25833",
                            dst_transform=src.transform, dst_crs=src.crs,
                            resampling=rasterio.enums.Resampling.average)
near_road = (pop * frac)[lived].sum()
print(f"residents living within 50 m of a major or secondary road, if spread evenly "
      f"within each cell: {near_road:,.0f} ({near_road / pop[lived].sum():.0%})")
residents living within 50 m of a major or secondary road, if spread evenly within each cell: 553,550 (15%)

About 550,000 people — 15% of residents — on this rough estimate. A 10 km model gives them the same exposure as their neighbours a street away. Road-side NO2 can be several µg/m³ higher — Berlin's seven kerbside stations averaged 24.6 µg/m³ in 2023 against 11.5 at its ten background stations, measured here — so the coarse field misclassifies a large, specific group — not at random, but everyone who lives on a main road.

Example 3 — the exposure distribution

order = np.argsort(no2[lived])
cum = np.cumsum(pop[lived][order]) / pop[lived].sum()
values = no2[lived][order]
for q in (0.1, 0.25, 0.5, 0.75, 0.9):
    print(f"{q:.0%} of residents live in cells at or below {values[np.searchsorted(cum, q)]:.1f} µg/m³")
print(f"spread between the 10th and 90th percentile resident: "
      f"{values[np.searchsorted(cum, 0.9)] - values[np.searchsorted(cum, 0.1)]:.1f} µg/m³")
10% of residents live in cells at or below 10.0 µg/m³
25% of residents live in cells at or below 11.0 µg/m³
50% of residents live in cells at or below 12.9 µg/m³
75% of residents live in cells at or below 13.5 µg/m³
90% of residents live in cells at or below 13.9 µg/m³
spread between the 10th and 90th percentile resident: 4.0 µg/m³

The whole exposure distribution on a 10 km field spans a few µg/m³ — less than the traffic increment at a single kerb. That narrowness is a property of the field, not of Berlin, and it is why the resolution of the concentration data has to be reported with every exposure figure.

Share of Berlin residents above three NO2 thresholds computed by 1 km cell and by district average. who counts as exposed> 12 µg/m³, by cell63%> 12 µg/m³, by district67%> 12.5 µg/m³, by cell52%> 12.5 µg/m³, by district34%> 13 µg/m³, by cell36%> 13 µg/m³, by district23%035%70%Compute on the finest grid you trust; aggregate results, not inputs spatialworkflow.io
The reporting scale decides who counts as exposed.

Residence, resolution and reporting scale

Why residential exposure dominates

Residence is known for everyone, from censuses and registers; activity patterns are known only from surveys of a few thousand people. Cohort studies link health records to home addresses, so the exposure metric that can be matched to outcomes is the residential one. Its errors are well studied: it tends to understate contrasts, which biases health effects towards zero rather than inventing them.

Why resolution creates misclassification

When the concentration field is smoother than reality, everyone in a cell gets the cell's mean. People by main roads are underestimated and people in quiet courtyards overestimated. In a regression of health on exposure, that is measurement error; in an equity analysis, it is the difference between seeing an injustice and missing it.

Why the reporting scale matters

Aggregation is a smoothing step of its own. District averages compress the exposure distribution; thresholds applied after aggregation count districts, not people. The rule is the same as for population weighting: compute on the finest grid you trust, aggregate the results, not the inputs.

The distribution of residential NO2 exposure in Berlin on a 10 km model field. a narrow distribution on a coarse fieldresidents at or belowNO210%10.025%11.050%12.975%13.590%13.910th to 90th percentile: 4.0 µg/m³ — less than one kerb’s incrementµg/m³; CAMS 2023 on GHSL 2020 population spatialworkflow.io
The narrowness belongs to the field, not to Berlin.

Indoor sources, mobility and ecological inference

  • Indoor sources (gas cooking for NO2, wood burning for PM2.5) can exceed outdoor contributions for some households.
  • Residential mobility: long-term studies need address histories, not a single address.
  • Night-time population grids place people at home; daytime population grids exist for some countries.
  • Uncertainty propagation: exposure errors should carry into health impact estimates as ranges.
  • Ecological inference: district-level correlations between exposure and deprivation do not describe individuals.

FAQ

What is exposure assessment?

Estimating how much of a pollutant people are exposed to, by combining concentration data with information about where people are.

What is residential exposure?

The outdoor concentration at a person's home, usually as an annual mean. It is the standard proxy in population studies.

Why use population-weighted rather than area-weighted means?

Because exposure is about people; area weighting counts empty land as heavily as dense districts.

Does the resolution of the pollution data matter?

Yes. Coarse fields misclassify people near roads, and the misclassification is not random.

Should I report exposure by district?

Compute on the finest grid you trust and aggregate the results. Thresholds applied to district averages count the wrong people.

Is residential exposure the same as what people breathe?

No. It ignores time indoors, at work and in transit, which can raise or lower personal exposure; time-activity models need survey data that population grids do not carry.