Spatiotemporal Data Foundations & Structures

Without a correct data model beneath them, even sophisticated movement analytics produce silent, hard-to-diagnose errors: projected velocities that are 40% too high, trajectory gaps that are silently interpolated as straight lines, timestamps that drift by hours across device fleets. Mobility data scientists, urban analysts, and logistics engineering teams lose weeks tracing these failures back to foundational misconfigurations in coordinate systems, sampling architecture, and storage layout.

Prerequisites & Scope

This guide assumes Python 3.10+ and familiarity with pandas DataFrames and basic GIS concepts. The primary libraries referenced are:

• geopandas ≥ 0.14 — spatial DataFrame operations and CRS management
• pyproj ≥ 3.6 — coordinate reference system transformations via PROJ 9
• movingpandas ≥ 0.17 — trajectory object construction and movement semantics
• shapely ≥ 2.0 — geometry operations (vectorized via GEOS)
• pyarrow / geoparquet — columnar storage for large trajectory datasets
• scipy — signal filtering (Savitzky-Golay, Kalman state estimation)

Data formats covered: raw GNSS NMEA, telematics JSON streams, AVL CSV exports, and GeoParquet archives. If you are new to the broader movement analytics stack, start with Movement Pattern Extraction & Trajectory Analysis to understand how the structures built here are consumed downstream.

Core Conceptual Model

Point Sequences and Trajectory Objects

Every movement dataset starts as a time-ordered sequence of coordinate-timestamp tuples (x, y, t). Raw sequences are analytically unusable at scale: they carry no entity context, bounding-box metadata, or schema guarantees. The first structural step is elevating them into trajectory objects, which encapsulate the entity identifier, temporal extent, spatial bounding box, transport mode, and confidence scores alongside the geometry. movingpandas.Trajectory is the standard Python representation; it wraps a GeoDataFrame with a DatetimeIndex and exposes vectorized spatial operations without leaving the pandas ecosystem.

Trajectory objects should support lazy evaluation and be backed by columnar storage (GeoParquet via pyarrow) to minimize memory overhead during batch processing. Crucially, the CRS must be embedded in the object schema and validated on construction—silent CRS mismatches are the single most common source of incorrect velocity calculations.

Space-Time Prisms and Reachability Envelopes

Between any two sampled positions, the true path is unknown. A space-time prism defines the feasible geographic region an entity could have occupied during that gap, bounded by maximum physically plausible velocity and the elapsed time. Rather than assuming a straight-line interpolation, the prism provides a mathematically rigorous uncertainty envelope—essential for privacy-preserving analytics, epidemiological exposure modelling, and probabilistic map-matching. When generating imputed segments after GPS dropout, represent the gap as a prism corridor, not a deterministic line.

Grid-Based and Network-Referenced Representations

Choosing between a grid and a network model is one of the most consequential architecture decisions in movement data engineering. Grid models discretize continuous coordinates into uniform cells—H3 hexagons at resolution 8 cover roughly 0.74 km² and enable O(1) spatial joins by cell key. This makes them ideal for density estimation, origin-destination matrices, and exposure analysis. Network models map coordinates to road segments using linear referencing systems, preserving the topological constraints—one-way streets, turn restrictions, grade separations—that grids cannot express. Use grids for aggregation and exposure; use networks for routing and capacity modelling.

Architecture Decision Map

Pipeline Integration

Spatiotemporal data foundations sit at the ingestion and cleaning layer of a broader movement analytics stack. The pipeline stages below are roughly sequential; each stage’s output is the next stage’s prerequisite.

Stage 1 — Ingestion. Raw streams arrive as NMEA sentences, JSON telemetry records, or CSV AVL exports. The schema contract (entity_id, timestamp, latitude, longitude, optional hdop, speed_ms, heading_deg) must be enforced at this boundary using Pydantic validators or great_expectations suites. Reject or quarantine records that fail rather than letting malformed data propagate silently.

Stage 2 — CRS normalization. Immediately after schema validation, reproject all coordinates to the analytical CRS using pyproj.Transformer. For regional datasets, UTM zone selection is deterministic from centroid longitude; for global datasets, use an equal-area projection like EPSG:6933 (WGS 84 / NSIDC EASE-Grid 2.0 Global). Store both the original WGS84 coordinates and the projected metric coordinates. Full guidance is in Coordinate Reference System Mapping.

Stage 3 — Temporal alignment. Normalize timestamps to UTC, detect and compensate for clock drift, and resample multi-source streams to a common interval. See Time-Series Synchronization Strategies for implementation patterns including NTP-based drift correction and event-driven bucketing.

