Exporting GTFS to Databases and Warehouses

A parsed GTFS feed only becomes useful once it lives somewhere queryable. The Python Parsing & Data Normalization overview covers turning raw .txt files into clean, typed DataFrames; this guide takes those normalized frames and lands them in two durable targets — a PostgreSQL/PostGIS database for transactional and spatial workloads, and partitioned Parquet for analytical, warehouse-style scans. Both destinations are driven from the same in-memory representation, so a single export pipeline can feed an application backend and a data-lake table without re-parsing the feed. The patterns here are production-ready for teams managing one agency or a rolling archive of hundreds of feed versions.

The Export Problem

GTFS arrives as a bundle of denormalized CSV files with permissive typing: every column is text on the wire, identifiers that look numeric ("0012") must stay strings, times can exceed 24:00:00, and referential integrity between trips.txt, stop_times.txt, and stops.txt is only a convention, never enforced. Loading that directly with a naive to_sql(if_exists="replace") produces tables with TEXT everywhere, no keys, no indexes, and no spatial columns — a store you cannot join efficiently or query geographically.

A real export layer has to do four things the CSV bundle does not. It must impose a typed schema so stop_sequence is an integer and stop_lat is a double, not text. It must materialize geometry so stops and route shapes are spatially indexable rather than loose coordinate pairs. It must enforce keys so a trip_id in stop_times provably points at a row in trips. And it must be idempotent across feed versions so that re-loading yesterday’s feed, or back-filling an archive, never duplicates rows or corrupts history. The two guides below implement the spatial database and the Parquet sides of this in full; this page frames the shared schema and the decisions that apply to both.

Prerequisites

Before running the code in this guide, ensure the following are in place:

• Python 3.9+ with the database and columnar libraries installed:

pip install pandas sqlalchemy psycopg2-binary geoalchemy2 pyarrow shapely

• A running PostgreSQL 14+ with the PostGIS extension enabled (CREATE EXTENSION postgis;), and connection credentials available as environment variables.
• A validated, normalized GTFS feed already loaded into DataFrames. If your ingestion layer is not yet producing typed frames, start with the step-by-step guide to parsing GTFS with partridge, which enforces foreign keys on read.
• For large operators, a working grasp of memory-efficient processing for large feeds, since stop_times.txt alone can exceed available RAM before it ever reaches the database.
• A convention for identifying feed versions — typically a feed_version string and a feed_date — following agency metadata and feed versioning practices.

Concept and Spec Background

The GTFS relational model

GTFS is already a relational schema in disguise. Each file is a table, and the specification defines the keys that connect them even though CSV cannot enforce them. Designing the target database is largely a matter of writing those implicit keys down explicitly.

The two composite keys — (trip_id, stop_sequence) and (shape_id, shape_pt_sequence) — are the ones most often missed, and they are exactly the keys a warehouse needs to deduplicate rows during an idempotent reload.

Typed columns vs. GTFS’s text-on-the-wire

Because every field arrives as text, the export layer is where typing decisions are made once and enforced forever. Identifiers stay TEXT/VARCHAR to preserve leading zeros. stop_sequence, shape_pt_sequence, direction_id, and route_type become integers. stop_lat, stop_lon, and shape_dist_traveled become DOUBLE PRECISION. Clock fields such as arrival_time and departure_time stay TEXT, because a value like 25:30:00 is spec-valid and would be rejected by a TIME column — normalize them to a UTC instant in a separate derived column only after reading timezone handling semantics, never by forcing the raw field into a temporal type on load.

Two storage targets, one model

The governing principle: normalize once into typed DataFrames, then fan out to both targets from that single source of truth — never re-parse the CSV bundle per destination.

Step-by-Step Implementation

Step 1: Declare a typed target schema

Write the schema as explicit DDL rather than letting to_sql infer it. Inference produces TEXT columns and no keys; explicit DDL bakes in the typing decisions from the concept section and gives you a place to attach a feed_version discriminator to every table.

import sqlalchemy as sa
from sqlalchemy import text
import logging

logging.basicConfig(level=logging.INFO, format="%(levelname)s: %(message)s")

engine = sa.create_engine(
    "postgresql+psycopg2://gtfs:gtfs@localhost:5432/transit",
    future=True,
)

