Predicting Delayed Train Trajectories on the Dutch Railway Network: Explainable AI Evaluation of Topological, Operational and Weather Features with Tree Based Ensemble Methods

Jia Long Bao; Ali Mohammed Mansoor Alsahag; Seyed Sahand Mohammadi Ziabari

Preprocessing Step	Rationale / Justification	Methodological Origin
Filtering non-domestic NS services (Freight & International)	Excludes distinct routing priorities and scheduling dynamics to strictly isolate domestic passenger patterns.	Novel exclusion criteria
Stop-level trajectory formulation (Graph structure)	Increases observation density and prevents spatial sparsity by treating all intermediate stops as discrete edge trajectories.	Inherited ( brakenhoff_dynamic_2026 )
Bidirectional spatial filtering	Prevents the generation of isolated, disconnected nodes and ensures topological validity within the graph structure.	Inherited ( brakenhoff_dynamic_2026 )
Terminal stop imputation	Aligns missing terminal timestamps with scheduled times, reflecting physical starting/ending realities.	Novel contribution
Handling of full/partial cancellations	Reflects true operational impact. Includes all types of cancellations in the definition for delayed trajectories.	Methodological extension ( kampere_predicting_2025 ; brakenhoff_dynamic_2026 )
Retention of October DST weather duplicates	Timestamps are duplicated but meteorological values are distinct and valid; retaining them prevents data loss and maintains unbiased monthly aggregates.	Novel justification

Appendix

Table 2 . Summary of Key Preprocessing Decisions and Methodological Origin

Split Type	Imbalance Handling	Purpose
Time-Based (All Models)	Time-Based Split	Primary temporal generalizability evaluation. Ensures 70%/30% train data split volume for each month.
Trajectory-Based (XGBoost)	Random Under-Sampling	Secondary baseline for strict literature comparability and standardized SHAP evaluation. Ensures 70%/30% train data split of trajectories across the whole dataset
Trajectory-Based Null Model	Random Under-Sampling	Establishes a theoretical performance floor (shuffled labels). If the Null Model performs similar to non-shuffled data, then the non-shuffled data is likely a result of overfitting.

Appendix

Table 3 . Summary of experimental configurations detailing split type, imbalance handling, and purpose.

Table 4 . Time-Based Split Results: F1 Score across Simultaneous and Non-Simultaneous Testing

Table 5 . Time-Based Split Results: ROC AUC across Simultaneous and Non-Simultaneous Testing

	Absolute Threshold (0.20)		Historical Seasonal Threshold
Feature Set	Simultaneous	Non-Simultaneous	Simultaneous	Non-Simultaneous
Topological Features	0.570	0.580	0.577	0.546
Weight Features	0.624	0.618	0.613	0.562
Weather Features	0.575	0.501	0.577	0.500
Operational Features	0.581	0.538	0.592	0.533
Topological + Weight	0.612	0.608	0.621	0.564

Appendix

Table 6 . Sensitivity Analysis: XGBoost Balanced Accuracy across Alternative Target Thresholds. Note: The average null baseline across all tested configurations remained mathematically stable at approximately 0.500 and is omitted for clarity.

Figure 3 . Sensitivity Analysis Boxplots: Distribution of XGBoost Balanced Accuracy across alternative target thresholds and temporal configurations. These distributions correspond to the mean values reported in Table 6 .

Feature Set	Direction	Mdn (Sim)	Mdn (Null)	Mdn Diff	W	p-value	Effect ( $r$ )
Topological Features	Sim $>$ Null	0.582	0.501	+0.075	2599.0	$<0.001$	0.601
Weight Features	Sim $>$ Null	0.601	0.496	+0.107	2538.0	$<0.001$	0.602
Weather Features	Sim $>$ Null	0.571	0.490	+0.082	2582.0	$<0.001$	0.593
Operational Features	Sim $>$ Null	0.584	0.496	+0.090	2520.0	$<0.001$	0.564
Topological + Weight	Sim $>$ Null	0.605	0.513	+0.104	2593.5	$<0.001$	0.598
Topological + Weight + Weather	Sim $>$ Null	0.630	0.512	+0.106	2627.0	$<0.001$	0.614