Stage 4 — Noise filtering. Apply GPS Precision & Error Handling techniques: HDOP threshold filtering (reject fixes with hdop > 4.0), kinematic plausibility checks (reject points implying speed > mode-specific ceiling), and Kalman or Savitzky-Golay smoothing for trajectory continuity.

Stage 5 — Analytical layer. Clean trajectory objects flow into Movement Pattern Extraction & Trajectory Analysis for segmentation, stay-point detection, and kinematic feature extraction.

Engineering Pitfalls and Production Gotchas

1. Distance calculations in geographic coordinates

Computing Euclidean distance directly on EPSG:4326 latitude/longitude values is the most widespread accuracy bug in movement data pipelines. At 52° N latitude, one degree of longitude spans roughly 69 km, while one degree of latitude spans 111 km—treating them as equal introduces an immediate ~38% east-west distortion. Any function that accepts a GeoDataFrame and computes distances, speeds, or areas must assert gdf.crs.is_projected before proceeding.

def assert_metric_crs(gdf: gpd.GeoDataFrame) -> None:
    """Raise ValueError if the GeoDataFrame is not in a projected (metric) CRS."""
    if gdf.crs is None:
        raise ValueError("GeoDataFrame has no CRS assigned.")
    if not gdf.crs.is_projected:
        raise ValueError(
            f"CRS {gdf.crs.srs!r} is geographic, not metric. "
            "Reproject to a local UTM zone or EPSG:6933 before computing distances."
        )

2. Projection mismatches across pipeline stages

Joining two GeoDataFrame objects in different CRSs silently returns geometrically incorrect results in geopandas ≤ 0.13 (later versions raise a CRSMismatchWarning). Always validate CRS equality before spatial joins:

if left.crs != right.crs:
    right = right.to_crs(left.crs)

3. Clock drift corrupting multi-entity correlation

Device clocks drift by 0.5–5 seconds per hour in unmanaged embedded systems. Over a 12-hour fleet shift, a 2 s/h drift accumulates 24 seconds of offset. When correlating two vehicles that were co-located at a junction, 24-second drift makes them appear 600 m apart at 90 km/h. Apply NTP-anchored drift correction at the ingestion layer; store the correction offset as a column for audit purposes.

4. Aliasing from under-sampling dynamic manoeuvres

A 30-second sampling interval cannot reconstruct a 10-second turn. The trajectory appears as a straight line through the intersection, which map-matching algorithms route incorrectly. For urban fleet analytics, a minimum of 5-second sampling is required to capture turn geometry. Use adaptive logging that increases frequency during deceleration events (detectable from speed_delta between consecutive fixes).

5. Silent trajectory splits from storage partitioning

Time-based table partitioning can split a single long trajectory across partition boundaries, causing queries filtered to one day to return incomplete entities. Always store a trajectory_id that spans partitions, and implement cross-partition trajectory reassembly in your query layer or materialise complete trajectories before downstream use.

Python Tooling Landscape

Sampling Rate Optimization

Sampling Rate Optimization is a first-order cost and quality decision. The trade-off is not simply “more points = better accuracy”:

• 1 Hz fixed: captures lane changes and tight turns; appropriate for insurance telematics and safety analytics; storage-intensive (~86,400 points/device/day).
• 0.1 Hz fixed: sufficient for origin-destination studies and daily route reconstruction; misses short stops and sharp manoeuvres.
• Adaptive: logs at 1 Hz during deceleration (Δspeed > 0.5 m/s² over 2 s) and 0.05 Hz during constant highway cruising. Reduces raw point count by 60–80% against fixed 1 Hz. Requires the receiver to buffer and timestamp events locally.

For Douglas-Peucker simplification applied post-collection, a spatial tolerance of 5–10 m preserves all analytically meaningful shape for urban traces while eliminating 50–70% of redundant collinear points. Implement with shapely.simplify(tolerance=8.0, preserve_topology=True).

Further detail, including downsampling high-frequency GPS tracks without losing path integrity, is covered in the Sampling Rate Optimization section.

Storage Formats and the Archive Decision

Format choice looks like an infrastructure question and behaves like a modelling one. By the time a movement archive is a year old, the format has decided which questions are cheap, which are merely slow, and which nobody will ever ask again — and none of that is recoverable by adding hardware. Spatial storage formats covers the comparison in full; the part that belongs at the foundations level is the ordering of the decisions.

The first decision is the row shape, and it comes before the format. One row per fix keeps every per-point attribute and lets a time-range predicate push down into the reader, at the cost of repeating the entity identifier a few thousand times per trip. One row per trip with nested arrays is roughly a third of the size and loads a whole journey in a single read, which is what a sequence model wants — but it cannot answer “every fix between 14:00 and 15:00” without exploding the arrays first. Most production stacks keep the long form as the source of truth and materialise the nested form for training, because that derivation runs in one direction only.