DDL = """
CREATE TABLE IF NOT EXISTS agency (
    feed_version    TEXT NOT NULL,
    agency_id       TEXT NOT NULL,
    agency_name     TEXT NOT NULL,
    agency_timezone TEXT NOT NULL,
    PRIMARY KEY (feed_version, agency_id)
);

CREATE TABLE IF NOT EXISTS routes (
    feed_version TEXT NOT NULL,
    route_id     TEXT NOT NULL,
    agency_id    TEXT,
    route_type   SMALLINT NOT NULL,
    PRIMARY KEY (feed_version, route_id)
);

CREATE TABLE IF NOT EXISTS trips (
    feed_version TEXT NOT NULL,
    trip_id      TEXT NOT NULL,
    route_id     TEXT NOT NULL,
    service_id   TEXT NOT NULL,
    shape_id     TEXT,
    direction_id SMALLINT,
    PRIMARY KEY (feed_version, trip_id)
);

CREATE TABLE IF NOT EXISTS stop_times (
    feed_version   TEXT NOT NULL,
    trip_id        TEXT NOT NULL,
    stop_id        TEXT NOT NULL,
    stop_sequence  INTEGER NOT NULL,
    arrival_time   TEXT,          -- kept TEXT: values may exceed 24:00:00
    departure_time TEXT,
    PRIMARY KEY (feed_version, trip_id, stop_sequence)
);
"""

with engine.begin() as conn:
    conn.execute(text("CREATE EXTENSION IF NOT EXISTS postgis;"))
    for statement in filter(str.strip, DDL.split(";")):
        conn.execute(text(statement))
logging.info("Target schema created")

Note that arrival_time and departure_time remain TEXT. Forcing them into a TIME column would reject the overnight 25:30:00 values that the GTFS spec explicitly permits.

Step 2: Bulk-load DataFrames with COPY

pandas.DataFrame.to_sql is convenient but issues row-by-row INSERTs. For anything larger than a small agency table, stream a CSV buffer through psycopg2’s COPY, which is dramatically faster because the server parses the payload once.

import io
import pandas as pd

def copy_dataframe(engine, frame: pd.DataFrame, table: str) -> int:
    """Bulk-load a DataFrame into an existing typed table via COPY."""
    buffer = io.StringIO()
    frame.to_csv(buffer, index=False, header=False, na_rep="")
    buffer.seek(0)
    columns = ", ".join(frame.columns)
    raw = engine.raw_connection()
    try:
        with raw.cursor() as cur:
            cur.copy_expert(
                f"COPY {table} ({columns}) FROM STDIN WITH (FORMAT csv, NULL '')",
                buffer,
            )
        raw.commit()
    finally:
        raw.close()
    logging.info("COPY loaded %d rows into %s", len(frame), table)
    return len(frame)

# stop_times is typically the largest table — COPY it rather than to_sql
copy_dataframe(engine, stop_times.assign(feed_version="mta_2024-04-20"), "stop_times")

Step 3: Materialize geometry columns

Stops become Point geometries and each shape becomes a LineString. Add the columns after the plain load, then populate them from the coordinate fields so the spatial index has typed geometry to work with.

from sqlalchemy import text

with engine.begin() as conn:
    conn.execute(text(
        "ALTER TABLE stops "
        "ADD COLUMN IF NOT EXISTS geom geometry(Point, 4326)"
    ))
    conn.execute(text(
        "UPDATE stops "
        "SET geom = ST_SetSRID(ST_MakePoint(stop_lon, stop_lat), 4326) "
        "WHERE geom IS NULL"
    ))
    # Assemble one LineString per shape_id, ordered by shape_pt_sequence
    conn.execute(text("""
        INSERT INTO shape_geom (feed_version, shape_id, geom)
        SELECT feed_version, shape_id,
               ST_MakeLine(
                   ST_SetSRID(ST_MakePoint(shape_pt_lon, shape_pt_lat), 4326)
                   ORDER BY shape_pt_sequence
               )
        FROM shapes
        GROUP BY feed_version, shape_id
    """))
logging.info("Geometry columns populated for stops and shapes")

The full end-to-end version of this — table creation, COPY, and geometry in one runnable script — lives in loading GTFS into PostGIS with Python. Keeping everything in EPSG:4326 matches the spec; see coordinate reference systems for transit data for when and how to reproject for metric work.