Appendix

Table 7 . Wilcoxon Signed-Rank Test results comparing simultaneous XGBoost performance against the null baseline (Trajectory-based split,

\alpha=0.05

,

N=72

per feature set).

Model	Direction	Mdn (Sim)	Mdn (Non-Sim)	Mdn Diff	W	p-value	Effect ( $r$ )
Topological Features	Sim $>$ Non-Sim	0.583	0.556	+0.033	2197.0	$<0.001$	0.442
Weight Features	Sim $>$ Non-Sim	0.601	0.569	+0.041	2305.0	$<0.001$	0.494
Weather Features	Sim $>$ Non-Sim	0.571	0.500	+0.069	2556.0	$<0.001$	0.615
Operational Features	Sim $>$ Non-Sim	0.584	0.540	+0.045	2394.0	$<0.001$	0.537
Topological + Weight	Sim $>$ Non-Sim	0.606	0.571	+0.037	2308.0	$<0.001$	0.495
Topological + Weight + Weather	Sim $>$ Non-Sim	0.630	0.553	+0.074	2539.0	$<0.001$	0.606

Appendix

Table 8 . Wilcoxon Signed-Rank Test results comparing simultaneous and non-simultaneous XGBoost performance (Trajectory-based split,

\alpha=0.05

,

N=71

per feature set).

Classifier	Feature Set	Direction	Mdn (Sim)	Mdn (Non-Sim)	Mdn Diff	W	p-value	Effect ( $r$ )
LGBM	Topological Features	Sim $<$ Non-Sim	0.576	0.609	-0.029	557.0	$1.000$	0.305
	Weight Features	Sim $>$ Non-Sim	0.619	0.601	+0.019	1363.0	$0.050$	0.143
	Weather Features	Sim $>$ Non-Sim	0.584	0.500	+0.079	2179.0	$<0.001$	0.597
	Operational Features	Sim $>$ Non-Sim	0.592	0.574	+0.014	1435.0	$0.018$	0.183
	Topological + Weight	Sim $>$ Non-Sim	0.616	0.616	+0.000	1058.0	$0.619$	0.026
	Top + Weight + Weather	Sim $>$ Non-Sim	0.632	0.568	+0.062	1998.0	$<0.001$	0.496

Appendix