The second decision is the partitioning key, and it is the one that decides whether row-group statistics are worth anything. A columnar file records per-group minima and maxima, so a reader can skip a group entirely when its bounding box or time range cannot satisfy the query. That mechanism is only as good as the sort order used when the file was written: rows arriving in ingestion order give every row group the full extent of the dataset and nothing can ever be skipped. Sorting by (date, spatial_cell) before writing is worth more than any compression setting, and it costs one shuffle at write time instead of a full scan on every read.

The third decision is retention, and it is where cost stops being theoretical. Raw 1 Hz fixes are irreplaceable for a few weeks — incident reconstruction, model retraining, disputed deliveries — and increasingly indistinguishable from a 5 s resample after that. Tiering by age rather than deleting outright keeps the questions that actually get asked answerable at a fraction of the storage.

Data Contracts Between Stages

Every failure mode listed above shares a shape: a downstream stage assumed something the upstream stage never promised. Distances were computed in degrees because nothing recorded that the frame was still geographic; a window closed early because nothing recorded that one partition was lagging; a stay point appeared in a car park because nothing recorded which points were interpolated. The fix is not more careful code in the consumer — it is making the promise explicit where it is produced.

A workable contract for a movement pipeline needs four things. The frame, meaning the CRS and, when the data spans years or crosses plate boundaries, the epoch as well — coordinate reference system mapping covers why the epoch matters as much as the code. The clock, meaning which of the timestamps in the provenance chain the t column actually holds, since a sensor stamp and a loader stamp can differ by half a minute and both look like ordinary UTC. The cadence, meaning whether rows are raw arrivals or lie on a regular grid, because every derivative calculation downstream depends on the answer. And the provenance of each value, meaning the quality bitfield that keeps observed, interpolated, masked and smoothed points separable after the fact.

Encode the contract as a schema the pipeline validates at the boundary rather than as documentation. A schema check that runs on every batch costs milliseconds and catches the class of bug that otherwise surfaces as a plausible-but-wrong number in a quarterly report. Version it, too: when the cadence changes from irregular arrivals to a 1 s grid, that is a new contract version, and the consumers that assumed the old one should fail loudly on the first batch rather than silently produce a different answer for the rest of the year.

The practical test of a contract is whether a new engineer can write a correct consumer without reading the producer’s source. If answering “are these coordinates in metres?” requires opening the ingest job, the contract is not written down — it is folklore, and it will be wrong the first time somebody changes the ingest job for an unrelated reason. Trajectory object design patterns sets out the concrete column-level schema this section argues for.

Validation and Testing Patterns

Spatiotemporal pipelines fail silently more often than they crash. Build these checks into your CI/CD data quality gates:

Inverse-transform round-trip. After projecting coordinates from EPSG:4326 to a metric CRS and back, the reprojected result should match the original to within 1 mm. Any larger residual indicates a datum shift misconfiguration.

import numpy as np
from pyproj import Transformer

def validate_roundtrip(lats: np.ndarray, lons: np.ndarray, target_epsg: int, tol_m: float = 0.001) -> None:
    """Assert that forward + inverse CRS transform returns to within tol_m metres."""
    fwd = Transformer.from_crs(4326, target_epsg, always_xy=True)
    inv = Transformer.from_crs(target_epsg, 4326, always_xy=True)
    xs, ys = fwd.transform(lons, lats)
    lons_rt, lats_rt = inv.transform(xs, ys)
    # Convert degree residual to approximate metres at equator (1 deg ≈ 111 km)
    residual_m = np.sqrt(((lons - lons_rt) * 111_320) ** 2 + ((lats - lats_rt) * 110_574) ** 2)
    if residual_m.max() > tol_m:
        raise AssertionError(f"CRS round-trip residual {residual_m.max():.6f} m exceeds tolerance {tol_m} m")

Velocity plausibility bounds. After noise filtering, assert that no trajectory segment exceeds a physically plausible ceiling for the declared transport mode.

MAX_SPEED_MS = {"pedestrian": 3.5, "cyclist": 12.0, "car": 55.6, "train": 97.2}

def assert_velocity_bounds(traj_gdf: gpd.GeoDataFrame, mode: str) -> None:
    ceiling = MAX_SPEED_MS.get(mode, 55.6)
    if (traj_gdf["speed_ms"] > ceiling).any():
        n_violations = (traj_gdf["speed_ms"] > ceiling).sum()
        raise AssertionError(f"{n_violations} points exceed {ceiling} m/s ceiling for mode {mode!r}")

