The spatial prediction of aftershocks following large earthquakes plays a critical role in seismic hazard assessment and emergency management. Traditional physical models, most notably the Coulomb Failure Stress change (ΔCFS) criterion (King et al., 1994), have been widely employed, but subsequent studies have reported significant inconsistencies in its predictive performance (e.g. Hardebeck et al., 1998; Sharma et al., 2020). Deep learning has offered a powerful data-driven alternative: DeVries et al. (2018) showed that a Deep Neural Network (DNN) trained on more than 131,000 mainshock-aftershock pairs from the SRCMOD database achieved an AUC of 0.849, far outperforming the Coulomb criterion (AUC = 0.583). However, a fundamental question remains insufficiently addressed: is the metric by which models are compared and selected valid? Most studies report performance based on cross-validation (CV), which rests on the independent and identically distributed data assumption (Vapnik, 1998). Seismic datasets violate this assumption at multiple levels. Shimodaira (2000) showed that under distribution shift, CV becomes biased; Sugiyama et al. (2007) proposed importance-weighted CV for covariate shift; but Ben-David et al. (2010) proved through impossibility theorems that under conditional shift, no domain adaptation method can guarantee success. Roberts et al. (2017) recommended block CV for structured data.

This study adopts a dual-pronged approach. First, systematic hyperparameter optimization of two neural network architectures is performed. The dataset comprises 199 slip distributions from 118 distinct mainshocks in the SRCMOD database. Following DeVries et al. (2018), static stress tensor components were computed at the centroids of 5×5×5 km cells within a volume extending 100×100×50 km around each mainshock, and aftershocks were labeled based on the ISC catalog within one second to one year after the main event. The final training set contained 164,293 samples (positive ratio 0.5225) and the test set contained 38,351 samples (positive ratio 0.4034), with zero scenario overlap verified.

Two architectures were investigated: a DNN with twelve static stress tensor features, and a shallow Artificial Neural Network (ANN) with two geometric-kinematic features (minimum distance to fault and mean slip). Bayesian Optimization with 30 trials for each architecture was executed. The optimized DNN selected 3 hidden layers with 100 neurons (21,601 parameters), ReLU activation, dropout of 0.35, Adam optimizer, and batch size of 128. It achieved a test AUC of 0.8406, consistent with DeVries et al. (2018). Threshold calibration on a held-out validation split raised F1 from 0.604 to 0.728 and Recall from 0.493 to 0.754. After calibration, the DNN-Bayes model outperformed the baseline in F1 (0.728 vs 0.720), Recall (0.754 vs 0.730), and Balanced Accuracy (0.769 vs 0.764), indicating that Bayesian Optimization ranked samples correctly and the apparent weakness at threshold 0.5 was purely a calibration issue.

The second part diagnoses the CV–Test gap through four controlled experiments. A systematic gap of +0.0638 AUC was observed between CV and test, independently confirmed by GroupKFold with five splits based on earthquake identifiers (AUC=0.7771±0.0415). Gap decomposition yielded three independent estimates of +0.0591, +0.0598, and +0.0599 AUC, all pointing to the intrinsic nature of the gap. Adversarial validation with three classifiers achieved a maximum AUC of 0.6830, indicating only moderate covariate shift. Density-ratio importance weighting (ESS = 67.8%) not only failed to improve test performance but actually reduced AUC from 0.8406 to 0.8270, fully consistent with the theoretical impossibility results of Ben-David et al. (2010) under conditional shift.

Analysis of feature-level class separability with three robust independent metrics, Cohen's d (ratio 1.4516), Mutual Information (ratio 1.4606), and Single-Feature AUC (ratio 1.2066), unanimously showed that all twelve features are more class-separable in the test set, with the largest improvements in shear components (ratios up to 1.87). Ensemble experiments provided an independent confirmation of the ranking inversion: the optimal weight from validation (w_base=0.00, only Bayes) differed completely from the optimal weight on test (w_base=1.00, only Baseline). The shallow ANN exhibited severe structural limitations because its two features contain only 199 unique values (one per earthquake); in GroupKFold it showed AUC=0.6196±0.0950 with instability more than twice that of DNN. Stability analysis across three seeds yielded AUC=0.8365±0.0040, confirming that the gap is approximately 16 times larger than training stochasticity.

These findings carry two key messages. First, Bayesian Optimization is computationally efficient, reducing evaluations by 89.6% compared to Grid Search. Second, CV-based model selection can be misleading under class-separability shift. We recommend cross-test evaluation, scenario-level GroupKFold, and threshold calibration as standard practice in aftershock prediction. The four-stage diagnostic framework developed here generalizes to any prediction task where training and test data are drawn from disjoint populations.