Step 4: Add foreign keys and indexes

Constraints turn the store from a pile of tables into an enforced model, and indexes make the joins and spatial queries fast.

with engine.begin() as conn:
    conn.execute(text(
        "ALTER TABLE stop_times "
        "ADD CONSTRAINT fk_stop_times_trip "
        "FOREIGN KEY (feed_version, trip_id) "
        "REFERENCES trips (feed_version, trip_id) NOT VALID"
    ))
    conn.execute(text(
        "CREATE INDEX IF NOT EXISTS ix_stop_times_stop "
        "ON stop_times (feed_version, stop_id)"
    ))
    # GiST spatial index drives ST_DWithin / ST_Intersects queries
    conn.execute(text(
        "CREATE INDEX IF NOT EXISTS ix_stops_geom "
        "ON stops USING GIST (geom)"
    ))
logging.info("Foreign keys and indexes created")

Declaring the foreign key NOT VALID lets the constraint apply to future writes without an upfront full-table scan — you can VALIDATE CONSTRAINT later during a maintenance window.

Step 5: Write partitioned Parquet from the same frames

For the warehouse side, the same DataFrames go to Parquet with an explicit Arrow schema and Hive-style partitioning. This keeps typing exact and lets query engines skip whole partitions.

import pyarrow as pa
import pyarrow.dataset as ds

stop_times_schema = pa.schema([
    ("trip_id", pa.string()),
    ("stop_id", pa.string()),
    ("stop_sequence", pa.int32()),
    ("arrival_time", pa.string()),
    ("departure_time", pa.string()),
    ("agency_id", pa.string()),
    ("feed_date", pa.string()),
])

table = pa.Table.from_pandas(
    stop_times.assign(agency_id="mta", feed_date="2024-04-20"),
    schema=stop_times_schema,
    preserve_index=False,
)

ds.write_dataset(
    table,
    base_dir="warehouse/stop_times",
    format="parquet",
    partitioning=["agency_id", "feed_date"],
    existing_data_behavior="delete_matching",  # overwrite the touched partition
)
logging.info("Wrote partitioned Parquet for stop_times")

The complete Parquet workflow — schema definition, partitioned write, and predicate-pushdown read-back — is covered in writing GTFS to partitioned Parquet.

Step 6: Make feed-version loads idempotent

Re-running a load must never duplicate rows. In PostgreSQL, stage the data and upsert on the natural key.

with engine.begin() as conn:
    conn.execute(text("""
        INSERT INTO trips AS t
            (feed_version, trip_id, route_id, service_id, shape_id, direction_id)
        SELECT feed_version, trip_id, route_id, service_id, shape_id, direction_id
        FROM trips_staging
        ON CONFLICT (feed_version, trip_id) DO UPDATE
        SET route_id     = EXCLUDED.route_id,
            service_id   = EXCLUDED.service_id,
            shape_id     = EXCLUDED.shape_id,
            direction_id = EXCLUDED.direction_id
    """))
logging.info("Idempotent upsert of trips complete")

On the Parquet side, idempotency is structural: because rows are partitioned by feed_date, existing_data_behavior="delete_matching" replaces exactly the partition being written and leaves every other feed version untouched.

Validation and Verification

Run these checks after every export to confirm the two targets agree with the source and with each other:

def verify_export(engine, source_frames: dict, feed_version: str) -> dict:
    """Confirm row counts and referential integrity after a GTFS export."""
    results = {}
    with engine.connect() as conn:
        for table, frame in source_frames.items():
            db_count = conn.execute(
                text(f"SELECT count(*) FROM {table} WHERE feed_version = :v"),
                {"v": feed_version},
            ).scalar_one()
            results[f"{table}_rows"] = (len(frame), db_count)
            assert db_count == len(frame), (
                f"{table}: loaded {db_count} rows, expected {len(frame)}"
            )

        # No stop_times may reference a trip_id absent from trips
        orphans = conn.execute(text("""
            SELECT count(*) FROM stop_times st
            LEFT JOIN trips t
              ON st.feed_version = t.feed_version AND st.trip_id = t.trip_id
            WHERE t.trip_id IS NULL AND st.feed_version = :v
        """), {"v": feed_version}).scalar_one()
        assert orphans == 0, f"{orphans} stop_times rows reference missing trips"
        results["orphaned_stop_times"] = orphans

        # Every stop must have a valid geometry
        null_geom = conn.execute(text(
            "SELECT count(*) FROM stops WHERE geom IS NULL AND feed_version = :v"
        ), {"v": feed_version}).scalar_one()
        assert null_geom == 0, f"{null_geom} stops have no geometry"
        results["stops_without_geometry"] = null_geom

    logging.info("Export verification passed: %s", results)
    return results

