diff --git a/entries/railway_delays/__init__.py b/entries/railway_delays/__init__.py
new file mode 100644
index 0000000..ef2f657
--- /dev/null
+++ b/entries/railway_delays/__init__.py
@@ -0,0 +1,221 @@
+#! /usr/bin/env python
+# -*- coding: utf-8 -*-
+#
+# Copyright (C) 2020-2026 Tiago de Paula Peixoto <tiago@skewed.de>
+#
+# This program is free software: you can redistribute it and/or modify it under
+# the terms of the GNU Affero General Public License as published by the Free
+# Software Foundation, either version 3 of the License, or (at your option) any
+# later version.
+#
+# This program is distributed in the hope that it will be useful, but WITHOUT
+# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
+# FOR A PARTICULAR PURPOSE.  See the GNU Affero General Public License for more
+# details.
+#
+# You should have received a copy of the GNU Affero General Public License along
+# with this program.  If not, see <http://www.gnu.org/licenses/>.
+
+from .. import *
+import json
+import math
+import sqlite3
+import struct
+
+title = "Italian and Chinese railway networks with delays (2020-2024)"
+description = r"""Railway networks of Italy and of the Chinese high-speed railway, compiled by Faccin, Scagliarini and De Domenico (CoMuNe Lab) to study how train delays respond to rainfall. Nodes are railway stations, and an undirected edge connects two stations between which trains run without intermediate stations of the network. The station locations and the shapes of the lines were extracted from OpenStreetMap, and combined with train schedules and reported delays: for Italy, from the TrainStats service (trainstats.altervista.org), with the years 2021-2023 used for model fitting and 2024 for testing; for China, from the dataset of train operation records and weather data of Zhang et al. (2022). Each station was assigned the population of its Voronoi cell (from WorldPop estimates), and a precipitation value from the ERA5-Land reanalysis of the Copernicus Climate Change Service. The coordinates are longitude and latitude in WGS 84 (EPSG:4326/CRS84).
+
+- `italy`: the Italian railway network. The upstream edge list contains 14,572 directed edges, which list both directions for 6,016 pairs of stations and a single one for 2,540 pairs; as in the authors' code, they are merged into 8,556 undirected edges. For 13 of the pairs listed in both directions the two geometries differ, and the `length` is that of the first one listed. Since edges join consecutive stops, some of them (e.g. of long-distance trains) connect stations hundreds of kilometers apart. The network has two connected components (the largest with 2,367 stations).
+- `china`: the Chinese high-speed railway network, from the GeoPackage of the repository (the station metadata from its separate GeoJSON file). The edges were obtained by merging the OpenStreetMap ways of the lines between stations; the edge properties below are OpenStreetMap tags of these ways, one value per edge as given in the data (other tags, mostly sparse names in various languages, are left out). Three pairs of homonymous stations were split by the authors into `<name>_0` and `<name>_1` (Jingzhou, Fengchengdong, Taian); for these, the `english_name`, `city` and `province` are those of only one of the two stations. The `degree` attribute of the node data (not equal to the degree in this network) is left out. The network has two connected components (the largest with 702 stations).
+
+The time series of the delays per station and the counts of train transitions per hour, also in the repositories, are not included.
+
+The networks have the following properties:
+
+- `name`: the station identifier (lower-case Italian name, for `italy`; pinyin name, for `china`).
+- `pos`: the position of the station (longitude, latitude).
+- `pop`: the population assigned to the station.
+- `capacity`: the maximum number of trains passing through the station in one hour.
+- `rain`: the average rain at the station (ERA5-Land precipitation, as given in the data; the units and the averaging period are not documented). For `china`, it is nonzero for only 6 stations.
+- `italy` only: `osm_name`, the name of the station in OpenStreetMap; `railway`, its OpenStreetMap `railway` tag (`station`, `stop`, `halt`, `manual`, `service_station`; empty for 4 stations); `region`, the Italian region; `delay_q10`, `delay_q50`, `delay_q90` and `delay_mean`, the 10%, 50% and 90% percentiles and the mean of the delays experienced by the station (in minutes).
+- `china` only: `osmid`, the OpenStreetMap node identifier; `english_name`, the English name; `name_zh`, the Chinese name (empty if missing); `osm_name` and `osm_name_en`, the `name` and `name:en` tags in OpenStreetMap (empty if missing); `city`; `province`.
+- Edge `length`: the length of the line between the two stations, along its geometry, in kilometers.
+- `china` only, edges: `line_name` and `line_name_en`, the `name` and `name:en` tags; `railway`, `highspeed`, `usage`, `service`, `electrified`, the OpenStreetMap tags of the same names; `maxspeed` and `maxspeed_designed`, the maximum speed and the designed maximum speed in km/h (-1 if missing). Empty strings mean missing tags."""
+tags = ['Transportation', 'Spatial', 'Unweighted', 'Metadata']
+url = ['https://github.com/CoMuNeLab/adaptsys_train-delays',
+       'https://github.com/CoMuNeLab/adaptsys_chinese_train_delays']
+citation = [('M. Faccin, T. Scagliarini, M. De Domenico, "Systemic collapse emergence from climate-induced perturbations in transportation networks." Proc. Natl. Acad. Sci. U.S.A. 123(42), e2532180123 (2026)', 'https://doi.org/10.1073/pnas.2532180123'),
+            ('D. Zhang, Y. Peng, Y. Xu, C. Du, Y. Zhang, N. Wang, Y. Chong, H. Wang, D. Wu, J. Liu, H. Zhang, L. Lu, J. Liu, "A high-speed railway network dataset from train operation records and weather data." Scientific Data 9, 244 (2022)', 'https://doi.org/10.1038/s41597-022-01349-8')]
+icon_hash = None
+upstream_license = ("MIT License (Italian network); none stated (Chinese network)",
+                    "https://github.com/CoMuNeLab/adaptsys_train-delays/blob/8104d02cd15560dd1eccdd572018ca5c6d28757d/LICENSE")
+# Both repositories at fixed commits; the repository and commit are part of
+# the file names, since both contain a data/graph_nodes_metadata.geojson
+upstream_prefix = 'https://raw.githubusercontent.com/CoMuNeLab'
+ita_repo = 'adaptsys_train-delays/8104d02cd15560dd1eccdd572018ca5c6d28757d/data'
+chn_repo = 'adaptsys_chinese_train_delays/eb1f7dc54d2b8fad583a0939625bf99b450818f6/data'
+files = [((f'{ita_repo}/graph_edges.geojson',
+           f'{ita_repo}/graph_nodes_metadata.geojson'), 'italy', None),
+         ((f'{chn_repo}/graph.gpkg',
+           f'{chn_repo}/graph_nodes_metadata.geojson'), 'china', None)]
+
+def fetch_upstream(force=False):
+    return fetch_upstream_files(__name__.split(".")[-1], upstream_prefix, files,
+                                force)
+
+def _line_length(coords):
+    """Length in km of a polyline of (longitude, latitude) points (haversine)."""
+    L = 0
+    for (lon1, lat1), (lon2, lat2) in zip(coords[:-1], coords[1:]):
+        lon1, lat1, lon2, lat2 = map(math.radians, (lon1, lat1, lon2, lat2))
+        a = (math.sin((lat2 - lat1) / 2) ** 2 +
+             math.cos(lat1) * math.cos(lat2) * math.sin((lon2 - lon1) / 2) ** 2)
+        L += 2 * 6371.0088 * math.asin(math.sqrt(min(a, 1)))
+    return L
+
+def _geojson_length(geom):
+    if geom["type"] == "LineString":
+        return _line_length(geom["coordinates"])
+    if geom["type"] == "MultiLineString":
+        return sum(_line_length(c) for c in geom["coordinates"])
+    raise ValueError(f"unexpected geometry: {geom['type']}")
+
+def _wkb_lines(b, pos=0):
+    """Parse a WKB LineString or MultiLineString (2D) at offset pos; returns
+    the list of the coordinate lists of its lines and the end offset."""
+    bo = "<" if b[pos] == 1 else ">"
+    gtype, = struct.unpack_from(bo + "I", b, pos + 1)
+    pos += 5
+    if gtype == 2:
+        n, = struct.unpack_from(bo + "I", b, pos)
+        pos += 4
+        xy = struct.unpack_from(bo + "%dd" % (2 * n), b, pos)
+        pos += 16 * n
+        return [list(zip(xy[::2], xy[1::2]))], pos
+    if gtype == 5:
+        n, = struct.unpack_from(bo + "I", b, pos)
+        pos += 4
+        lines = []
+        for i in range(n):
+            ls, pos = _wkb_lines(b, pos)
+            lines += ls
+        return lines, pos
+    raise ValueError(f"unexpected WKB geometry type: {gtype}")
+
+def _gpkg_length(blob):
+    """Length in km of a GeoPackage (Multi)LineString geometry blob."""
+    assert blob[:2] == b"GP"
+    flags = blob[3]
+    env = {0: 0, 1: 32, 2: 48, 3: 48, 4: 64}[(flags >> 1) & 7]
+    lines, _ = _wkb_lines(blob, 8 + env)
+    return sum(_line_length(c) for c in lines)
+
+def _str(x):
+    return "" if x is None else str(x)
+
+def _speed(x):
+    try:
+        return int(x)
+    except (TypeError, ValueError):
+        return -1
+
+@cache_network()
+@coerce_props()
+@annotate()
+def parse(alts=None):
+    name = __name__.split(".")[-1]
+    for fnames, alt, fmt in files:
+        if alts is not None and alt not in alts:
+            continue
+
+        with open_upstream_file(name, fnames[1], "r") as f:
+            nodes = json.load(f)["features"]
+
+        g = Graph(directed=False)
+        g.add_vertex(len(nodes))
+        key = "index" if alt == "italy" else "station_name_original"
+        idx = {f["properties"][key]: i for i, f in enumerate(nodes)}
+        assert len(idx) == len(nodes)
+
+        g.vp.name = g.new_vp("string", vals=[f["properties"][key] for f in nodes])
+        g.vp.pos = g.new_vp("vector<double>")
+        for v, f in zip(g.vertices(), nodes):
+            g.vp.pos[v] = f["geometry"]["coordinates"][:2]
+
+        if alt == "italy":
+            sprops = ["osm_name", "railway", "region"]
+            fprops = ["pop", "capacity", "rain", "delay_q10", "delay_q50",
+                      "delay_q90", "delay_mean"]
+            iprops = []
+            renames = {}
+        else:
+            sprops = ["english_name", "zh_name", "name", "name:en", "city",
+                      "province"]
+            fprops = ["pop", "rain"]
+            iprops = ["osmid", "capacity"]
+            renames = {"zh_name": "name_zh", "name": "osm_name",
+                       "name:en": "osm_name_en"}
+        for p in sprops:
+            g.vp[renames.get(p, p)] = g.new_vp("string",
+                                               vals=[_str(f["properties"][p])
+                                                     for f in nodes])
+        for p in fprops:
+            g.vp[p] = g.new_vp("double", vals=[float(f["properties"][p])
+                                               for f in nodes])
+        for p in iprops:
+            g.vp[p] = g.new_vp("int64_t", vals=[int(f["properties"][p])
+                                                for f in nodes])
+
+        g.ep.length = g.new_ep("double")
+
+        if alt == "italy":
+            with open_upstream_file(name, fnames[0], "r") as f:
+                edges = json.load(f)["features"]
+            seen = set()
+            for f in edges:
+                s, t = f["properties"]["source"], f["properties"]["target"]
+                if s == t or (s, t) in seen or (t, s) in seen:
+                    continue
+                seen.add((s, t))
+                e = g.add_edge(idx[s], idx[t])
+                g.ep.length[e] = _geojson_length(f["geometry"])
+        else:
+            eprops = {"name": "line_name", "name:en": "line_name_en",
+                      "railway": "railway", "highspeed": "highspeed",
+                      "usage": "usage", "service": "service",
+                      "electrified": "electrified"}
+            for p in eprops.values():
+                g.ep[p] = g.new_ep("string")
+            g.ep.maxspeed = g.new_ep("int")
+            g.ep.maxspeed_designed = g.new_ep("int")
+            with open_upstream_file(name, fnames[0], "rb") as f:
+                data = f.read()
+            con = sqlite3.connect(":memory:")
+            con.deserialize(data)
+            del data
+            cols = list(eprops.keys()) + ["maxspeed", "maxspeed:designed"]
+            q = ("SELECT source, target, geom, " +
+                 ", ".join(f'"{c}"' for c in cols) +
+                 " FROM edges ORDER BY fid")
+            seen = set()
+            for row in con.execute(q):
+                s, t, geom = row[:3]
+                vals = dict(zip(cols, row[3:]))
+                if s == t or (s, t) in seen or (t, s) in seen:
+                    continue
+                seen.add((s, t))
+                e = g.add_edge(idx[s], idx[t])
+                g.ep.length[e] = _gpkg_length(geom)
+                for p, pn in eprops.items():
+                    g.ep[pn][e] = _str(vals[p])
+                g.ep.maxspeed[e] = _speed(vals["maxspeed"])
+                g.ep.maxspeed_designed[e] = _speed(vals["maxspeed:designed"])
+            con.close()
+
+        yield alt, g
+
+
+def get_pos(g, sbm_fit):
+    # Spatial network: laid out with the default parameters of cfd_layout(),
+    # without the partition of the SBM fit
+    return cfd_layout(g)