Table 9 . Wilcoxon Signed-Rank Test results comparing simultaneous and non-simultaneous performance across all models (Time-based split,

\alpha=0.05

, N=66 per feature set.

Model	Direction	Mdn (Sim)	Mdn (Non-Sim)	Mdn Diff	W	p-value	Effect ( $r$ )
Topological Features	Sim $<$ Non-Sim	0.574	0.570	-0.000	1212.0	$0.647$	0.032
Weight Features	Sim $>$ Non-Sim	0.620	0.607	+0.024	1876.0	$<0.001$	0.288
Weather Features	Sim $>$ Non-Sim	0.569	0.502	+0.068	2538.0	$<0.001$	0.606
Operational Features	Sim $>$ Non-Sim	0.581	0.543	+0.039	2318.0	$<0.001$	0.500
Topological + Weight	Sim $>$ Non-Sim	0.606	0.595	+0.010	1695.0	$0.008$	0.201
Topological + Weight + Weather	Sim $>$ Non-Sim	0.619	0.588	+0.032	2314.0	$<0.001$	0.498

Appendix

Table 10 . Wilcoxon Signed-Rank Test results comparing simultaneous and non-simultaneous XGBoost performance under the Absolute Threshold (Trajectory-based split,

\alpha=0.05

,

N=71

per feature set).

Model	Direction	Mdn (Sim)	Mdn (Non-Sim)	Mdn Diff	W	p-value	Effect ( $r$ )
Topological Features	Sim $>$ Non-Sim	0.580	0.547	+0.034	2178.0	$<0.001$	0.433
Weight Features	Sim $>$ Non-Sim	0.618	0.562	+0.053	2400.0	$<0.001$	0.539
Weather Features	Sim $>$ Non-Sim	0.580	0.501	+0.080	2495.0	$<0.001$	0.585
Operational Features	Sim $>$ Non-Sim	0.582	0.532	+0.056	2505.0	$<0.001$	0.590
Topological + Weight	Sim $>$ Non-Sim	0.619	0.565	+0.055	2406.0	$<0.001$	0.542
Topological + Weight + Weather	Sim $>$ Non-Sim	0.628	0.556	+0.075	2489.0	$<0.001$	0.582

Appendix

Table 11 . Wilcoxon Signed-Rank Test results comparing simultaneous and non-simultaneous XGBoost performance under the Month-Relative Threshold (Trajectory-based split,

\alpha=0.05

,

N=71

per feature set).

Step	Base Model	Enhanced Model	Direction	Mdn (Base)	Mdn (Enh)	Mdn Diff	W	p-value	Effect ( $r$ )
Adding Weights	Topological Features	Topological + Weight	Enh $>$ Base	0.547	0.565	+0.018	368.0	$<0.001$	0.611
Adding Weather	Topological + Weight	Top + Weight + Weather	Enh $<$ Base	0.565	0.556	-0.006	1931.0	$1.000$	0.314
Adding Operational	Topological + Weight	Top + Weight + Op	Enh $>$ Base	0.565	0.565	+0.003	699.0	$<0.001$	0.278
Adding Op. to Weather	Top + Weight + Weather	All Features	Enh $<$ Base	0.556	0.553	-0.006	1831.0	$0.999$	0.266
Adding Wth. to Op.	Top + Weight + Op	All Features	Enh $<$ Base	0.565	0.553	-0.014	2528.0	$1.000$	0.601
Base vs. Best	Topological Features	Top + Weight + Op	Enh $>$ Base	0.547	0.565	+0.023	21.0	$<0.001$	0.604

Appendix

Table 12 . Wilcoxon Signed-Rank Test results evaluating incremental feature engineering configurations (XGBoost Trajectory based split, Non-simultaneous testing,

\alpha=0.05

).

Model	Direction	Mdn (Sim)	Mdn (Non-Sim)	Mdn Diff	W	p-value	Effect ( $r$ )
All Features	Sim $>$ Non-Sim	0.650	0.572	+0.071	2555.0	$<0.001$	0.614

Appendix

Table 13 . Wilcoxon Signed-Rank Test results comparing simultaneous and non-simultaneous tuned XGBoost performance (Trajectory-based split,

\alpha=0.05

,

N=71

).

Feature Set	MAE	RMSE	MAE (Null)	RMSE (Null)
(a) Simultaneous Testing
Topological Features	0.0758	0.1326	0.0961	0.1574
Weight Features	0.0600	0.1196	0.0942	0.1542
Weather Features	0.0771	0.1428	0.0978	0.1622
Operational Features	0.0553	0.1145	0.0936	0.1548
Topological + Weight	0.0585	0.1141	0.0934	0.1529

Appendix

Table 14 . Regression Performance: Regulatory Punctuality Threshold (3-minute delay margin). Performance metrics are evaluated via Mean Absolute Error (MAE) and Root Mean Squared Error (RMSE) across both testing paradigms.

Feature Set	MAE	RMSE	MAE (Null)	RMSE (Null)
(a) Simultaneous Testing
Topological Features	0.1366	0.1864	0.1645	0.2200
Weight Features	0.1228	0.1758	0.1617	0.2169
Weather Features	0.1486	0.2028	0.1673	0.2241
Operational Features	0.1290	0.1827	0.1588	0.2123
Topological + Weight	0.1209	0.1721	0.1620	0.2185

Appendix

Table 15 . Regression Performance: Schedule Adherence Threshold (0-minute delay margin). Performance metrics are evaluated via Mean Absolute Error (MAE) and Root Mean Squared Error (RMSE) across both testing paradigms.

Baseline Model	Statistic	Balanced Accuracy	F1 Score	ROC AUC
Naive ( $t-1$ )	Mean	0.8080	0.7798	0.8080
	Median	0.8246	0.8147	0.8246
	Q1	0.7807	0.7172	0.7807
	Q3	0.8459	0.8703	0.8459
Seasonal Naive ( $t-12$ )	Mean	0.7190	0.6531	0.7190
	Median	0.7257	0.6612	0.7257

Appendix

Table 16 . Temporal Baselines: Predictive Performance Summary. Evaluated over

N=60

truncated months to ensure an equal comparison window between Naive (

t-1

) and Seasonal Naive (

t-12

) models.

Baseline Model	Metric	Default Rate	Default Delay Rate	Non-Def. Delay Rate	Delay Rate Gap	Non-Def. Bal. Acc.
Naive	Mean	0.0356	0.6983	0.4953	0.2030	0.8274
	Median	0.0329	0.7208	0.5355	0.2155	0.8453
Seasonal Naive	Mean	0.0525	0.6215	0.4962	0.1254	0.7377
	Median	0.0466	0.6250	0.5396	0.1760	0.7396

Appendix

Table 17 . Temporal Baselines: Survivorship-Bias and Imputation Summary (

N=60

). This table details the percentage of predictions relying on zero imputation (Default Rate) due to missing historical trajectories, alongside the resulting Delay Rate Gap and Non-Defaulted Accuracy.

Figure 4 . Temporal Baselines: Balanced Accuracy by Calendar Month. Illustrates the aggregated seasonal variance and interquartile range for the Naive and Seasonal Naive baselines.

Figure 5 . Temporal Baselines: Balanced Accuracy Over Time. Displays the chronological predictive stability of the persistence baselines over the 60-month evaluation period.

Figure 6 . Survivorship-Bias Assessment. The top panel shows the Fallback Usage (default rate) over time due to missing historical trajectories. The bottom panel displays the Delay-Rate Gap between defaulted and non-defaulted predictions.

Feature Set	$N$ (Features)	Mean $\tau$	Min $\tau$	Max $\tau$	Significant Transitions ( $p<0.05$ )
Topological Features	13	0.540	0.128	0.769	59 / 71
Weight Features	7	0.421	-0.048	0.905	12 / 71
Weather Features	14	0.449	-0.209	0.846	48 / 71
Operational Features	9	0.831	0.556	1.000	71 / 71
Topological + Weight	20	0.479	-0.158	0.705	63 / 71
Topological + Weight + Weather	34	0.367	-0.005	0.647	58 / 71

Appendix

Table 18 . Kendall rank correlation (

\tau

) results evaluating sequential SHAP feature importance stability for the XGBoost base model (Trajectory-based split,

\alpha=0.05

,

T=71

sequential transitions).

Figure 7 . Temporal trajectory of Kendall rank correlation coefficients across sequential transitions for all XGBoost feature configurations.

Figure 8 . Non-simultaneous balanced accuracy performance across months for all feature sets (XGBoost, trajectory0based split.

Category	Drift Frequency (%)	Mean D-Stat	Max D-Stat
Operational (OPS)	3.29	0.028	0.994
Topological (TOP)	6.94	0.045	0.190
Weight (WGT)	1.41	0.043	1.000
Weather (WTH)	84.81	0.645	1.000

Appendix

Table 19 . Macro-Categorical Kolmogorov-Smirnov (K-S) Drift Summary. This table quantifies the frequency and magnitude (D-Statistic) of statistically significant distributional shifts across the primary feature categories.

Topological Features (TOP)
Feature	Category	Significance Rate (%)	Mean D-Stat	Max D-Stat
src_degree	TOP	4.23	0.045	0.079
tgt_degree	TOP	1.41	0.044	0.073
src_avg_distance_in	TOP	4.23	0.047	0.162
common_neighbors	TOP	12.68	0.040	0.095
jaccard_coefficient	TOP	12.68	0.045	0.100

Appendix

Table 20 . Feature-Level Kolmogorov-Smirnov (K-S) Diagnostic Summary (Part 1: Topological and Weight Features). Details the significance rate of temporal drift across all 71 sequential transitions.

Weather Features (WTH)
Feature	Category	Significance Rate (%)	Mean D-Stat	Max D-Stat
src_avg_temperature_2m	WTH	100.00	0.878	1.000
tgt_avg_temperature_2m	WTH	100.00	0.878	1.000
src_avg_soil_temperature	WTH	100.00	0.913	1.000
tgt_avg_soil_temperature	WTH	100.00	0.914	1.000
src_avg_rain	WTH	100.00	0.784	1.000

Appendix

Table 21 . Feature-Level Kolmogorov-Smirnov (K-S) Diagnostic Summary (Part 2: Weather and Operational Features). Details the significance rate of temporal drift across all 71 sequential transitions.

Figure 9 . Spatial overview of the Dutch railway network (April 2024) infrastructure and historical passenger service coverage.

Figure 10 . Binary Class Distribution with Target Label ’Is Significantly Delayed’

Figure 11 . Monthly class balance illustrating the proportion of on-time versus significantly delayed edges over the observation period.

Figure 12 . Total volume of scheduled train services aggregated per month across the entire railway network.

Figure 13 . Mean proportion of delayed arrivals per edge over time, measured against the historical global median.

Figure 14 . Distribution of the binary target label across unique edges within the network infrastructure of NS data.

Figure 15 . Distribution of the binary target label across unique edges within the network infrastructure utilizing only origin-destination as trajectories of NS data (Kämper’s).

Figure 16 . Distribution over time for selected network and environmental features. Shaded regions denote the interquartile range (25th–75th percentiles). Note that variations in src_avg_soil_temperature are present but extremely minor, causing the interquartile range to collapse and appear as a single line.

Figure 17 . Comparative analysis of balanced accuracy distributions across all evaluated feature sets under Non-Simultaneous and Simultaneous testing conditions.

Figure 18 . Overview of class balance and XGBoost balanced accuracy performance utilizing the absolute delay threshold.

Figure 19 . Overview of class balance and XGBoost balanced accuracy performance utilizing the month relative delay threshold.

Figure 20 . Time-Series of Monthly Operational Metrics: Pure Delays and Cancellations Compared Against Total Scheduled Service Volume.

Topological Features (Unweighted)
Feature Name	Measured In
[src/tgt] degree	Count
[src/tgt] avg distance in	Kilometers
[src/tgt] avg distance out	Kilometers
[src/tgt] avg distance total	Kilometers
common neighbors	Count

Appendix

Table 22 . Overview of engineered feature sets utilized for trajectory delay prediction and their corresponding measurement units. Prefixes [src/tgt] denote source or target station calculations.

Feature	Type	Formula / Description
Node Degree	Node-level	$deg(v)=\|\{u\in V:(v,u)\in E\}\|$
Out Weighted Degree	Node-level	$wdeg_{out}(v)=\sum_{u\in N(v)}w(v,u)$ , where $w(v,u)$ is the edge weight.
In Weighted Degree	Node-level	$wdeg_{in}(v)=\sum_{u\in N(v)}w(u,v)$ , where $w(u,v)$ is the edge weight.
Total Weighted Degree	Node-level	$wdeg_{total}(v)=wdeg_{out}(v)+wdeg_{in}(v)$
Out Average Distance	Node-level	$\frac{1}{\|N(v)\|}\sum_{u\in N(v)}d(v,u)$
In Average Distance	Node-level	$\frac{1}{\|N(v)\|}\sum_{u\in N(v)}d(u,v)$

Appendix

Table 23 . Topological Feature Definitions and Mathematical Formulations. The monthly railway network is formally defined as a directed graph

G_{t}=(V_{t},E_{t})

, where the nodes

V_{t}

represent individual train stations and the directed edges

E_{t}

denote active train trajectories. For any given node

v\in V_{t}

,

deg(v)

and

wdeg(v)

denote the standard degree and weighted degree, respectively. The neighborhood set is defined as

N(v)

, and the physical spatial distance between stations is denoted by

d(u,v)

. Within this framework, the physical distance between stations was utilized as the edge distance, while the decayed operational edge weight (

W_{\text{decayed}}

) previously defined in Section 3.2.1 was formally assigned as the structural edge weight

w(u,v)

to map historical service density directly onto the network topology.

Classifier	Hyperparameters
LogisticRegression (scikit-learn 1.8.0)	penalty: {’l2’}, C: {1.0}, l1_ratio: {0.0}, fit_intercept: {True}, intercept_scaling: {1}, solver: {’lbfgs’}, max_iter: {1000}, tol: {0.0001}, dual: {False}, class_weight: {None}, random_state: {88}, warm_start: {False}, n_jobs: {None}
RandomForestClassifier (scikit-learn 1.8.0)	n_estimators: {100}, criterion: {’gini’}, max_depth: {None}, min_samples_split: {2}, min_samples_leaf: {1}, min_weight_fraction_leaf: {0.0}, max_features: {’sqrt’}, max_leaf_nodes: {None}, min_impurity_decrease: {0.0}, bootstrap: {True}, oob_score: {False}, n_jobs: {None}, random_state: {88}, class_weight: {None}, ccp_alpha: {0.0}, max_samples: {None}, monotonic_cst: {None}
XGBClassifier (xgboost 3.2.0)	objective: {’binary:logistic’}, n_estimators: {100}, learning_rate: {0.3}, max_depth: {6}, booster: {’gbtree’}, tree_method: {’auto’}, gamma: {0}, min_child_weight: {1}, max_delta_step: {0}, subsample: {1}, colsample_bytree: {1}, colsample_bylevel: {1}, colsample_bynode: {1}, reg_alpha: {0}, reg_lambda: {1}, scale_pos_weight: {1}, base_score: {auto}, random_state: {88}
LGBMClassifier (lightgbm 4.6.0)	objective: {’binary’}, boosting_type: {’gbdt’}, num_leaves: {31}, max_depth: {-1}, learning_rate: {0.1}, n_estimators: {100}, subsample_for_bin: {200000}, scale_pos_weight: {1}, min_split_gain: {0.0}, min_child_weight: {0.001}, min_child_samples: {20}, subsample: {1.0}, subsample_freq: {0}, colsample_bytree: {1.0}, reg_alpha: {0.0}, reg_lambda: {0.0}, random_state: {88}

Appendix

Table 24 . Utilized Default Hyperparameter Configurations by Classifier

Figure 21 . Monthly Year-over-Year temporal feature importance for the XGBoost model evaluating isolated baseline feature sets.

Figure 22 . Monthly Year-over-Year temporal feature importance for the XGBoost model evaluating initial combined topological + weight and weather feature sets.

Figure 23 . Monthly Year-over-Year temporal feature importance for the XGBoost model evaluating combined topological + weight + operational feature and all feature set.

Figure 24 . Aggregated (2019-2024) Month-over-Month SHAP feature importance for the XGBoost model evaluating isolated baseline feature sets.

Figure 25 . Aggregated (2019-2024) Month-over-Month SHAP feature importance for the XGBoost model evaluating initial combined topological + weight and weather feature sets.

Figure 26 . Aggregated (2019-2024) Month-over-Month SHAP feature importance for the XGBoost model evaluating combined operational features and all features set.

Figure 27 . SHAP summary plots illustrating feature impact magnitudes for the XGBoost model utilizing individual feature sets of period 2024-01.

Figure 28 . SHAP summary plots illustrating feature impact magnitudes for the XGBoost model utilizing combined feature sets of period 2024-01. Color indicates feature value (red = high, blue = low).

Figure 29 . Side-by-side confusion matrix comparisons for XGBoost models utilizing individual Topological, Weight, Weather, and Operational feature sets of period 2024-01.

Figure 30 . Side-by-side confusion matrix comparisons for XGBoost models utilizing progressively combined feature sets and the complete feature space of period 2024-01.

Figure 31 . Side-by-side confusion matrix comparisons for Tuned XGBoost models utilizing the complete feature space (Topological + Weight + Weather + Operational) of period 2024-01.

Predicting Delayed Train Trajectories on the Dutch Railway Network: Explainable AI Evaluation of Topological, Operational and Weather Features with Tree Based Ensemble Methods

Abstract

Figures & tables

Appendix

Explore similar work

RIDE: An Open Dataset and Benchmark for Train Delay Prediction

Predicting Train Delays in Finland Using Machine Learning and Weather Data

FDN: Interpretable Spatiotemporal Forecasting with Future Decomposition Networks