Failure Modes and Edge Cases

• to_sql(if_exists="replace") drops your schema. The replace mode issues a DROP TABLE, destroying your typed DDL, keys, indexes, and geometry columns, then recreates everything as TEXT. Always load into pre-created tables with if_exists="append" or COPY, never replace.

• Leading-zero identifiers coerced to integers. If a stop_id like "007" is read without dtype=str, pandas stores 7, and the COPY writes 7 into a TEXT column — now it no longer joins against stop_times. Enforce string dtypes on every ID column at read time, well before the export layer sees the frame.

• Overnight times rejected by a temporal column. A TIME or INTERVAL column will reject 25:30:00, the spec-valid encoding for 1:30 AM on the following service day. Keep the raw fields TEXT and derive normalized instants separately, following timezone handling and schedule normalization.

• Shapes assembled in the wrong order. ST_MakeLine without an ORDER BY shape_pt_sequence produces a LineString that zig-zags in row-insertion order. Always aggregate shape points ordered by their sequence, or the route geometry is geometrically wrong even though the load “succeeds”. See spatial analysis and route geometry for downstream consumers that depend on correct shape ordering.

• Reloading without a version key duplicates everything. If tables lack a feed_version column and you append a new feed, stop_times doubles. The version discriminator is not optional; it is the mechanism that makes reloads safe. Anchor it to your agency metadata and feed versioning convention.

• Parquet partition explosion. Partitioning by a high-cardinality column such as trip_id creates millions of tiny files and cripples reads. Partition only on low-cardinality keys — agency_id and feed_date — and let row groups handle the rest.

Performance and Scale Notes

For national or multi-agency archives, stop_times dominates everything else — tens of millions of rows per feed version. Three tactics keep the export tractable:

# Stream stop_times to COPY in chunks so the full frame never sits in RAM twice
import io

def copy_in_chunks(engine, path: str, table: str, chunksize: int = 500_000) -> int:
    total = 0
    reader = pd.read_csv(
        path,
        dtype={"trip_id": str, "stop_id": str, "stop_sequence": "Int32",
               "arrival_time": str, "departure_time": str},
        chunksize=chunksize,
    )
    raw = engine.raw_connection()
    try:
        with raw.cursor() as cur:
            for chunk in reader:
                buffer = io.StringIO()
                chunk.to_csv(buffer, index=False, header=False, na_rep="")
                buffer.seek(0)
                cur.copy_expert(
                    f"COPY {table} FROM STDIN WITH (FORMAT csv, NULL '')", buffer
                )
                total += len(chunk)
        raw.commit()
    finally:
        raw.close()
    logging.info("Chunked COPY loaded %d rows into %s", total, table)
    return total

First, drop indexes and constraints before a bulk load and rebuild them afterward — maintaining a GiST index during a multi-million-row COPY is far slower than building it once at the end. Second, chunk the read as above so a 5 GB stop_times.txt never fully materializes; the broader techniques are in optimizing pandas memory usage for transit feeds. Third, prefer Parquet for the archive tier: keeping every historical feed version in PostgreSQL is expensive, whereas partitioned Parquet on object storage costs almost nothing and still answers analytical queries through DuckDB or a cloud warehouse.

Frequently Asked Questions

Related

• Loading GTFS into PostGIS with Python — the complete spatial-database load: typed tables, COPY, Point and LineString geometry, and GiST indexes
• Writing GTFS to Partitioned Parquet — explicit Arrow schemas, partitioning by agency_id/feed_date, and predicate-pushdown reads
• Memory-Efficient Processing for Large Feeds — chunked reading so stop_times.txt never overruns RAM during export
• Coordinate Reference Systems for Transit Data — keeping geometry in EPSG:4326 and projecting only for metric analysis
• Agency Metadata and Feed Versioning Practices — the versioning convention that makes idempotent reloads possible

Up: Python Parsing & Data Normalization | Home