Temporal monotonicity. Timestamps must be strictly increasing within each trajectory. Backward jumps indicate either duplicate records or clock-reset events; both corrupt kinematic derivations.

def assert_temporal_monotonicity(ts: pd.DatetimeIndex, entity_id: str) -> None:
    deltas = ts[1:] - ts[:-1]
    if (deltas.total_seconds() <= 0).any():
        raise AssertionError(f"Non-monotonic timestamps detected for entity {entity_id!r}")

Bounding-box containment. All points in a trajectory should fall within the expected geographic bounding box for the deployment region. Points outside it are almost always projection failures or hardware faults, not genuine observations.

Frequently Asked Questions

Why can’t I compute distances directly in EPSG:4326?

EPSG:4326 measures angular degrees, not metres, so Euclidean distance in degree-space is only meaningful near the equator. At 52° latitude a degree of longitude is about 68 km against 111 km at the equator, so an unprojected east-west distance is understated by roughly 38%. The failure is worse than a constant scale error, because it varies with latitude: a fleet operating across Europe accumulates a different bias in every city, and no single correction factor fixes it. Reproject to a metric CRS at ingestion, before anything computes a speed, a distance or an area.

What sampling rate should I use for urban fleet telematics?

1 Hz captures lane changes and turn behaviour and is the right default for insurance telematics and safety analytics. 0.1 Hz is adequate for origin-destination work and daily route reconstruction, and will miss short stops entirely. Adaptive logging — high frequency during acceleration and deceleration, low on constant-speed segments — cuts raw point counts by 60–80% against fixed 1 Hz with no loss of analytical fidelity, at the cost of requiring the device to buffer and stamp events locally. Choose the rate against the shortest manoeuvre you intend to detect, not against a storage budget.

How do I handle GPS dropouts without inventing false positions?

Distinguish the two cases first. A gap of a few seconds is a dropout, and dead-reckoning forward from the last known velocity and heading gives a corridor rather than a point; represent it as an uncertainty envelope so downstream models can weight it. A gap of minutes is a trajectory boundary, and interpolating across it fabricates a journey. The rule that survives contact with production is to split at a threshold — typically 30 s — and to record every imputed value in the quality bitfield so a later analysis can exclude them.

When should I use H3 cells versus a road-network model?

Grids are right when topology is irrelevant: density estimation, exposure analysis, heat maps, any aggregation where “how much activity happened near here” is the question. A network model becomes mandatory as soon as feasibility matters — routing, travel-time estimation, capacity, one-way streets, turn restrictions. A grid cannot represent a road you cannot legally take, so any analysis that implies a route needs the network. Many stacks run both, using the grid for aggregation and map matching against a road network for anything that must correspond to a real path.

How do I know whether my foundations are actually sound?

Test the invariants, not the outputs. A pipeline with sound foundations passes four assertions on every batch: a CRS round-trip that returns to within a millimetre, velocities inside a plausible ceiling for the declared transport mode, strictly increasing timestamps within each entity, and every point inside the expected bounding box. All four are cheap, all four run in CI against a fixed fixture, and each of them catches a class of error that otherwise reaches a dashboard as a number nobody questions. If a pipeline has no such fixture, the honest answer is that its foundations are untested rather than sound.

Should the archive keep raw fixes forever?

Almost never, and the interesting question is what to keep instead. Raw 1 Hz data is irreplaceable for incident reconstruction and model retraining for a few weeks, then becomes hard to distinguish from a 5 s resample for every question that actually gets asked. Tier by age rather than deleting, and hold back one un-demoted month of a sampled cohort permanently — without it, no future change to the resampler or the smoother can ever be evaluated against real raw input again.

Related

• Movement Data Privacy & Anonymization — why four points re-identify most people, and the controls that actually reduce that risk.
• GPS Precision & Error Handling — noise filtering, HDOP thresholding, and Kalman smoothing for raw GNSS streams
• Coordinate Reference System Mapping — CRS selection, PROJ transformation pipelines, and spatial join optimization
• Time-Series Synchronization Strategies — aligning asynchronous device clocks, UTC normalization, and resampling
• Sampling Rate Optimization — adaptive logging strategies and Douglas-Peucker simplification
• Trajectory Object Design Patterns — movingpandas schema design, serialization, and GeoParquet storage
• Spatial Storage Formats — row shape, partitioning keys, and the format comparison behind the archive decision
• Movement Pattern Extraction & Trajectory Analysis — the downstream stage that consumes the structures built here

Back to Home