Conceptio › Archive › arXiv CS

arXiv CSopen access

Reducing cross-sample prediction churn in scientific machine learning

2026 · arxiv_cs

arXiv CS · Papers · License: Open Access · 2026

Open Source ↗Direct PDF ↓

neural-networks

machine learning, deep learning, neural networks

Reducing cross-sample prediction churn in scientific machine learning

arXiv:2605.13826v1 [cs.LG] 13 May 2026

Gordan Prastalo1,2 [email protected]

Kevin Maik Jablonka2,3,4,5 [email protected]

Helmholtz-Zentrum Berlin für Materialien und Energie GmbH, Hahn-Meitner-Platz 1, 14109 Berlin, Germany HIPOLE Jena (Helmholtz Institute for Polymers in Energy Applications Jena), Lessingstrasse 12–14, 07743 Jena, Germany 3 Laboratory of Organic and Macromolecular Chemistry (IOMC), Friedrich Schiller University Jena, Humboldtstrasse 10, 07743 Jena, Germany 4 Center for Energy and Environmental Chemistry Jena (CEEC Jena), Friedrich Schiller University Jena, Philosophenweg 7a, 07743 Jena, Germany 5 Jena Center for Soft Matter (JCSM), Friedrich Schiller University Jena, Philosophenweg 7, 07743 Jena, Germany

Abstract Scientific machine learning reports predictive performance. It does not report whether the same prediction would survive a different draw of training data. Across 9 chemistry benchmarks, two classifiers trained on independent bootstraps of the same training set agree on aggregate accuracy to within 1.3–4.2 percentage points but disagree on the class label of 8.0–21.8% of test molecules. We call this gap cross-sample prediction churn. The standard parameter-side techniques (deep ensembles, MC dropout, stochastic weight averaging) do not reduce this gap; two data-side methods do. The first is K-bootstrap bagging, which cuts the rate 40–54% on every dataset at no accuracy cost (K×-ERM compute). The second is twinbootstrap, our proposal: two networks trained jointly on independent bootstraps with a sym-KL consistency loss between their predictions, which at matched 2×-ERM compute reduces churn a further median 45% beyond bagging-K=2. Cross-sample prediction churn deserves a column alongside predictive performance in scientific-ML benchmark reports, because without it the parameter-side and data-side methods are indistinguishable on the metric they actually differ on.

Introduction

A bioactivity classifier is retrained on a new batch of assay results. Aggregate accuracy moves by 1.8 percentage points; the class label flips on 16.1% of test molecules (Figure 1). We call this cross-sample prediction churn — the fraction of test predictions that change class between two models trained on independent samples of the same training population. This is a per-prediction stability gap that aggregate-performance reporting in scientific-ML benchmarks systematically hides. In closed-loop laboratories, Bayesian-optimisation campaigns, and virtual-screening pipelines, model predictions feed directly into experimental decisions: the molecule selected for synthesis, the candidate prioritised in a screening campaign, the hit returned by a virtual screen. Churn is the rate at which retraining changes which molecules the model selects, and is therefore the operational stability the practitioner cares about. Chemistry models are designed to be invariant to atomic permutation, rotation, and translation; E(3)-equivariant message passing and permutation-equivariant aggregators are standard. Sample invariance — invariance to which independent draw of the training population the model was fit Preprint.

Twin-bootstrap ( = 300)

Top 80 test molecules (by ERM argmax churn)

ERM

retraining 1 10

class 0

class 1

retraining 1 10

Figure 1: Twin-bootstrap eliminates most of ERM’s retraining-induced prediction flips on BACE. Each row is one of the 80 test molecules with the largest cross-sample contrast; each column is one of ten retrainings on an independent bootstrap of the BACE training pool; cells are coloured by predicted class. Visible vertical stripes in the left (ERM) panel are predictions that flip class across retrainings; under twin-bootstrap (right; λ=300, matched compute against bagging-K=2) the same molecules become near-uniform. The mean class-flip rate over the full test set drops from 16.1% (ERM) to 5.7% (twin-bootstrap). on — is the analogue on the data axis, but is not routinely measured in scientific ML. However, a practitioner who updates a screening set, retrains, and reroutes synthesis based on the new top predictions confronts the challenge of cross-sample churn on every dataset update. Two existing approaches are conceptually nearby but do not measure or reduce this gap. Prediction churn [Cormier et al., 2016, Bhojanapalli et al., 2021, Jiang et al., 2022] measures retraining-induced disagreement on web-scale data; prior work varies training conditions on a fixed dataset or compares a candidate model against a deployed incumbent. In scientific ML there is often no incumbent, and the dominant source of variance is the small training set itself, redrawn between retrainings. Epistemic uncertainty tools — deep ensembles [Lakshminarayanan et al., 2017], Monte Carlo (MC) dropout [Gal and Ghahramani, 2015], and stochastic weight averaging (SWA) [Izmailov et al., 2018] — sample over weights at fixed data and so capture the parameter-side slice of model variance, but not the data-sampling slice we measure: across our 9 chemistry benchmarks, the three together shift the class-flip rate by −22.3% to +12.5% relative to empirical risk minimisation (ERM), with no consistent sign. Methods that vary the data — bagging [Breiman, 1996] and the twin-bootstrap procedure we propose — do reduce it. Contributions. • Cross-sample prediction churn. The data-sampling analogue of model variance — the slice of uncertainty that deep ensembles, MC dropout, and SWA all miss because they sample over weights at fixed data (Section 3). • Magnitudes. On 9 chemistry benchmarks (MoleculeNet, TDC ADME and Tox, materialsscience), per-prediction churn flips 8.0–21.8% of test predictions while aggregate accuracy moves only 1.3–4.2 pp (Table 1). • Bagging. K-bootstrap bagging [Breiman, 1996] cuts churn 40–54% on every dataset at no accuracy cost, at K×-ERM compute (Table 2). • Twin-bootstrap. Two networks are trained jointly on independent bootstraps with a symmetric-KL consistency loss between their predictions (∼ 40% inter-network data overlap by construction). At matched 2×-ERM compute it cuts churn a further median 45% beyond bagging-K=2 and matches 5×-compute bagging-K=5 in mean rank (Table 2, Figure 2); the same reduction shows up at the Bayesian-optimisation acquisition layer (Appendices P and Q).

• Per-prediction triage. On top of any deployed model, sorting by per-example churn from a single extra retraining captures 58–100% of class flips at the top-30% review fraction; better than predictive entropy on every dataset (Section 6).

Related work

Prediction churn. Cormier et al. [2016] introduced prediction churn as the disagreement between two classifiers trained on the same task with different runs or data updates; their MCMC stabiliser and Jiang et al. [2022]’s distillation-from-incumbent are deployment-loop methods that compare a candidate model against a frozen incumbent on web-scale data. Bhojanapalli et al. [2021] instead hold the training data fixed and isolate parameter-side variance (initialisation, mini-batch order, hardware) as the source of churn on CIFAR/ImageNet. Codistillation operating points on the data-overlap axis. Two-network sym-KL consistency has been applied at two operating points on the inter-network data-overlap axis: Anil et al. [2018] train on disjoint shards (0% overlap) with a KL agreement loss for distributed-training throughput; Bhojanapalli et al. [2021] share the entire training set (100% overlap) for parameter-side reproducibility. Our twin-bootstrap inherits the sym-KL loss form and applies it at a third operating point: independent bootstraps of the canonical training pool, sharing ∼40% of indices in expectation (Appendix I). Predictive multiplicity and the Rashomon set. The Rashomon set [Breiman, 2001] is the set of near-loss-equivalent models on a fixed dataset; Marx et al. [2020] operationalise this as predictive multiplicity, D’Amour et al. [2020] attribute it to under-specified pipelines, and Black et al. [2022] build selective ensembles that abstain when ensemble consistency cannot be certified by a hypothesis test. These works hold the dataset fixed and vary models within the loss-equivalent class — orthogonal to cross-sample churn we discuss here, where the dataset varies and each retraining is one point in the resulting Rashomon set. Parameter-side uncertainty and bagging. Deep ensembles [Lakshminarayanan et al., 2017] and Monte Carlo dropout [Gal and Ghahramani, 2015] sample over model parameters at fixed data. Bagging [Breiman, 1996] samples over the data — an axis that has been studied for accuracy in deep learning but, to our knowledge, not for cross-sample stability. Clinical-biostatistics precedent. Riley and Collins [2023] define four levels of stability in clinical prediction-model risk estimates (mean, distribution, subgroups, individuals) and assess the individual level by retraining on 1000 bootstrap samples of the development data, reporting prediction-instability plots and a mean absolute prediction error. They demonstrate considerable instability of individualised risk estimates and argue that stability checks should be a routine part of model development — the closest published precedent for the argument we make in Section 8, in a different domain and without a deep-learning treatment.

Cross-sample prediction churn

Setting. A learner A maps a training set S = {(xi , yi )}N i=1 drawn iid from population D to a |C| classifier fS = A(S) whose output fS (x) ∈ R scores the classes c ∈ C; the predicted label at x is ŷ(x) = arg maxc fS (x)c . Cross-sample churn quantifies how ŷ(x) changes when S is replaced by an independent iid draw from D. Definition. Given two iid training samples SA , SB ∼ DN and a test example x, the cross-sample churn of A at x is ρ(A, x) = PSA ,SB arg max fSA (x) ̸= arg max fSB (x) . (1) A churn-aware report supplements the standard predictive-performance columns with the expected per-example argmax disagreement ρ̄(A) = Ex [ρ(A, x)], and optionally with a distributional analogue ρ̄KL (A) = Ex 12 (KL(fSA (x)∥fSB (x)) + KL(fSB (x)∥fSA (x))) . We use “cross-sample churn” (or just “churn” when context is clear) for ρ̄ and “class-flip rate” for the same quantity reported as a percentage; the two names refer to the same quantity. 3

Bootstrap estimator. A practitioner has one observed training set S, not a population D. We estimate ρ by drawing two independent bootstraps from S (size N , with replacement); each bootstrap is a draw from the empirical distribution D̂S , which converges to D as N → ∞. All numerical magnitudes we report are with respect to this bootstrap proxy, which is computable from a single observed training set. Measurement protocol. A canonical seed fixes one train/test split. Within a canonical seed, the training pool and test set are reused across every method and train seed, so test-side variance is zero by construction. For each train seed s ∈ {1, . . . , 10} we draw a bootstrap S (s) from the canonical training pool and train f (s) = A(S (s) ) on it. Cross-sample churn is the average of ′ ′ 1[arg max f (s) (x) ̸= arg max f (s ) (x)] over all 10 2 =45 pairs (s, s ) and over the canonical test set. Confidence intervals come from 10,000-sample bootstrap on the seed-pair distribution; cross-method deltas are paired. For the chemistry-MLP comparison we run 3 replicates; quoted magnitudes are across-replicate means and the λ-selection rule picks the same value on every replicate (Appendix F). Datasets. We evaluate on seventeen chemistry benchmarks (N ∈ [ 304, 4658 ]): MoleculeNet [Wu et al., 2018] (BACE, BBBP, ClinTox); the TDC [Huang et al., 2021] ADME and Tox suites (HIA_Hou, Bioavailability_Ma, Pgp_Broccatelli, BBB_Martins, CYP{2C9,2D6,3A4}-Sub, hERG, DILI, AMES, Skin_Reaction); and three materials-science benchmarks (TADF [Jablonka, 2026, Huang and Cole, 2024], MOF-thermal-stability, and MOF-solvent-removal [Jablonka, 2026, Nandy et al., 2022]). We require ERM accuracy to exceed the majority-class baseline by at least 5 pp on a canonical test set of at least 60 examples (threshold fixed before any cross-sample analysis); this keeps nine datasets for the method comparison, flags three borderline datasets, and excludes five (per-dataset outcomes in Appendix E). Featurisation follows standard chemistry conventions for each dataset: 2048-bit Morgan radius-2 fingerprints for MoleculeNet/TDC; +217 RDKit descriptors for TADF; 174 RAC descriptors plus geometric solvent-accessibility features for MOFs [Moosavi et al., 2020]. The default architecture is a 256-unit MLP trained from scratch on these features. For each MoleculeNet/TDC dataset we group molecules by Bemis–Murcko scaffold and partition the within-scaffold-group pool 80/20 into train and id-test under the canonical seed; TADF and the MOF benchmarks use a random 80/20 split.

Cross-sample churn flips 8.0–21.8% of test predictions on 9 chemistry benchmarks

Table 1 reports the cross-bootstrap class-flip rate and the corresponding distributional gap (sym-KL) for each of the 9 chemistry datasets that pass the +5pp ERM-vs-majority filter – on datasets where ERM fails to beat majority-class prediction by at least 5pp, cross-sample churn conflates the method shifting its decision boundary with the majority class itself shuffling under sampling noise. Excluded datasets are listed in Appendix E. Alongside the churn columns we report training-set size N , ERM accuracy, and the mean pairwise accuracy difference |∆acc| between two retrainings. Class-flip rate ranges from 8.0% (BBB-Martins) to 21.8% (MOF-thermal). Aggregate accuracy moves 1.3–4.2 pp between two retrainings while individual predictions disagree at 3–14× that rate. Sym-KL varies by ∼3× across datasets, and the variation is task-dependent rather than a function of training-set size alone (e.g. TADF and BACE both have N ≈ 1000 but differ by ∼2× in sym-KL). At these magnitudes, 8.0–21.8% of the per-molecule decisions a downstream pipeline would route to synthesis or screening change between two retrainings on the same training population — a reshuffling that is invisible in aggregate metrics such as mean accuracy. Two acquisition-layer ablations confirm this carries through to Bayesian optimisation: the top-10 Jaccard overlap of ERM’s acquisition shortlist across retrainings sits at 0.03–0.56 on chemistry benchmarks (Appendix P), and on regression the cross-trajectory standard deviation of the BO final-best y under ERM reaches 1.27–2.53% of the response range (Appendix Q). The within-task N -scaling on BACE (Appendix A) shows sym-KL trending downward with a log-log slope of −0.20 in M ; argmax churn is noisier and does not decrease monotonically with M .

Table 1: Two retrainings on independent bootstraps differ in aggregate accuracy by 1.3 – 4.2 pp on average, but disagree on 8–22% of individual test predictions. Cross-bootstrap class-flip rate on the canonical id-test of the nine chemistry datasets that pass a +5pp ERM-vs-majority filter on test sets of at least 60 examples (BACE is the development dataset; the other eight are held-out). Three datasets that pass the filter only marginally are reported in Appendix D. ERM id-acc is the mean across 10 retrainings; |∆acc| is the mean absolute accuracy difference between two retrainings, averaged over the same 45 pairs of 10 retrainings as the churn column. Class-flip rate is the perexample argmax-disagreement rate (cross-sample churn); sym-KL is the corresponding distributional gap. All paired columns report mean with 95% paired-bootstrap CIs (10,000 resamples). Aggregate accuracy

Per-prediction disagreement

Dataset

Ntrain

Nid-test

ERM id-acc

|∆acc| (pp)

Argmax churn (%)

Sym-KL (nats)

DILI CYP2D6-Sub Pgp BACE (dev) TADF MOF-thermal BBB-Martins BBBP AMES

304 427 780 968 1007 1251 1300 1305 4658

76 106 194 242 428 627 324 326 1164

0.724 [0.671, 0.776] 0.752 [0.698, 0.811] 0.842 [0.820, 0.871] 0.779 [0.756, 0.810] 0.806 [0.787, 0.822] 0.718 [0.703, 0.745] 0.867 [0.849, 0.895] 0.838 [0.819, 0.868] 0.779 [0.756, 0.799]

4.1 [3.4, 4.8] 4.2 [3.4, 4.9] 2.0 [1.6, 2.4] 1.8 [1.5, 2.2] 1.3 [1.1, 1.5] 1.5 [1.2, 1.9] 1.7 [1.4, 2.1] 1.5 [1.1, 1.8] 1.6 [1.3, 1.9]

16.8 [15.8, 17.7] 13.3 [12.6, 14.1] 10.3 [9.7, 10.9] 16.1 [15.6, 16.7] 12.7 [12.3, 13.1] 21.8 [21.2, 22.4] 8.0 [7.6, 8.4] 8.5 [8.1, 8.9] 15.2 [14.8, 15.5]

0.751 [0.692, 0.807] 0.640 [0.582, 0.701] 0.512 [0.478, 0.547] 0.753 [0.712, 0.796] 0.401 [0.384, 0.418] 0.383 [0.368, 0.398] 0.500 [0.477, 0.524] 0.470 [0.444, 0.496] 1.121 [1.092, 1.153]

Methods

We compare six approaches under the protocol of Section 3. All methods share the same architecture and optimiser unless stated otherwise: a 256-unit, two-hidden-layer MLP with ReLU activations on the dataset-specific input features (Section 3), trained from scratch with AdamW [Loshchilov and Hutter, 2019] (weight decay 10−4 , gradient clipping at 1.0), learning rate 10−3 , batch size 64, and 30 epochs. ERM.

A single model trained with AdamW on one bootstrap of the canonical training set.

Stochastic weight averaging (SWA) [Izmailov et al., 2018]. Trained as ERM, but at the end of every epoch in the second half of training we snapshot model weights and accumulate a running average. The averaged weights replace the final-epoch weights at inference time; prediction is a single forward pass through the averaged-weight model. SWA is a parameter-side smoother: it modifies the weight trajectory without varying the data sample. Cross-sample churn is measured between two such weight-averaged models trained on independent bootstraps. MC dropout. A single model with dropout layers (p = 0.2) inserted after each hidden ReLU, trained by SGD on one bootstrap. At inference time, T = 20 stochastic forward passes are averaged with dropout active [Gal and Ghahramani, 2015]. This captures dropout-mask variance at fixed data; cross-sample churn is then measured between the T -averaged predictions of two models trained on independent bootstraps. Deep ensemble. K models trained on the same bootstrap with different initialisation seeds; predictions are averaged at inference time. This captures parameter-distribution variance only. K-bootstrap bagging. The classical estimator of Breiman [1996]: K models trained on K independent bootstraps of the canonical training set, predictions averaged at inference time. No auxiliary loss; the data axis is varied along with parameters. Twin-bootstrap. Two networks θA , θB are trained jointly on independent bootstraps SA , SB . At each step we draw one mini-batch from each loader, forward both networks on both mini-batches, and minimise (t) (t) L = LCE (θA ; SA ) + LCE (θB ; SB ) + λ · Lcons , (2) where LCE is standard cross-entropy on each network’s own bootstrap, and the consistency term Lcons is the symmetric KL between the two networks’ softmax distributions evaluated on the union of the two mini-batches. At inference time, predictions are averaged. 5

Mechanism. At λ=0 the procedure is K=2 bagging. At λ > 0 the consistency term penalises divergence on every example seen in either mini-batch, including those that appear in both bootstraps. Bootstrap-with-replacement at size N covers ∼63% of unique indices in expectation, so two independent bootstraps share ∼0.632 ≈ 40% of indices. The consistency loss therefore acts on the overlap while the cross-entropy losses specialise on the non-shared remainder. Algorithm 1 gives the per-epoch procedure. Algorithm 1 Twin-bootstrap training (one epoch). Require: canonical training set S, networks θA , θB , optimisers optA , optB , weight λ 1: Sample SA , SB ∼ Bootstrap(S) ▷ N draws with replacement, independent 2: Build dataloaders LA over SA , LB over SB 3: for (BA , BB ) zipped from (LA , LB ) do (A) (A) 4: p̂A ← fθA (BA ); p̂B ← fθB (BA ) (B) (B) 5: p̂A ← fθA (BB ); p̂B ← fθB (BB ) (A) (B) 6: LCE ← CE(p̂A , yA ) + CE(p̂B , yB ) (A) (A) (B) (B) 7: Lcons ← 12 symKL(p̂A , p̂B ) + symKL(p̂A , p̂B ) 8: L ← LCE + λ Lcons 9: Backprop L through both networks; step optA , optB 10: end for Hyperparameter selection. A single development dataset (BACE) is used to select λ via the pre-registered rule: largest λ in {1, 3, 10, 30, 100, 300} such that BACE id-accuracy is within 0.02 of ERM id-accuracy. The rule yields λ=300 (Pareto curve in Appendix G). We apply λ=300 unchanged to every held-out dataset; per-dataset Bayesian optimisation of λ with a cross-sample-churn objective matches or beats this frozen choice on 6/9 datasets and strictly improves on 4/9 (Appendix H). K for ensemble methods is set a priori to 5 following Lakshminarayanan et al. [2017]; we additionally report K=2 for compute-matched comparison against twin-bootstrap.

Empirical evaluation

Cross-sample reduction across benchmarks. Across the 9 chemistry benchmarks (BACE plus 8 held-out), bagging-K=5 and twin-bootstrap both beat ERM on the cross-sample class-flip rate with 95% paired-bootstrap CIs that exclude zero (Table 2, Figure 2). The three parameter-side techniques — MC dropout, deep ensembles, and stochastic weight averaging (SWA) — do not reduce churn consistently. Their paired ∆ class-flip rate vs. ERM ranges −22.3% to +12.5% across datasets, with no consistent sign. Averaging T =20 stochastic forward passes (MC dropout), averaging K=5 networks at fixed data (deep ensembles), or averaging weight snapshots from a single training trajectory (SWA) all keep the data-resampling axis constant. At matched compute against bagging-K=2, twin-bootstrap wins on 7/8 held-out datasets (paired CI excludes zero); the comparison ties on DILI, whose canonical id-test holds only 76 examples. Against the 5×-compute bagging-K=5, twin-bootstrap wins on 5/8 held-out datasets and trails on 3 (DILI, Pgp, TADF). The aggregate Friedman test rejects equal-ranks at χ2 = 44.8, p = 5.2×10−8 ; twin-bootstrap ranks first, indistinguishable from bagging-K=5 at less than half the compute (Appendix M). What reduces churn turns out to be different from what improves accuracy: deep ensembles average five independently initialised networks at fixed data and raise mean accuracy by +0.1 pp without changing churn, while twin-bootstrap varies the data sample and cuts churn a median 65% at a small accuracy cost (−0.7 pp on average). Distributional disagreement. Argmax churn counts only whether the top-predicted class changes; it ignores how the rest of the probability distribution shifts. Symmetric KL between two bootstraps captures that and distinguishes methods more finely: bagging-K=5 reduces it 6–9× vs. ERM, and twin-bootstrap reduces it by another factor of ∼9 (per-dataset paired CIs in Table 11, Appendix L). Pipelines that consume the softmax — active-learning acquisition, Bayesian-optimisation, virtualscreening top-K ranking — see a larger gap between methods than argmax churn reflects. Appendix P quantifies this directly: twin-bootstrap raises the Jaccard overlap of the top-10 predicted-active sets between two retrainings on every chemistry dataset (paired ∆J10 ≥ +0.12 on every cell), with the 6

DILI CYP2D6-Sub Pgp BBB-Martins MOF-thermal BBBP TADF AMES -15.0 -12.5 -10.0

Paired

MC dropout

-7.5

-5.0

-2.5

0.0

id-churn vs. ERM (pp; <0 better)

SWA

Deep Ens. K=5

2.5

Bagging K=2

-5

vs. Bagging K=2 (pp) Bagging K=5

Twin-bootstrap

= 300

Figure 2: Bagging and twin-bootstrap beat ERM on every chemistry benchmark; MC dropout, deep ensembles, and SWA do not. Left: paired ∆ id-churn vs. ERM for six methods, one row per dataset (smallest N at top), 95% paired-bootstrap CIs across the 45 seed pairs. Vertical reference lines mark each method’s across-dataset mean (twin-bootstrap solid, others dashed; colours match markers); the solid black line is parity with ERM. Right: ∆ vs. matched-compute bagging-K=2, twin-bootstrap only; CI excludes zero on 7/8 datasets and ties on DILI. Frozen λ=300 is selected on the BACE development dataset alone; no held-out tuning. largest gains where ERM is most unstable. In a full BO-loop simulation on the three regression benchmarks, twin-bootstrap reduces the cross-trajectory standard deviation of the final-best acquired y by 34%–100% on every dataset; bagging-K=5 helps on 2/3 but is comparable to ERM on the third (Appendix Q). Churn-ranked triage. Per-example churn, computed across pairs of ERM bootstraps, ranks test examples by how likely they are to flip class on retraining. Figure 3 reports the cumulative fraction of total flip-mass captured if the practitioner sorts predictions by churn and reviews the top fraction. Reviewing the top 30% captures 58–100% of all flips when churn is scored from all 10 bootstraps. With a single extra retraining (K=2, the one-bootstrap practitioner workflow), the top-30% recall drops only to 48–83% – a gap of just 10–24 pp from the K=10 gold standard on every dataset (Appendix N). Predictive entropy from a single ERM model ranks worse on this task: churn captures 25–60% of flips in its top decile vs. 21–47% for entropy (per-dataset gap 1–13 pp, Appendix O). A practitioner can therefore identify the 10% most fragile predictions from a single additional retraining: one extra bootstrap on top of the deployed pipeline. A per-molecule view of these flips on six BACE id-test compounds is in Appendix K (Figure 5). 6.1

The recipe generalises across architectures and tasks

The main result on small-N from-scratch MLPs holds across three other architectures and on regression. Different architectures require different λ values, but the same accuracy-tolerance rule on the development set (Section 5, “Hyperparameter selection”) picks each one.

Table 2: Bagging and twin-bootstrap reduce the class-flip rate on every chemistry benchmark; MC dropout and deep ensembles do not. Paired ∆ id-churn vs. ERM in percentage points (negative is better). Each cell is the mean across 3 canonical-seed replicates; each replicate-mean averages over the 45 pairs that 10 retrainings produce, with paired-bootstrap 95% CIs from 10,000 resamples. Per-replicate values in Appendix F. Best per dataset in bold; entries whose CI excludes zero are significant at α=0.05. Twin-bootstrap λ=300 is selected on BACE only and applied unchanged to every held-out benchmark. Parameter-side Dataset

N 304 427 780 968 1007 1251 1300 1305 4658

Cumulative % of flip-mass captured

DILI CYP2D6-Sub Pgp BACE (dev) TADF MOF-thermal BBB-Martins BBBP AMES

Data-side

MC dropout

Deep Ens. K=5

Bagging K=2

Bagging K=5

Twin-bootstrap λ=300

+1.1 [+0.3, +1.8] +0.0 [-0.4, +0.4] -0.2 [-0.6, +0.1] +0.3 [-0.2, +0.7] -0.5 [-0.8, -0.1] -1.6 [-2.4, -0.9] +1.0 [+0.7, +1.3] +0.6 [+0.3, +0.9] +1.8 [+1.5, +2.1]

+0.2 [-0.6, +1.1] +0.1 [-0.5, +0.7] +0.3 [-0.1, +0.6] -0.3 [-0.7, +0.1] -0.5 [-1.0, +0.0] -4.6 [-5.3, -3.9] -0.4 [-0.7, -0.2] -0.6 [-0.9, -0.4] -0.5 [-0.7, -0.3]

-3.9 [-5.4, -2.5] -6.0 [-7.0, -5.1] -2.3 [-3.2, -1.5] -2.8 [-3.6, -2.0] -3.0 [-3.5, -2.5] -5.8 [-6.6, -5.0] -2.7 [-3.2, -2.3] -2.0 [-2.5, -1.4] -3.9 [-4.3, -3.6]

-8.4 [-9.8, -7.0] -7.2 [-8.2, -6.2] -4.4 [-5.1, -3.7] -6.5 [-7.2, -5.7] -6.0 [-6.5, -5.4] -10.5 [-11.2, -9.7] -4.1 [-4.5, -3.7] -4.0 [-4.5, -3.6] -7.2 [-7.6, -6.8]

-2.7 [-4.2, -1.1] -9.3 [-10.2, -8.5] -3.6 [-4.3, -3.0] -10.5 [-11.0, -9.9] -4.7 [-5.4, -4.0] -14.1 [-14.9, -13.3] -7.4 [-7.8, -7.0] -7.0 [-7.5, -6.6] -8.3 [-8.6, -7.9]

100 75 50 Review top 30%: catch 58 100% of flips

25 0

100

Test predictions reviewed, ranked by fragility (%)

BACE (dev) held-out: DILI, CYP2D6-Sub, Pgp, BBB-Martins, MOF-thermal, BBBP, TADF, AMES random ranking

Figure 3: Routing the top 30% of test predictions (ranked by per-example churn from one extra retraining) to a human reviewer captures 58–100% of all retraining-induced class flips on the 9 chemistry datasets. Each curve is one dataset (BACE in red, held-out in grey); x-axis is the fraction of test predictions reviewed in churn-rank order, y-axis is the cumulative fraction of total flip-mass captured. The diagonal marks random ranking; the shaded band marks the practical review fraction. The triage workflow is a one-bootstrap operation on top of an existing pipeline. Pretrained backbones. The λ=300 chosen on the BACE MLP over-constrains pretrained representations: twin-bootstrap collapses accuracy by 26pp on Waterbirds (ImageNet ResNet-50) and by 11–18pp on ChemBERTa-77M-MTR. Re-applying the 0.02-tolerance selection rule on each backbone’s development set picks λ=10, which preserves accuracy and cuts the class-flip rate on every ChemBERTa dataset (Appendix S) and on Waterbirds (Appendix T). Bagging is independent of λ and works on every backbone we tested. Graph networks. On a 3-layer GIN on RDKit graphs of BACE, bagging-K=5 transfers cleanly and improves id-accuracy by +4.5pp. Twin-bootstrap at λ=300 collapses accuracy and is rejected by the same 0.02-tolerance rule; re-running the rule on BACE-GIN selects λ=10, which cuts the class-flip rate 52% at no accuracy cost (Appendix R). 8

Regression. On three regression benchmarks (ESOL, FreeSolv, Lipophilicity; Appendix U) the methods extend to continuous targets and rank the same way as on classification. Cross-sample churn for regression is the per-example absolute prediction difference between two retrainings; on ERM, that churn is 42–57% of the deployed model’s mean absolute error — a single retraining moves predictions by roughly half the test-set MAE on average. At matched 2×-ERM compute, twin-bootstrap cuts prediction-difference churn 39–42% vs. ERM while bagging-K=2 cuts 29–32% (twin-bootstrap wins on every dataset). Bagging-K=5 at 5×-ERM compute is strongest at 55–58%.

Limitations

The main magnitudes and method comparison are on chemistry binary classification (MoleculeNet, TDC ADME and Tox, materials-science). Regression is covered for three MoleculeNet benchmarks (Appendix U); Waterbirds is a single vision cross-check; multi-class classification and structured outputs are not tested. We use K=2 throughout twin-bootstrap and do not study larger head counts. We measure cross-sample churn at the prediction layer, at the top-10 ranking layer that gates a single Bayesian-optimisation acquisition step (Appendix P), and at the trajectory level on the three regression benchmarks (Appendix Q: twin-bootstrap reduces the cross-trajectory std of the final-best y by 34%–100% in a 10-step greedy BO loop on every dataset; bagging on 2/3). We use greedy top-1 acquisition only; UCB, expected improvement, and Thompson sampling are not tested.

Discussion

Cross-sample churn is a metric scientific-ML benchmarks have been missing. Parameter-side uncertainty methods that the field treats as adequate — deep ensembles, MC dropout, stochastic weight averaging — do not reduce cross-sample churn on these benchmarks, and on several they make it worse. This follows from how they are constructed: they sample over weights at fixed data, while the variance practitioners act on lives on the data axis they hold constant. A performance-only benchmark cannot separate them from data-side methods that do reduce the rate, because the perprediction signal they differ on is not part of the standard report. Adding cross-sample churn to the standard report exposes the methods that produce different per-prediction calls between retrainings even when aggregate accuracy looks the same — the regime where the model’s choice of next molecule is bootstrap-sensitive and the workflow built on top of it is not reproducible. The data-resampling protocol is the design lever. Cross-sample stability is set less by the model class than by how the training procedure resamples data between independent runs. The codistillationvs. twin-bootstrap comparison (Appendix I) isolates the inter-network bootstrap overlap as the operative parameter; the consistency weight λ that controls how strongly the procedure exploits the overlap takes a different numerical value on different architectures, but the development-set rule that selects it transfers unchanged across MLP, GIN, ChemBERTa, and ImageNet-pretrained ResNet-50 (Section 6.1). Cross-sample stability can therefore be designed in at training time, without changing the model class or the deployment pipeline, by tuning a single hyperparameter on the development set with the same rule. Broader impacts. Reducing cross-sample churn has operational consequences because deploymentdriven domains route substantial wet-lab effort to the molecules a model ranks highest at any given retraining; reducing the rate at which that ranking flips between dataset updates cuts wasted experimental work and makes computational triage decisions reproducible. However, low churn is not correctness: a stably wrong model looks identical to a stably right one under this metric, and the churn-ranked triage workflow only flags predictions that disagree across retrainings, not those that agree but are jointly miscalibrated. Cross-sample churn is also a bootstrap-variance quantity — it does not characterise stability under distribution shift between the training population and the deployment domain, which requires separate measurement.

Acknowledgements K.M.J.’s work is supported by the Carl Zeiss Foundation. G.P.’s work was supported by the HPC Gateway measure of the Helmholtz Association. The authors thank Martiño Ríos-García for feedback on an early draft of the manuscript.

Declaration of generative AI and AI-assisted technologies in the research and writing process In addition to using Anthropic’s Claude models to write this declaration, we used Anthropic’s Claude models and OpenAI’s Codex as copilots during code development and cluster experiment submission, and Claude additionally to improve language and readability of the manuscript. After using these services we reviewed and edited the content as needed and take full responsibility for the content of the publication.

Code availability Code is available at https://github.com/lamalab-org/data-invariance.

References Rohan Anil, Gabriel Pereyra, Alexandre Passos, Róbert Ormándi, George E. Dahl, and Geoffrey E. Hinton. Large scale distributed neural network training through online distillation. International Conference on Learning Representations, 2018. Srinadh Bhojanapalli, Kimberly Wilber, Andreas Veit, Ankit Singh Rawat, Seungyeon Kim, Aditya Menon, and Sanjiv Kumar. On the reproducibility of neural network predictions. arXiv preprint arXiv: 2102.03349, 2021. Emily Black, Klas Leino, and Matt Fredrikson. Selective ensembles for consistent predictions. In The Tenth International Conference on Learning Representations, ICLR 2022, Virtual Event, April 2529, 2022. OpenReview.net, 2022. URL https://openreview.net/forum?id=HfUyCRBeQc. Leo Breiman. Bagging predictors. Machine Learning, 24(2):123–140, August 1996. ISSN 1573-0565. doi: 10.1007/bf00058655. URL http://dx.doi.org/10.1007/BF00058655. Leo Breiman. Statistical modeling: The two cultures (with comments and a rejoinder by the author). Statistical Science, 16(3), August 2001. ISSN 0883-4237. doi: 10.1214/ss/1009213726. URL http://dx.doi.org/10.1214/ss/1009213726. Seyone Chithrananda, Gabriel Grand, and Bharath Ramsundar. Chemberta: Large-scale selfsupervised pretraining for molecular property prediction. arXiv preprint arXiv: 2010.09885, 2020. Q. Cormier, M. Milani Fard, K. Canini, and M. R. Gupta. Launch and iterate: reducing prediction churn. In Proceedings of the 30th International Conference on Neural Information Processing Systems, NIPS’16, page 3179–3187, Red Hook, NY, USA, 2016. Curran Associates Inc. ISBN 9781510838819. Alexander D’Amour, Katherine Heller, Dan Moldovan, Ben Adlam, Babak Alipanahi, Alex Beutel, Christina Chen, Jonathan Deaton, Jacob Eisenstein, Matthew D. Hoffman, Farhad Hormozdiari, Neil Houlsby, Shaobo Hou, Ghassen Jerfel, Alan Karthikesalingam, Mario Lucic, Yian Ma, Cory McLean, Diana Mincu, Akinori Mitani, Andrea Montanari, Zachary Nado, Vivek Natarajan, Christopher Nielson, Thomas F. Osborne, Rajiv Raman, Kim Ramasamy, Rory Sayres, Jessica Schrouff, Martin Seneviratne, Shannon Sequeira, Harini Suresh, Victor Veitch, Max Vladymyrov, Xuezhi Wang, Kellie Webster, Steve Yadlowsky, Taedong Yun, Xiaohua Zhai, and D. Sculley. Underspecification presents challenges for credibility in modern machine learning. arXiv preprint arXiv: 2011.03395, 2020.

Y. Gal and Zoubin Ghahramani. Dropout as a bayesian approximation: Representing model uncertainty in deep learning. International Conference on Machine Learning, 2015. Dingyun Huang and Jacqueline M. Cole. A database of thermally activated delayed fluorescent molecules auto-generated from scientific literature with chemdataextractor. Scientific Data, 11(1), January 2024. ISSN 2052-4463. doi: 10.1038/s41597-023-02897-3. URL http://dx.doi.org/ 10.1038/s41597-023-02897-3. Kexin Huang, Tianfan Fu, Wenhao Gao, Yue Zhao, Yusuf H. Roohani, J. Leskovec, Connor W. Coley, Cao Xiao, Jimeng Sun, and M. Zitnik. Therapeutics data commons: Machine learning datasets and tasks for drug discovery and development. NeurIPS Datasets and Benchmarks, 2021. Pavel Izmailov, Dmitrii Podoprikhin, Timur Garipov, Dmitry Vetrov, and Andrew Gordon Wilson. Averaging weights leads to wider optima and better generalization. In Conference on Uncertainty in Artificial Intelligence (UAI), 2018. Kevin Maik Jablonka. Clever materials: When models identify good materials for the wrong reasons. arXiv preprint arXiv: 2602.17730, 2026. Heinrich Jiang, Harikrishna Narasimhan, Dara Bahri, Andrew Cotter, and Afshin Rostamizadeh. Churn reduction via distillation. In The Tenth International Conference on Learning Representations, ICLR 2022, Virtual Event, April 25-29, 2022. OpenReview.net, 2022. URL https://openreview.net/forum?id=HbtFCX2PLq0. Balaji Lakshminarayanan, Alexander Pritzel, and Charles Blundell. Simple and scalable predictive uncertainty estimation using deep ensembles. Advances in neural information processing systems, 30, 2017. Ilya Loshchilov and Frank Hutter. Decoupled weight decay regularization. In International Conference on Learning Representations, 2019. URL https://openreview.net/forum?id= Bkg6RiCqY7. Charles T. Marx, Flávio P. Calmon, and Berk Ustun. Predictive multiplicity in classification. In Proceedings of the 37th International Conference on Machine Learning, ICML 2020, 13-18 July 2020, Virtual Event, Proceedings of Machine Learning Research, pages 6765–6774. PMLR, 2020. URL http://proceedings.mlr.press/v119/marx20a.html. Seyed Mohamad Moosavi, Aditya Nandy, Kevin Maik Jablonka, Daniele Ongari, Jon Paul Janet, Peter G. Boyd, Yongjin Lee, Berend Smit, and Heather J. Kulik. Understanding the diversity of the metal-organic framework ecosystem. Nature Communications, 11(1), August 2020. ISSN 2041-1723. doi: 10.1038/s41467-020-17755-8. URL http://dx.doi.org/10.1038/ s41467-020-17755-8. Aditya Nandy, Gianmarco Terrones, Naveen Arunachalam, Chenru Duan, David W. Kastner, and Heather J. Kulik. Mofsimplify, machine learning models with extracted stability data of three thousand metal–organic frameworks. Scientific Data, 9(1), March 2022. ISSN 2052-4463. doi: 10.1038/s41597-022-01181-0. URL http://dx.doi.org/10.1038/s41597-022-01181-0. Richard D. Riley and Gary S. Collins. Stability of clinical prediction models developed using statistical or machine learning methods. Biometrical Journal, 65(8), 2023. ISSN 1521-4036. doi: 10.1002/bimj.202200302. URL http://dx.doi.org/10.1002/bimj.202200302. Zhenqin Wu, Bharath Ramsundar, Evan N. Feinberg, Joseph Gomes, Caleb Geniesse, Aneesh S. Pappu, Karl Leswing, and Vijay Pande. Moleculenet: a benchmark for molecular machine learning. Chemical Science, 9(2):513–530, 2018. ISSN 2041-6539. doi: 10.1039/c7sc02664a. URL http://dx.doi.org/10.1039/C7SC02664A. Keyulu Xu, Weihua Hu, Jure Leskovec, and Stefanie Jegelka. How powerful are graph neural networks? In 7th International Conference on Learning Representations, ICLR 2019, New Orleans, LA, USA, May 6-9, 2019. OpenReview.net, 2019. URL https://openreview.net/forum?id= ryGs6iA5Km.

Appendix Appendix roadmap. Appendices A–E support the measurement framework (N -scaling, per-class breakdown, aggregate-metric drift, borderline datasets, and the ERM-vs-majority filter outcomes). Appendix F reports the canonical-seed sensitivity of every magnitude in the main table. Appendices G–J support the methods (Pareto curve, twin-bootstrap bayesian optimisation, overlap spectrum, ablations). Appendices K–Q support the empirical results (per-molecule case study, distributional disagreement, Friedman test, triage convergence, entropy baseline, top-K ranking stability and a full BO-loop simulation on regression as Bayesian-optimisation analogues). Appendices R–U document the cross-architecture and cross-task generalisation evidence. Appendix V reports the per-method compute footprint.

N -scaling on BACE

The cross-dataset variation in Table 1 confounds training-set size with task difficulty. We isolate N by subsampling a single dataset (BACE) to M ∈ {200, 300, 400, 500, 600, 700, 800, 900, 968} and re-running ERM at each M for 10 train-seeds; the pool at each M is the deterministic prefix of a fixed canonical-seed shuffle, so all 10 train-seeds at the same M draw their bootstraps from the same underlying pool. Sym-KL trends downward as M grows: the seed-averaged value falls from 0.83 at M =200 to 0.62 at M =800 before ticking back up slightly to 0.69 at M =968 (Table 3; log-log slope −0.20, computed on canonical-seed 99 alone for protocol consistency). More data reduces churn but does not drive it to zero, motivating the methods we propose for the regime where simply collecting more data is not available. Table 3: Cross-sample churn decreases with training-pool size on BACE, with a visible plateau between M =600 and M =800. Within-dataset N -scaling under the cross-sample bootstrap protocol. Each row is mean cross-pair sym-KL and argmax churn at training-pool size M , computed over all 45 pairs of 10 retrainings with paired-bootstrap 95% CIs (10,000 resamples). The within-dataset log-log slope of sym-KL vs M is −0.20; the plateau reflects an irreducible boundary-disagreement floor (churn decreases with N but not to zero). M

sym-KL (mean [95% CI])

Argmax churn (%, [95% CI])

200 300 400 500 600 700 800 900 968

0.886 [0.839, 0.934] 0.838 [0.792, 0.885] 0.796 [0.756, 0.838] 0.603 [0.565, 0.641] 0.570 [0.549, 0.592] 0.628 [0.592, 0.665] 0.549 [0.511, 0.588] 0.725 [0.681, 0.771] 0.753 [0.712, 0.796]

18.2 [17.5, 19.0] 19.6 [18.7, 20.5] 18.7 [18.1, 19.2] 15.3 [14.6, 16.1] 15.0 [14.5, 15.4] 15.4 [14.7, 16.1] 14.9 [14.2, 15.6] 17.4 [16.7, 18.2] 16.1 [15.6, 16.7]

Sensitivity to the canonical-pool draw. The canonical-seed shuffle determines which examples enter the pool at each M , so the per-M sym-KL is itself a single-seed estimate. We re-ran the entire M -grid on two additional canonical seeds (7 and 42) for 10 train-seeds each; the per-seed sym-KL at fixed M ranges by up to ∼30% in absolute terms, but the averaged trajectory preserves the broaddecreasing shape described above. The slight uptick from M =800 to M =968 is canonical-seed sensitive (visible on canonical seeds 99 and 7, absent on 42) and should be read as a soft floor near M ≈800 rather than a sharp inflection.

Per-class breakdown: minority-class predictions are more unstable

The overall cross-bootstrap argmax-churn rate is an average over the canonical id-test set. On imbalanced datasets, that average can hide class-conditional structure: are the disagreements concentrated on majority-class examples (where the model is essentially predicting the dominant prior) or on minority-class examples (where the model has to learn structure)? Table 4 restricts churn to each true-label subset. On the three most imbalanced chemistry datasets in our suite (BBB-Martins, 12

BBBP at positive fraction 0.78; CYP2D6-Sub at 0.30), minority-class predictions are 2–4× more unstable than majority-class predictions. The disagreement is therefore concentrated on the class practitioners care most about: the predicted-positive (active, permeable, substrate) calls that drive next-step decisions. Table 4: On imbalanced datasets, minority-class predictions are 2–4× more unstable across retrainings than majority-class predictions. For each chemistry dataset: overall cross-bootstrap argmax-churn and its restriction to y=0 and y=1 subsets of the canonical id-test, across the 10 2 = 45 ERM seed pairs (mean [ 95% CI ], 10,000 resamples). Pos. frac. is the fraction of y=1 examples; the minority class is bolded. On the most imbalanced datasets (BBB-Martins, BBBP at 0.78 pos-frac; CYP2D6-Sub at 0.30) the minority-class churn rate is 2–4× the majority-class rate, so per-example disagreement is concentrated on the rarer class — exactly the predictions practitioners care most about (active toxicity, BBB-permeable, substrate). On balanced datasets the rates are comparable. Dataset DILI CYP2D6-Sub Pgp BACE (dev) TADF MOF-thermal BBB-Martins BBBP AMES

Nid-test

Pos. frac.

Overall (%)

churn|y=0 (%)

churn|y=1 (%)

304 427 780 968 1007 1251 1300 1305 4658

76 106 194 242 428 627 324 326 1164

0.47 0.30 0.53 0.44 0.50 0.48 0.78 0.78 0.51

16.8 [15.8, 17.7] 13.3 [12.6, 14.1] 10.3 [9.7, 10.9] 16.1 [15.6, 16.7] 12.7 [12.3, 13.1] 21.8 [21.2, 22.4] 8.0 [7.6, 8.4] 8.5 [8.1, 8.9] 15.2 [14.8, 15.5]

13.6 [12.3, 14.9] 8.2 [7.6, 8.9] 13.5 [12.5, 14.6] 15.8 [15.2, 16.4] 13.5 [12.9, 14.1] 20.7 [19.9, 21.4] 18.0 [16.4, 19.7] 15.7 [14.5, 17.0] 13.9 [13.5, 14.3]

20.3 [18.7, 21.9] 25.1 [22.8, 27.4] 7.5 [6.9, 8.0] 16.6 [15.6, 17.5] 11.9 [11.3, 12.4] 23.0 [22.1, 23.9] 5.2 [4.9, 5.5] 6.4 [5.9, 6.8] 16.4 [16.0, 16.8]

Aggregate-metric drift under precision, recall, F1 , AP

Table 1 reports |∆accuracy| between two retrainings. Accuracy on imbalanced binary tasks can mask class-specific drift, so we additionally compute paired |∆| for precision, recall, F1 , and average precision (AP) on each of the 9 chemistry datasets, using the same 45 seed pairs of ERM bootstraps. Table 5 reports the result. On the imbalanced ADME datasets (CYP2D6-Sub, DILI), recall and precision drifts are 2–3× accuracy drift, but the per-example argmax-churn rate (13.3% on CYP2D6Sub, 16.8% on DILI) still exceeds every aggregate-metric drift on every dataset. Table 5: Per-example argmax-disagreement (right column) dominates aggregate-metric drift on every dataset, regardless of which summary statistic the aggregate-metric column uses. Paired |∆| of five aggregate metrics and the per-example argmax-churn rate, computed across the 10 = 45 seed pairs of ERM bootstraps; mean [ 95% CI ] over 10,000 resamples, in percentage 2 points. Per-example argmax-churn ranges 8.0–21.8%; the strongest aggregate-metric drift on any cell is |∆recall| = 10.5 pp on the imbalanced CYP2D6-Sub dataset. Aggregate-metric drift (pp) Dataset DILI CYP2D6-Sub Pgp BACE (dev) TADF MOF-thermal BBB-Martins BBBP AMES

N |∆acc| 304 4.1 [3.4, 4.8] 427 4.2 [3.4, 4.9] 780 2.0 [1.6, 2.4] 968 1.8 [1.5, 2.2] 1007 1.3 [1.1, 1.5] 1251 1.5 [1.2, 1.9] 1300 1.7 [1.4, 2.1] 1305 1.5 [1.1, 1.8] 4658 1.6 [1.3, 1.9]

Per-example

|∆prec|

|∆rec|

|∆F1 |

|∆AP|

argmax churn (%)

5.4 [4.4, 6.4] 9.9 [8.2, 11.7] 3.6 [2.9, 4.4] 2.6 [2.1, 3.2] 2.1 [1.7, 2.5] 1.6 [1.3, 1.9] 1.4 [1.1, 1.7] 0.9 [0.7, 1.0] 1.8 [1.5, 2.2]

7.1 [5.6, 8.6] 10.5 [8.5, 12.6] 2.2 [1.7, 2.7] 0.5 [0.4, 0.7] 1.5 [1.2, 1.8] 4.8 [3.9, 5.8] 1.6 [1.3, 1.9] 2.0 [1.7, 2.4] 1.4 [1.2, 1.7]

4.2 [3.4, 4.9] 6.1 [4.9, 7.4] 2.0 [1.6, 2.4] 1.8 [1.4, 2.1] 1.3 [1.1, 1.5] 1.5 [1.2, 1.9] 3.0 [2.4, 3.6] 2.0 [1.6, 2.4] 1.6 [1.3, 2.0]

3.2 [2.6, 3.9] 7.4 [6.0, 8.7] 1.1 [0.9, 1.4] 2.7 [2.2, 3.3] 0.9 [0.7, 1.2] 2.4 [1.9, 2.9] 0.9 [0.7, 1.1] 1.6 [1.3, 1.9] 0.9 [0.7, 1.0]

16.8 [15.8, 17.7] 13.3 [12.6, 14.1] 10.3 [9.7, 10.9] 16.1 [15.6, 16.7] 12.7 [12.3, 13.1] 21.8 [21.2, 22.4] 8.0 [7.6, 8.4] 8.5 [8.1, 8.9] 15.2 [14.8, 15.5]

Borderline datasets

hERG, HIA_Hou, Skin_Reaction. The headline analysis requires ERM to clear majority-class accuracy by at least 5 pp, otherwise cross-sample churn would conflate the model shifting its decision boundary with the majority class itself shuffling under sampling noise. These three datasets pass that 13

filter only marginally (+3 to +4 pp on test sets of 57–104 examples). Their cross-sample magnitudes are reported below for transparency; the method comparison is not run on them because the small test sets do not give enough statistical power. Table 6: Cross-sample magnitudes for the three borderline datasets. These pass the ERM-vsmajority filter only marginally (+3 to +4 pp on test sets of 57–104 examples) and are reported here for transparency; method comparisons are not run on them because the small test sets do not give enough statistical power. Columns and CI conventions match Table 1. Aggregate accuracy

Per-prediction disagreement

Dataset

Ntrain

Nid-test

ERM id-acc

|∆acc| (pp)

Argmax churn (%)

Sym-KL (nats)

SkinReact HIA hERG

232 370 420

57 92 104

0.675 [0.596, 0.737] 0.908 [0.891, 0.924] 0.751 [0.702, 0.798]

4.3 [3.4, 5.3] 1.0 [0.8, 1.3] 3.1 [2.5, 3.8]

20.3 [18.8, 21.8] 2.1 [1.8, 2.5] 15.9 [14.9, 16.9]

1.035 [0.931, 1.140] 0.106 [0.086, 0.127] 0.724 [0.668, 0.779]

Filter outcomes for excluded datasets

Table 7 reports the per-dataset ERM id-acc and the canonical majority-class baseline on the five chemistry datasets that fail the +5 pp filter and are therefore excluded from the main analysis. Table 7: Five datasets fail the +5pp ERM-vs-majority filter and are excluded from the main analysis. ERM id-acc is the mean across 10 retrainings; majority is the largest class proportion on the canonical id-test set. The filter requires ERM to exceed majority by at least 5pp; on each of the five rows below it does not, so cross-sample churn would conflate “method shifts the decision boundary” with “majority-class shuffling under noise”. Reported here for transparency.

Dataset

Ntrain

Nid-test

Majority

ERM id-acc

Gap (pp)

CYP2C9-Sub CYP3A4-Sub ClinTox Bioavailability MOF-solvent

428 429 948 410 849

107 107 236 102 436

0.776 0.617 0.928 0.794 0.589

0.718 0.596 0.918 0.788 0.709

-5.8 -2.1 -1.0 -0.6 +11.9

Canonical-seed sensitivity

The canonical seed determines the train/test split. All single-canonical-seed estimates of cross-sample churn therefore inherit a sensitivity to which examples ended up in the test set. We replicate the entire main-table protocol on 3 canonical seeds (99, 7, 42) for each method-dataset cell (10 train-seeds each), and report the per-cell across-seed mean (Table 2, Table 1) with the per-seed values tabulated below. What changes across canonical seeds. The largest across-seed range observed for any single (dataset, method) cell is 5.6 pp on the absolute class-flip rate, attained on BACE, ERM. Magnitudes are seed-sensitive because the test set composition changes between canonical seeds. Paired comparisons — the main ∆-churn columns of Table 2 — are computed within seed and so factor out the test-set composition; the largest across-seed range of any paired ∆ in our suite is 3.9 pp (BACE, twinbootstrap λ=300), typically a small fraction of the magnitude of the reduction. What does not change. The qualitative ranking “data-side > parameter-side” (bagging and twinbootstrap reduce churn vs. ERM on every seed; deep ensemble / MC dropout / SWA do not consistently reduce churn) holds on every canonical seed individually. Wins counts and Friedman ranks computed per-seed agree with the across-seed average reported in Section 6. The 0.02-tolerance λ-selection rule on BACE selects λ=300 on every canonical seed.

Table 8: Per-canonical-seed values for the main table. Each method-dataset cell is reported on three independent canonical splits (99, 7, 42). Top number per cell: id-churn rate (%); bottom: paired ∆ id-churn vs. ERM (pp). Both columns are aggregated over the 45 paired-bootstrap seed pairs at fixed canonical seed. ERM rows have no paired-∆ entry. Canonical seed Dataset BACE (dev)

Method ERM SWA MC dropout Deep Ens. K=5 Bag K=2 Bag K=5 Twin λ=300

DILI

ERM SWA MC dropout Deep Ens. K=5 Bag K=2 Bag K=5 Twin λ=300

CYP2D6-Sub

ERM SWA MC dropout Deep Ens. K=5 Bag K=2 Bag K=5 Twin λ=300

Pgp

ERM SWA MC dropout Deep Ens. K=5 Bag K=2 Bag K=5 Twin λ=300

BBB-Martins

ERM SWA MC dropout Deep Ens. K=5 Bag K=2 Bag K=5

18.1

12.5

16.1

17.1 -1.0 18.1 -0.0 17.4 -0.8 13.3 -4.8 9.1 -9.1 7.2 -11.0

11.9 -0.6 13.0 +0.5 12.8 +0.2 8.3 -4.2 5.8 -6.7 5.5 -7.1

15.6 -0.6 16.4 +0.3 15.8 -0.3 13.4 -2.8 9.7 -6.5 5.7 -10.5

16.5

17.1

16.8

16.7 +0.2 16.2 -0.3 16.3 -0.2 10.8 -5.7 7.8 -8.7 10.8 -5.6

17.0 -0.1 16.3 -0.8 16.3 -0.8 14.6 -2.5 11.1 -6.0 11.3 -5.8

17.3 +0.5 17.8 +1.1 17.0 +0.2 12.8 -3.9 8.4 -8.4 14.1 -2.7

12.9

13.8

13.3

12.9 +0.0 14.2 +1.3 12.6 -0.3 9.7 -3.2 7.0 -5.9 6.8 -6.1

13.6 -0.2 13.9 +0.2 13.9 +0.1 9.2 -4.5 6.5 -7.2 4.5 -9.2

13.3 +0.0 13.3 +0.0 13.4 +0.1 7.3 -6.0 6.1 -7.2 4.0 -9.3

11.9

10.7

10.3

12.2 +0.3 13.8 +1.9 12.2 +0.3 7.9 -4.0 6.8 -5.1 7.3 -4.6

10.6 -0.0 10.4 -0.2 10.4 -0.3 6.9 -3.7 6.0 -4.7 6.4 -4.2

10.2 -0.1 10.1 -0.2 10.5 +0.3 7.9 -2.3 5.9 -4.4 6.7 -3.6

7.4

9.4

8.0

7.3 -0.1 7.4 -0.0 7.1 -0.4 5.0 -2.4 2.8 -4.6

9.2 -0.2 10.3 +0.9 9.1 -0.3 6.5 -2.9 5.1 -4.3

8.0 +0.0 9.0 +1.0 7.6 -0.4 5.3 -2.7 3.9 -4.1

continued on next page

Table 8 continued from previous page Canonical seed Dataset

Method Twin λ=300

MOF-thermal

ERM SWA MC dropout Deep Ens. K=5 Bag K=2 Bag K=5 Twin λ=300

BBBP

ERM SWA MC dropout Deep Ens. K=5 Bag K=2 Bag K=5 Twin λ=300

TADF

ERM SWA MC dropout Deep Ens. K=5 Bag K=2 Bag K=5 Twin λ=300

AMES

ERM SWA MC dropout Deep Ens. K=5 Bag K=2 Bag K=5 Twin λ=300

2.1 -5.4

1.8 -7.6

0.6 -7.4

21.6

22.1

21.8

18.6 -3.0 20.1 -1.5 17.5 -4.1 16.6 -5.1 11.5 -10.1 7.7 -13.9

18.1 -4.0 21.3 -0.7 18.2 -3.9 15.3 -6.8 10.3 -11.8 6.7 -15.4

16.9 -4.9 20.1 -1.6 17.2 -4.6 16.0 -5.8 11.3 -10.5 7.7 -14.1

7.9

8.6

8.5

7.8 -0.1 8.5 +0.5 8.0 +0.0 5.2 -2.7 4.1 -3.8 1.3 -6.6

8.7 +0.0 8.9 +0.3 8.4 -0.2 6.1 -2.5 3.9 -4.8 1.5 -7.2

8.3 -0.2 9.1 +0.6 7.8 -0.6 6.5 -2.0 4.4 -4.0 1.4 -7.0

13.8

12.7

13.7 -0.1 14.4 +0.6 13.0 -0.7 9.7 -4.0 6.3 -7.4 8.2 -5.6

13.4 -0.4 14.0 +0.2 13.5 -0.3 10.6 -3.2 7.4 -6.4 8.6 -5.2

12.2 -0.5 12.2 -0.5 12.2 -0.5 9.6 -3.0 6.7 -6.0 8.0 -4.7

15.8

16.1

15.2

15.6 -0.2 16.7 +0.9 14.7 -1.1 11.1 -4.7 8.2 -7.6 7.0 -8.8

15.9 -0.2 16.9 +0.8 15.6 -0.5 11.5 -4.6 8.3 -7.8 6.8 -9.3

15.1 -0.1 17.0 +1.8 14.6 -0.5 11.2 -3.9 8.0 -7.2 6.9 -8.3

Pretrained-backbone scope studies. For the pretrained-backbone scope studies (GIN, ChemBERTa, Waterbirds: Appendices R–T) we replicate ERM and the rule-selected λ on all 3 canonical seeds; the development λ-Pareto sweeps within those studies are run on canonical seed 99 only, since the goal there is architecture-feasibility evidence and the rule-selected λ does not change across seeds.

id-accuracy (%)

= 30 = 1

= 100

80 79

= 300

78 6 ERM

8 Bagging K = 5

id-churn (%) Twin-bootstrap (

{1, 3, 10, 30, 100, 300})

Figure 4: Twin-bootstrap λ=300 sits at the accuracy-preserving end of the BACE Pareto frontier. Twin-bootstrap at λ ∈ {1, 3, 10, 30, 100, 300} traces an accuracy-vs.-churn trajectory; the pre-registered selection rule (largest λ with id-acc ≥ ERM-id-acc −0.02) picks λ=300 (filled red). Bagging-K=5 (blue circle) achieves a similar accuracy at higher churn; ERM (open square, upper-right) at the highest churn. Error bars are bootstrap 95% CIs.

Development Pareto curve on BACE

A plateau of configurations λ ∈ {1, 3, 10, 30} all preserve accuracy and reduce churn from ERM’s 16.1% to a low-single-digit range; the curve drops further at higher λ until λ=300 at 5.7% within the 0.02 accuracy tolerance (Figure 4). The selected λ sits at the extreme of the feasible interval that satisfies the 0.02 accuracy tolerance, and the choice is robust to the exact tolerance value: relaxing the tolerance to 0.01 would still pick the same λ, while tightening to 0.005 would pick λ=100 (still on the knee). Accuracy at the frozen λ across held-out datasets. At the frozen λ=300, twin-bootstrap improves id-accuracy on 3/9 chemistry datasets (BACE +1.2 pp, DILI +2.2 pp, TADF +0.5 pp) and is within paired-bootstrap CIs on AMES, BBB-Martins, CYP2D6-Sub, and BBBP. Pgp shows a small loss (−1.5 pp), and MOF-thermal collapses (−7.1 pp, well outside the 0.02 selection-rule tolerance and the only dataset on which a per-dataset re-application of the rule would reject the BACE-frozen λ).

Per-dataset Bayesian optimisation of λ

The pre-registered rule (Section 5) freezes λ=300 on BACE and applies it to every other dataset. Does per-dataset tuning win on the metric the paper reports (cross-sample churn between independent retrainings)? We test this with a Bayesian optimisation that minimises the argmax disagreement between two independent twin-bootstrap ensembles on a held-out validation fold, subject to an accuracy constraint relative to the unregularised baseline – the val-side analogue of the cross-sample churn the paper reports on test. Protocol. For each (dataset, training seed), 50 BO trials over log10 λ ∈ [log10 10−3 , log10 104 ]. Each trial trains two independent twin-bootstrap pairs on the same fold-train pool with different training seeds, and measures the argmax disagreement between their two ensembles on the held-out fold (val_cross_churn). The first trial is forced to λ=0 so the unregularised baseline accuracy a0 is set before the constraint kicks in; the next 4 trials are log-uniform random, after which BO acquisition 17

uses GP-EI (Matérn-2.5 kernel) on the scalar score(λ) = −val_cross_churn(λ) − C · max 0, a0 − 0.02 − val_acc(λ) . With val_cross_churn ∈ [0, 1] and C=100, an accuracy shortfall of 0.01 subtracts 1.0 from the score, i.e. the full feasible-region range of the first term – so BO is pulled back inside the accuracy constraint within a few trials of any violation. Val scores average over k=3 cross-validation folds carved deterministically from the canonical training pool with seed 99; id_test and ood_test are never seen during selection. Algorithm 2 gives the loop; we denote a twin-bootstrap pair as θA , θB (two networks trained jointly on independent bootstraps with the sym-KL consistency loss; deployment averages their softmax outputs). Algorithm 2 Per-dataset BO of λ with cross-sample-churn objective. Stage 1 (the outer for-loop over s) runs once per (dataset, training seed) and returns λ⋆s ; stage 2 (the two trailing lines) takes the median across seeds and retrains. Require: fold-train indices {If }kf =1 and val loaders {Vf }kf =1 ; trial count T ; init-trial count T0 ; accuracy tolerance δ; penalty C; training seeds {s1 , . . . , sS } 1: for s ∈ {s1 , . . . , sS } do ▷ stage 1: per-seed BO 2: λ1 ← 0 ▷ forced baseline trial 3: for t = 1, . . . , T do 4: ā, c̄ ← 0, 0 5: for f = 1, . . . , k do 6: (θA , θB ) ← twin-bootstrap on If at λt ′ ′ 7: (θA , θB ) ← twin-bootstrap on If at λt with an independent training seed 1 8: p ← 2 (fθA + fθB ), p′ ← 21 (fθA′ + fθB′ ) on Vf 9: ā += acc(p, Vf )/k 10: c̄ += Prx∈Vf [arg max p(x) ̸= arg max p′ (x)]/k 11: end for 12: if t = 1 then a0 ← ā 13: end if 14: scoret ← −c̄ − C · max(0, a0 − δ − ā) 15: λt+1 ← log-uniform random if t ≤ T0 , else GP-EI argmax over {(log10 λi , scorei )}i≤t 16: end for 17: λ⋆s ← λarg maxt scoret 18: end for ⋆ 19: λdataset ← median λ⋆ ▷ stage 2: aggregate s1 , . . . , λsS 20: for s ∈ {s1 , . . . , sS } do 21: retrain twin-bootstrap on the full canonical pool at λdataset with training seed s 22: end for 23: return the S retrained twins Cross-sample churn between any two retrainings within a dataset is then fixed-λ variance on both arms of the comparison vs. the frozen λ=300 baseline. The full protocol (per-seed BO sweep, median aggregation, fixed-λ retraining, and table regeneration) is reproducible with the repository’s make bayes-twin target. Findings. The per-dataset tuned λ varies by roughly two orders of magnitude and lands on both sides of 300 (Table 9, “Mean λ” column). On cross-sample ID-churn – the metric Table 2 reports – per-dataset BO beats the rule on 4/9 datasets (CI excludes zero), ties on 2/9, and loses on 3/9. Against bootstrap-matched ERM (Table 2), per-dataset BO cuts cross-sample ID-churn on every dataset, so per-dataset tuning preserves and modestly extends the twin-bootstrap effect rather than erasing it. The frozen λ=300 is therefore a robust default – it matches or beats per-dataset BO on 5/9 datasets – but not the per-dataset optimum everywhere: per-dataset tuning is worth the BO cost when the dataset’s class structure or size diverges from BACE.

Table 9: Twin-bootstrap bayesian optimization discovers dataset-specific operating points on the accuracy–stability frontier. Cells report the first method named in each block minus the second with paired bootstrap 95% confidence intervals. Positive accuracy deltas indicate improved performance, while negative churn deltas indicate improved stability. Bayesian optimization (BO) generally identifies intermediate λ values that improve ID accuracy and often improve OOD accuracy relative to λ=0, while trading off some stability relative to strongly regularized models (λ=300). The BO-vs.-ERM block checks whether optimized λ values retain a churn reduction relative to ordinary retraining. The λ column reports the mean optimized λ for the first method in each comparison block. Bolded entries indicate cases where the point estimate favors the first method. BO vs. λ=300

Dataset

N Mean λ ID acc

BACE (dev) 968 93.03 DILI 304 4.78 CYP2D6-Sub 427 256.3 Pgp 780 84.66 BBB-Martins 1300 362.3 MOF-thermal 1251 14.91 BBBP 1305 373.7 TADF 1007 74.6 AMES 4658 424.1

OOD acc

ID churn

OOD churn

+1.3 [+0.5,+2.2] +1.6 [+0.8,+2.4] +2.1 [+1.7,+2.4] -0.2 [-0.6,+0.2] -1.1 [-3.2,+0.8] +6.0 [+3.9,+8.0] -3.6 [-5.0,-2.3] -1.7 [-2.7,-0.7] +0.6 [+0.2,+1.0] +0.5 [+0.2,+0.8] -0.9 [-1.2,-0.6] +0.3 [-0.0,+0.7] +2.8 [+1.9,+3.8] +4.1 [+2.8,+5.5] -0.5 [-1.1,+0.0] +0.5 [-0.1,+1.1] -0.0 [-0.2,+0.2] -0.2 [-0.5,+0.1] +0.1 [+0.0,+0.3] +0.1 [-0.0,+0.3] +8.1 [+7.3,+8.8] +8.1 [+7.3,+8.8] +4.8 [+4.3,+5.4] +4.8 [+4.3,+5.4] -0.3 [-0.6,-0.1] -0.4 [-0.8,-0.0] -0.1 [-0.2,+0.1] -0.3 [-0.6,+0.0] +0.9 [+0.1,+1.6] +0.9 [+0.1,+1.6] -0.8 [-1.4,-0.3] -0.8 [-1.4,-0.3] -0.7 [-0.9,-0.5] -0.2 [-0.4,+0.1] -0.7 [-0.9,-0.5] -0.2 [-0.4,-0.1] BO vs. ERM

Dataset

N Mean λ ID acc

BACE (dev) 968 93.03 DILI 304 4.78 CYP2D6-Sub 427 256.3 Pgp 780 84.66 BBB-Martins 1300 362.3 MOF-thermal 1251 14.91 BBBP 1305 373.7 TADF 1007 74.6 AMES 4658 424.1

OOD acc

ID churn

OOD churn

+2.0 [+0.9,+3.1] +4.7 [+3.3,+6.1] -8.4 [-9.0,-7.9] -6.7 [-7.4,-6.1] +2.4 [+0.5,+4.2] +1.7 [-0.0,+3.3] -6.3 [-7.6,-5.0] -4.4 [-5.4,-3.5] +0.8 [-1.0,+2.4] +1.5 [+0.2,+2.9] -10.2 [-11.1,-9.4] -9.4 [-10.5,-8.4] -0.2 [-0.9,+0.6] -1.5 [-2.3,-0.8] -4.2 [-4.8,-3.6] -3.2 [-3.9,-2.6] +0.4 [-0.5,+1.1] +1.4 [+0.9,+1.9] -7.2 [-7.7,-6.8] -11.7 [-12.1,-11.4] +0.3 [-0.7,+1.2] +0.3 [-0.7,+1.2] -9.2 [-9.8,-8.6] -9.2 [-9.8,-8.6] +1.7 [+1.0,+2.5] +2.8 [+1.3,+4.3] -7.1 [-7.5,-6.7] -11.6 [-12.2,-11.0] +0.9 [+0.0,+1.8] +0.9 [+0.0,+1.8] -5.5 [-6.0,-5.0] -5.5 [-6.0,-5.0] -3.1 [-4.1,-2.2] +0.1 [-0.7,+1.0] -8.9 [-9.3,-8.6] -10.6 [-10.9,-10.3]

Overlap spectrum: ∼40% balances churn reduction and accuracy preservation

The codistillation [Anil et al., 2018], twin-bootstrap, and twin-shared variants share the same twonetwork sym-KL consistency objective and differ only in the inter-network bootstrap-overlap: 0% (disjoint shards), ∼40% (independent bootstraps with replacement), and 100% (the same bootstrap to both networks). Table 10 reports paired ∆ class-flip rate vs. ERM at each operating point on every dataset where all three methods completed. Codistillation has the largest mean churn reduction on most datasets but drops id-accuracy by more than 5 pp on Pgp and MOF-thermal; twin-shared drops accuracy on MOF-thermal and increases churn on DILI; twin-bootstrap at ∼40% matches codistillation on never increasing churn while only collapsing accuracy on a single dataset (MOF-thermal, where all three operating points collapse). Among the three points, only twin-bootstrap at ∼40% overlap satisfies both criteria simultaneously, which is what motivates its use as the default operating point.

Weight-decay ablation

Weight decay alone is insufficient. Sweeping AdamW weight decay over five orders of magnitude (10−5 to 10−1 ) on BACE, BBBP, and TADF produced no detectable movement in cross-sample churn (paired CIs vs. ERM include zero at every value tested). The effective per-step regularisation under AdamW does not constrain the data-resampling sensitivity that drives churn.

Table 10: Twin-bootstrap (∼40% overlap) reduces churn on every dataset without an accuracy collapse; codistillation (0%) collapses accuracy on two datasets; twin-shared (100%) increases churn on DILI. Cells: paired ∆ id-churn vs. ERM in percentage points, [ 95% CI ], † marks idaccuracy drop > 5 pp. Dataset SkinReact DILI HIA hERG Pgp BACE TADF MOF-thermal BBB-Martins BBBP AMES Mean across 11 datasets Acc-dagger count Churn-positive count

Codistillation 0%

Twin-bootstrap ∼40%

Twin-shared 100%

232 304 370 420 780 968 1007 1251 1300 1305 4658

-17.0 [-18.8, -15.2] -10.3 [-11.5, -9.2] -2.1 [-2.5, -1.8] -11.7 [-12.9, -10.5] -4.8 [-5.7, -4.0]† -14.1 [-14.7, -13.5] -5.4 [-6.4, -4.3] -17.7 [-18.4, -17.0]† -7.4 [-7.8, -7.0] -7.9 [-8.3, -7.5] -12.8 [-13.2, -12.5]

-15.0 [-16.6, -13.5] -2.7 [-4.2, -1.1] -1.7 [-2.1, -1.3] -9.7 [-10.7, -8.7] -3.6 [-4.3, -3.0] -10.5 [-11.0, -9.9] -4.7 [-5.4, -4.0] -14.1 [-14.9, -13.3]† -7.4 [-7.8, -7.0] -7.0 [-7.5, -6.6] -8.3 [-8.6, -7.9]

-7.4 [-9.0, -5.8] +0.4 [-1.2, +2.2] -1.2 [-1.7, -0.6] -6.8 [-8.0, -5.6] -1.4 [-2.1, -0.7] -7.1 [-7.8, -6.3] -3.4 [-4.1, -2.8] -12.5 [-13.3, -11.7]† -6.5 [-6.9, -6.1] -6.8 [-7.2, -6.4] -6.2 [-6.7, -5.8]

-10.1 2/11 0/11

-7.7 1/11 0/11

-5.3 1/11 1/11

— — —

Per-molecule case study on BACE

Figure 5 shows six BACE id-test molecules where ERM flips class on ≥ 36% of seed pairs and twin-bootstrap flips on 0% over the same ten retrainings. Both methods see the same canonical training pool and test set; the difference is the consistency loss.

ERM 56%

Twin-bootstrap 0%

ERM 53%

Twin-bootstrap 0%

ERM 53%

Twin-bootstrap 0%

ERM 47%

Twin-bootstrap 0%

ERM 47%

Twin-bootstrap 0%

ERM 36%

Twin-bootstrap 0%

class 0

class 1

Figure 5: Twin-bootstrap eliminates per-molecule retraining flips on BACE. Six id-test molecules where ERM flips class on ≥ 36% of seed pairs and twin-bootstrap flips on 0%. Below each structure: ten ERM and ten twin-bootstrap seed predictions, each cell is the predicted class at one retraining; the gap separates the two methods. ERM stripes flip; twin-bootstrap collapses to a single colour.

Distributional disagreement (sym-KL) per dataset

Table 11 reports paired ∆ sym-KL vs. ERM per dataset and per method, in the same format as the argmax-churn main result (Table 2). Aggregating over datasets, bagging-K=5 cuts sym-KL by 6–9× vs. ERM and twin-bootstrap by a further median ∼9×; the table below shows the per-method, per-dataset paired CIs behind those aggregates. Table 11: Twin-bootstrap reduces distributional disagreement (sym-KL) by an additional factor of ∼8 beyond the strongest argmax-churn reducer. Paired ∆ sym-KL vs. ERM on the canonical id-test (in nats; ∆ < 0 better). Cells show mean ∆ and the relative reduction vs. ERM in parentheses, over all 45 pairs of 10 retrainings; “∗ ” marks cells whose 95% paired-bootstrap CI excludes zero. Parameter-side Dataset DILI CYP2D6-Sub Pgp BACE (dev) TADF MOF-thermal BBB-Martins BBBP AMES

ERM 0.75 0.64 0.51 0.75 0.40 0.38 0.50 0.47 1.12

Data-side

Deep Ens. K=5

MC dropout ∗

-0.10 (-13%) -0.08 (-13%)∗ -0.03 (-5%)∗ -0.06 (-8%)∗ -0.09 (-21%)∗ -0.21 (-55%)∗ -0.06 (-11%)∗ -0.05 (-11%)∗ -0.26 (-23%)∗

∗

-0.07 (-9%) -0.03 (-5%)∗ -0.02 (-5%)∗ -0.09 (-12%)∗ -0.03 (-8%)∗ -0.17 (-45%)∗ -0.05 (-10%)∗ -0.06 (-13%)∗ -0.25 (-23%)∗

Bagging K=2

Bagging K=5

∗

-0.50 (-67%) -0.41 (-65%)∗ -0.30 (-58%)∗ -0.40 (-54%)∗ -0.21 (-54%)∗ -0.21 (-56%)∗ -0.29 (-59%)∗ -0.25 (-54%)∗ -0.63 (-56%)∗

-0.67 (-89%) -0.55 (-85%)∗ -0.43 (-84%)∗ -0.64 (-85%)∗ -0.33 (-83%)∗ -0.32 (-83%)∗ -0.43 (-86%)∗ -0.40 (-85%)∗ -0.96 (-86%)∗

Twin-bootstrap λ=300 -0.74 (-98%)∗ -0.63 (-98%)∗ -0.50 (-98%)∗ -0.74 (-98%)∗ -0.39 (-98%)∗ -0.38 (-99%)∗ -0.49 (-98%)∗ -0.46 (-98%)∗ -1.11 (-99%)∗

Friedman test: twin-bootstrap matches bagging-K=5 in mean rank

The pairwise CIs reported in Table 2 test each (method, dataset) cell separately. As an aggregate cross-dataset test, we apply the Friedman non-parametric rank test on mean cross-sample class-flip rate over the nine datasets and seven methods (ERM, MC dropout, SWA, Deep Ensemble K=5, Bagging K=2, Bagging K=5, Twin-bootstrap λ=300). The test rejects the null of equal ranks at χ2 = 44.8, p = 5.2×10−8 . Mean ranks (lower is better, 1 = best on cross-sample class-flip rate) are twin-bootstrap 1.44, bagging K=5 1.67, bagging K=2 2.89, SWA 4.94, deep-ensemble K=5 5.17, ERM 5.78, MC dropout 6.11. The Nemenyi critical difference at α=0.05 is CD=3.0; rank gaps exceeding CD identify pairs that differ significantly. Twin-bootstrap is significantly better in rank than ERM, MC dropout, SWA, and deep ensemble; it is not significantly different from bagging K=5 (rank gap 0.22). Twin-bootstrap therefore matches bagging K=5 in mean-rank at roughly 2× ERM compute (Appendix V) versus the 5× that bagging K=5 requires.

Convergence: how many bootstraps suffice for triage?

The triage workflow asks the practitioner to train one extra ERM bootstrap and rank predictions by churn. We verify that K=2 (one extra bootstrap) is sufficient by comparing the top-30% cumulative recall when the churn ranking is computed from K ∈ {2, 3, 5, 10} ERM bootstraps.

Predictive entropy: a weaker flip-predictor than churn on every dataset

A natural deployment-cheaper alternative to per-example churn is the predictive entropy of a single trained model: the practitioner already has access to a softmax distribution at deployment time and need not train a second model. We test whether the per-example entropy ranks test predictions for retraining-flip risk as well as churn does. Across the eight held-out datasets, churn captures more of the total flip-mass in its top decile and has a larger area under the precision-vs-coverage curve than predictive entropy. The gap is consistent and ranges from 1pp to 13pp on top-decile recall. The practitioner cost of computing churn is one additional ERM training run, which the speedup over training a full K-bootstrap-bagging ensemble may amortise; entropy is free at deployment time but is a meaningfully weaker signal.

Top-30% cumulative recall (%)

100 75 50 25 0

ERM bootstraps used to score fragility (K)

BACE (dev) held-out: DILI, CYP2D6-Sub, Pgp, BBB-Martins, MOF-thermal, BBBP, TADF, AMES

Figure 6: One extra bootstrap is enough: the top-30% recall at K=2 is within 10–24 pp of the K=10 gold standard on every dataset. Mean recall across 30 random K-subsets per K<10. Per-dataset at K=2: 48–83%. At K=10 (using all ten bootstraps to score churn): 58–100%. MOFthermal is the floor on every K (48–58%); BBBP and BBB-Martins reach 100% recall by K=5. Table 12: Per-example churn beats single-model predictive entropy as a retraining-flip predictor on every chemistry dataset. Churn is computed from one extra bootstrap pair on the canonical id-test. “Top-10% recall” is the fraction of all retraining-induced flips captured by the top decile of the score; “AuPC” is the area under the precision-vs-coverage curve. Higher is better for both columns. Dataset DILI CYP2D6-Sub Pgp BACE TADF MOF-thermal BBB-Martins BBBP AMES

Top-10% recall (%) Churn Entropy 30.0 38.6 50.2 33.2 40.8 24.6 59.8 59.4 35.5

29.1 29.3 39.0 26.7 29.0 21.1 46.4 46.7 27.6

AuPC Churn Entropy 0.355 0.306 0.263 0.343 0.298 0.397 0.223 0.233 0.331

0.323 0.269 0.233 0.303 0.247 0.351 0.193 0.209 0.288

Bayesian-optimisation analogue: top-K ranking stability

A virtual-screening or Bayesian-optimisation pipeline ranks the candidate pool by predicted P (active) and acquires the top-K for the next round of evaluation. The BO-relevant question is then how often the top-K set would change between two retrainings on independent bootstraps of the same training pool — the fraction of candidates that would be re-routed to different downstream evaluations depending on which bootstrap the surrogate happened to land on. Measurement. For each chemistry dataset and each method (ERM, bagging-K=5, twin-bootstrap λ=300), we use the same 10 retrainings as the main argmax-churn protocol (Section 3). The canonical seed pins the train/test split once: every retraining is trained on its own bootstrap of the canonical training pool and is then evaluated on the same canonical id-test pool, which the surrogate has never seen during training. The candidate library is therefore fixed across retrainings; only the surrogate varies (its training-data bootstrap and initialisation).

For every retraining we rank the candidate library by the deployment-time prediction P (class=1) — a single forward pass for ERM, the K=5-model averaged softmax for bagging, and the 2-network averaged softmax for twin-bootstrap — and take the top-K molecules (K=10). K is an absolute batch size, mirroring a typical single-round wet-lab acquisition; the corresponding fraction of the canonical test pool ranges from ∼1% (AMES, Ntest =1164) to ∼13% (DILI, Ntest =76). Between each pair of retrainings (s, s′ ) of a given method we compute the Jaccard overlap (s) (s′ ) | TK ∩ TK | (s,s′ ) JK = , (s) (s′ ) | TK ∪ TK |

(s) where TK is the method’s top-K set under retraining s. We report the mean JK over the 10 2 =45 pairs with paired-bootstrap 95% CIs over the seed-pair distribution and the paired ∆JK vs. ERM (positive = stabler ranking). We additionally report the ERM top-K hit rate — the fraction of the surrogate’s top-K predicted-active molecules that are actually y=1 in the held-out test set — alongside the dataset’s positive-class prior. The class prior is the chance baseline for hit rate (a random ranker scores its positive fraction); ERM hit rate substantially above the prior confirms the surrogate has non-trivial signal and the Jaccard comparison is over informative top-K sets. This measures the static input-stability of a single BO acquisition step under the same candidate library, not the dynamics of an iterative BO loop where two divergent trajectories would expose their surrogates to different candidate pools downstream. Table 13: Top-K ranking stability (10 molecules) between independent retrainings: a Bayesianoptimisation analogue of cross-sample churn. Jaccard overlap of the top-K predicted-active sets across the same 45 pairs of 10 retrainings as the main table; 1.0 = identical sets, 0.0 = disjoint. Paired ∆ vs. ERM in the right two columns (positive = stabler ranking). All cells: mean [ 95% CI ] from 10,000 paired-bootstrap resamples. Class prior is the positive-class fraction on the canonical id-test, the chance baseline for top-K hit rate; the surrogate has signal whenever the ERM hit rate exceeds the prior. The Jaccard difference is the BO-relevant consequence of cross-sample churn.

Dataset BACE (dev) DILI CYP2D6-Sub Pgp BBB-Martins MOF-thermal BBBP TADF AMES

class

ERM

prior (%)

hit rate (%)

Jaccard

44 47 30 53 78 48 78 50 51

Bagging-K=5 Jaccard

92 0.24 [0.21, 0.26] 0.49 [0.45, 0.52] 90 0.56 [0.53, 0.59] 0.67 [0.64, 0.70] 67 0.44 [0.40, 0.48] 0.57 [0.54, 0.61] 100 0.24 [0.20, 0.27] 0.23 [0.19, 0.26] 96 0.08 [0.05, 0.11] 0.15 [0.12, 0.18] 73 0.33 [0.29, 0.36] 0.42 [0.37, 0.47] 98 0.11 [0.08, 0.15] 0.15 [0.12, 0.19] 99 0.24 [0.21, 0.26] 0.33 [0.29, 0.37] 89 0.03 [0.02, 0.05] 0.09 [0.06, 0.13]

∆ vs ERM

Twin-bootstrap λ=300 Jaccard

+0.25 [+0.20, +0.29] 0.68 [0.65, 0.71] +0.11 [+0.06, +0.16] 0.69 [0.65, 0.73] +0.14 [+0.09, +0.18] 0.65 [0.62, 0.68] -0.01 [-0.06, +0.04] 0.63 [0.59, 0.67] +0.07 [+0.03, +0.10] 0.65 [0.60, 0.69] +0.09 [+0.03, +0.16] 0.45 [0.40, 0.50] +0.04 [-0.01, +0.08] 0.64 [0.61, 0.67] +0.09 [+0.05, +0.14] 0.46 [0.41, 0.51] +0.06 [+0.02, +0.10] 0.49 [0.45, 0.54]

∆ vs ERM +0.44 [+0.40, +0.49] +0.13 [+0.08, +0.18] +0.21 [+0.17, +0.26] +0.39 [+0.35, +0.44] +0.56 [+0.51, +0.62] +0.12 [+0.06, +0.19] +0.53 [+0.47, +0.58] +0.22 [+0.16, +0.28] +0.46 [+0.42, +0.51]

Findings. ERM top-10 Jaccard is low across the 9 chemistry datasets, ranging from 0.03 (AMES) to 0.56 (DILI). On AMES, BBB-Martins, and BBBP two retrainings on independent bootstraps overlap on a single-digit fraction of their top-10 picks: a BO loop on these datasets would route nearly every wet-lab evaluation to a different molecule depending on the bootstrap. Twin-bootstrap raises Jaccard on every dataset (paired ∆J10 ∈ [+0.12, +0.56], lifting the lower end of the range to 0.45 and the upper end to 0.69), with the largest absolute gains where ERM is most unstable. The Jaccard gain exceeds what the argmax churn reduction alone would predict: top-K ranking depends on the full softmax distribution, not just the argmax, and that is the quantity sym-KL — and twin-bootstrap — stabilise (Section 6, “Distributional disagreement”). Scope. The Jaccard quantity above measures input-stability of a single BO acquisition step under the identical candidate library; it does not simulate the full BO loop, where two trajectories that acquire different molecules at step 1 subsequently see different candidate pools and may diverge further. What the Jaccard does establish is the necessary condition for any BO loop on these datasets to make a reproducible first-batch decision: with ERM, two equally plausible bootstrap surrogates would route the next wet-lab batch to almost disjoint molecule sets on AMES, BBB-Martins, and BBBP; with twin-bootstrap, the routing is largely bootstrap-invariant. Appendix Q extends the test to a full BO loop on the regression benchmarks.

Q BO trajectory variance on regression: cross-sample churn drives outcome and identity spread Section P measures the input-stability of a single BO acquisition step on classification. This appendix complements that result with a full BO-loop simulation on the three regression benchmarks (ESOL, FreeSolv, Lipo): are independent BO trajectories of the same method more reproducible at the trajectory level under bagging-K=5 and twin-bootstrap than under ERM? Protocol. For each regression dataset and each of {ERM, bagging-K=5, twin-bootstrap λ=3} we run T =20 independent BO trajectories (Algorithm 3). All T trajectories of a given method share the same random initial labelled subset of 50 molecules; trajectories diverge only in the in-loop training-data bootstraps, the same source of variance the measured elsewhere. Acquisition is greedy top-1 maximisation of the predicted regression target (ESOL log-solubility, FreeSolv hydration free energy, Lipo log-octanol/water partition); higher y is the BO objective on all three datasets and “final-best y” below refers to the largest y acquired along a trajectory. ERM has no natural per-prediction variance estimate, while bagging and twin-bootstrap do (head-spread); implementing UCB / EI / Thompson acquisition asymmetrically across methods would not be a fair comparison, so we use greedy for every method. The budget is 10 acquisitions on top of the initial subset; the candidate library is the canonical training pool with the oracle y values held back until acquisition. Algorithm 3 BO trajectory (one trajectory of one method on one dataset). Require: candidate pool (X, y), training method M, trajectory index k, initial labelled subset L0 , budget B 1: L ← L0 2: for t = 1, . . . , B do 3: seed ← k · 106 + t ▷ step- and trajectory-keyed bootstrap seed 4: f ← M.train(XL , yL ; seed) ▷ ERM / bagging-K=5 / twin-bootstrap 5: ŷ ← f.predict(X[ 1..N ]\L ) 6: i∗ ← arg max ŷ ▷ greedy top-1 acquisition 7: L ← L ∪ { i∗ } ▷ oracle reveals yi∗ 8: end for 9: return L \ L0 ▷ the B acquired molecules of trajectory k

Table 14: BO trajectory variance on the three regression benchmarks: bagging-K=5 and twin-bootstrap reduce the cross-trajectory standard deviation of the final-best y on every dataset. For each (dataset, method) we run T =20 BO trajectories sharing the same random initial subset of 50 labelled molecules; trajectories diverge only in the in-loop training-data bootstraps. At each step the surrogate is retrained from scratch, predicts ŷ on the unlabelled remainder, and acquires the arg max ŷ. Final best reports cross-trajectory mean and std with 95% bootstrap CIs over the T trajectories (10,000 resamples). std/range is the std as a percentage of each dataset’s y range, anchoring its absolute scale. Acquired Jaccard is the mean overlap of per-trajectory acquiredmolecule sets across all T2 trajectory pairs. Dataset

Method

ESOL

ERM Bagging-K=5 Twin-λ=3 ERM Bagging-K=5 Twin-λ=3 ERM Bagging-K=5 Twin-λ=3

FreeSolv

Lipo

Final best mean [95% CI]

Final best std [95% CI]

std/range (%)

Acquired Jaccard

1.50 [1.43, 1.57] 1.57 [1.57, 1.57] 1.57 [1.57, 1.57] 2.50 [2.27, 2.74] 3.11 [3.03, 3.18] 2.74 [2.58, 2.88] 3.94 [3.91, 3.98] 3.95 [3.92, 3.99] 3.91 [3.90, 3.93]

0.167 [0.000, 0.218] 0.000 [0.000, 0.000] 0.000 [0.000, 0.000] 0.547 [0.406, 0.665] 0.178 [0.042, 0.268] 0.362 [0.260, 0.425] 0.082 [0.045, 0.101] 0.089 [0.062, 0.102] 0.045 [0.000, 0.073]

1.3 0.0 0.0 2.5 0.8 1.7 1.4 1.5 0.7

0.48 [0.47, 0.50] 0.65 [0.63, 0.67] 0.61 [0.59, 0.62] 0.38 [0.36, 0.40] 0.45 [0.43, 0.47] 0.39 [0.37, 0.42] 0.37 [0.36, 0.39] 0.52 [0.50, 0.53] 0.47 [0.45, 0.49]

Findings. Twin-bootstrap reduces the cross-trajectory standard deviation of the final-best y on every regression dataset (Table 14): the per-dataset std drops from 0.082–0.547 under ERM to 0.000–0.362 under twin-bootstrap, a relative reduction of 34%–100%. Bagging-K=5 reduces std on 2/3 datasets 24

(down to 0.000–0.178 across the three) but does not improve on ERM on the third (Lipo), where ERM’s std is already small (0.082). In absolute terms the std is 1.27%–2.53% of each dataset’s y range under ERM. What matters for BO reproducibility is whether independent trajectories converge to the same molecule, not whether their final-y spread is fractionally small: on ESOL every bagging and twin-bootstrap trajectory converges to the same final molecule (std = 0), while ERM trajectories scatter across final-best y values from 1.10 to 1.57. The identity-level Jaccard of acquired-molecule sequences also rises under both methods: mean 0.37– 0.48 under ERM, 0.45–0.65 under bagging, 0.39–0.61 under twin-bootstrap. The spread reduction therefore reflects more consistent molecule selection, not just lucky convergence to molecules with similar y values. Caveats. Greedy top-1 acquisition isolates surrogate stability but is narrower than the UCB / EI / Thompson policies that real BO campaigns typically run; whether the reproducibility gain transfers to those acquisition functions is untested. Cells in the table report 95% bootstrap CIs over the K=20 trajectories (10,000 resamples), but K is small for tight std-of-std estimates — the std CIs widen accordingly. We treat this experiment as a focused check that the cross-sample-churn reduction translates to BO-loop reproducibility on small-N regression, not a full benchmark of BO algorithms.

GIN: bagging transfers; rule selects λ=10

The main comparison uses a 256-unit MLP on Morgan fingerprints (2048 bits, radius 2). As a non-fingerprint architecture cross-check, we train a graph isomorphism network (GIN, Xu et al., 2019) on RDKit-derived atom-and-bond graphs of BACE: 3 GINConv layers, hidden dimension 128, mean-pool readout, 50 epochs. All ten train-seeds and the canonical-data-seed protocol are unchanged. Table 15 reports per-method id-accuracy, class-flip rate, and sym-KL with 95% paired-bootstrap CIs, plus paired ∆ vs. ERM on the same 45 seed-pairs. Table 15: GIN on BACE: bagging transfers cleanly at the same K=5; twin-bootstrap requires re-running the λ-selection rule on the GIN backbone. ERM-GIN is more fragile than ERM-MLP (23.0% vs. 16.1% argmax churn), making the methods more rather than less relevant on this backbone. Bagging-K=5 cuts churn 54% and improves id-accuracy by +4.5pp. Twin-bootstrap at the λ=300 chosen on the BACE MLP reduces sym-KL by −99% but drops id-accuracy by 16pp — well outside the 0.02 selection-rule tolerance ERM-GIN id-acc would impose. Bold cells mark the best mean per column among the three methods. Method ERM Bagging-K=5 Twin-bootstrap λ=300

id-acc

id-churn (%)

sym-KL

0.742 [0.725, 0.757] 0.786 [0.778, 0.794] 0.582 [0.568, 0.599]

23.0 [21.9, 24.2] 10.6 [10.2, 11.0] 7.1 [5.1, 9.1]

0.486 [0.451, 0.524] 0.068 [0.065, 0.070] 0.004 [0.003, 0.005]

-12.4 [-13.6, -11.2] (−54%) -16.0 [-17.9, -14.0] (−69%)

-0.418 [-0.456, -0.383] (−86%) -0.482 [-0.519, -0.447] (−99%)

Paired ∆ vs. ERM (same 45 seed-pairs) Bagging-K=5 +4.46 pp Twin-bootstrap λ=300 −15.95 pp

The consistency loss is active and suppresses inter-network disagreement on GIN (99% sym-KL reduction), but the fixed λ from the MLP development run is too strong for GIN. We re-applied the same selection rule on BACE-GIN, sweeping λ ∈ {1, 3, 10, 30, 100, 300}: λ

id-acc

id-churn (%)

sym-KL

1 3 10 30 100 300

0.781 [0.775, 0.786] 0.769 [0.760, 0.777] 0.747 [0.736, 0.755] 0.706 [0.699, 0.713] 0.613 [0.602, 0.624] 0.582 [0.568, 0.599]

13.3 [12.6, 14.0] 10.6 [10.1, 11.1] 11.1 [10.5, 11.7] 16.1 [14.9, 17.3] 16.4 [14.0, 18.7] 7.1 [5.1, 9.1]

0.144 [0.133, 0.157] 0.091 [0.087, 0.096] 0.058 [0.054, 0.063] 0.028 [0.025, 0.032] 0.007 [0.006, 0.008] 0.004 [0.003, 0.005]

ERM-GIN id-acc is 0.742, so the 0.02-tolerance rule admits any λ keeping id-acc ≥ 0.722. The largest λ satisfying this is λ=10 (bold). At λ=10 twin-bootstrap cuts the class-flip rate −51.8% 25

vs. ERM-GIN (paired CI [−13.2, −10.6] pp), close to the bagging-K=5 reduction of −54%, and improves id-accuracy by +0.5 pp. The selection rule transfers across architectures unchanged; the numerical value does not. For practitioners changing the architecture, the practical recipe is a single λ sweep on the development dataset before transferring to held-out data.

ChemBERTa: rule-selected λ=10 preserves accuracy and cuts churn 9–76%

We fine-tuned ChemBERTa-77M-MTR [Chithrananda et al., 2020] on six main-comparison datasets (BACE, BBBP, Pgp, BBB-Martins, AMES, DILI), running ERM and a λ sweep for twin-bootstrap on BACE-ChemBERTa, then twin-bootstrap at the rule-selected λ on the five held-out datasets. All runs use five train-seeds, giving 52 = 10 seed-pairs per dataset. The λ chosen on the BACE MLP does not transfer. At λ=300 (the value selected on BACE-MLP) twin-bootstrap drops ChemBERTa accuracy by 11–18pp on every dataset and fails to reduce the classflip rate on 5/6 datasets: the across-seed-mean paired ∆ vs. ERM is positive on BACE, Pgp, and DILI (strict increase, all three canonical seeds agree), and the one apparent decrease (BBBP) is a model collapse to the majority-class predictor, not a real cut. AMES and BBB-Martins are within paired CI of zero on average across seeds. The underlying cause matches the GIN architecture cross-check (Appendix R): the λ tuned on a from-scratch MLP is too strong for pretrained representations. Rule-selected λ on BACE-ChemBERTa is λ=10. On BACE-ChemBERTa we sweep λ ∈ {1, 3, 10, 30, 100, 300} and apply the same 0.02 id-acc tolerance rule used on BACE-MLP. The first three values keep id-acc within tolerance; λ ≥ 30 fall below. The rule picks the largest admissible, λ=10. Table 16: The rule transfers; the value λ takes does not. At the λ=300 chosen on the BACE MLP, twin-bootstrap over-regularises ChemBERTa (accuracy drops 9–17 pp; churn rises on 5/6 datasets, BBBP collapses to majority). Re-applying the same 0.02-tolerance rule on BACE-ChemBERTa picks λ=10, at which twin-bootstrap preserves accuracy (within 2 pp of ERM) and cuts churn 15–82% on every dataset. Paired ∆ churn columns report mean [ 95% CI ] in percentage points over 52 =10 seed pairs. ERM Dataset BACE BBBP Pgp BBB-Martins AMES DILI

Twin-bootstrap λ=300

Twin-bootstrap λ=10 (rule)

churn (%)

acc

∆ churn (pp)

acc

∆ churn (pp)

15.9 4.7 8.1 3.3 8.8 27.8

0.60 0.78 0.63 0.77 0.58 0.50

+10.2 [+8.4, +12.1] -4.5 [-4.8, -4.1]∗ +10.2 [+8.9, +11.4] -0.0 [-0.7, +0.6] +3.4 [+1.9, +5.1] +16.2 [+12.0, +20.4]

0.70 0.88 0.78 0.87 0.75 0.67

-1.4 [-2.5, -0.2] -1.3 [-1.7, -0.8] -1.6 [-2.2, -1.0] -2.5 [-2.8, -2.2] -2.1 [-2.4, -1.7] -5.4 [-7.9, -3.0]

∗ BBBP at λ=300 collapses to the majority-class predictor (acc 0.78= majority), so the churn drop is meaningless.

At rule-selected λ=10, twin-bootstrap works on every ChemBERTa dataset. The failure documented at λ=300 on ChemBERTa is a λ-transfer failure, not a method failure: the rule transfers across architectures (MLP → GIN, MLP → ChemBERTa) and modalities (Morgan fingerprints → SMILES tokens) unchanged; the value λ takes does not. Appendix T reports the matching result on the vision pretrained backbone.

Waterbirds (ResNet-50): rule-selected λ=10 cuts churn 55%

The same protocol is applied to Waterbirds (ImageNet-pretrained ResNet-50, N =4795, single-task binary classification) with five train-seeds. We sweep λ ∈ {1, 3, 10, 30, 60, 100, 300} and apply the same 0.02-tolerance selection rule. The rule transfers across all three architectures and modalities we tested (MLP, GIN, ResNet50/ChemBERTa); the value λ takes does not. Both pretrained backbones (ChemBERTa on SMILES, 26

Table 17: Waterbirds (ImageNet-ResNet50): the 0.02-tolerance rule picks λ=10, recovering the closed-loop result on a vision pretrained backbone. ERM id-acc 0.875 (rule threshold ≥ 0.855). Twin-bootstrap at the λ=300 chosen on BACE collapses accuracy by 27 pp; at the rule-selected λ=10, twin-bootstrap preserves accuracy and cuts argmax churn 52%. All cells report mean [ 95% CI ] over 5 five train-seeds ( 2 =10 pairs for paired quantities). id-acc

id-churn (%)

∆ id-churn vs ERM (pp)

0.876 [0.859, 0.890] 0.872 [0.863, 0.880] 0.618 [0.612, 0.626]

11.1 [10.5, 11.7] 5.0 [4.8, 5.2] 3.5 [3.0, 4.0]

— -6.0 [-6.7, -5.4] (-55%) -7.6 [-8.5, -6.7] (-69%)⋆

Method ERM Twin-bootstrap λ=10 Twin-bootstrap λ=300

⋆ λ=300 fails the rule’s 0.02 accuracy tolerance (id-acc 0.618 is 27 pp below ERM’s 0.876); shown as a Pareto-sweep diagnostic, not a viable operating po

ResNet-50 on Waterbirds) and the from-scratch GIN pick λ=10; only the from-scratch fingerprintMLP picks λ=300.

Regression: methods rank the same way as on classification

The main experiments are binary classification. To test whether the methods extend to continuous targets we run the same protocol on three small-N MoleculeNet regression benchmarks: ESOL (N =1128, log aqueous solubility), FreeSolv (N =642, hydration free energy in kcal/mol), and Lipophilicity (N =4200, log D7.4 ). Cross-sample churn for regression is the per-example absolute prediction difference between two retrainings on independent bootstraps: ρreg (A, x) = ESA ,SB |fSA (x) − fSB (x)|. The twin-bootstrap consistency loss is MSE between the two networks’ predictions (instead of sym-KL on softmax outputs); the rule-selected λ is the largest value in {1, 3} keeping id-MAE within 0.04 of ERM id-MAE, which is λ=3 on all three datasets (both candidates improve id-MAE, so the tolerance constraint is not active). All other protocol choices (canonical-data seed 99, 10 train-seeds, paired-bootstrap CIs over 45 seed-pairs) are unchanged. Table 18: At matched 2×-ERM compute, twin-bootstrap beats bagging-K=2 on every regression dataset; the methods rank identically to the classification main-table result. Per-method id-MAE and cross-sample churn (mean |fA − fB | between bootstrap retrainings) on three MoleculeNet regression benchmarks. All reported quantities are mean [ 95% CI ] over the 45 pairs of 10 retrainings (or over the 10 retrainings themselves for id-MAE). Paired ∆ churn vs. ERM in the bottom rows; bold marks the best matched-compute method per dataset. Bagging-K=5 (5×-ERM compute) included as a stronger no-cost reference. ESOL (N =1128) Method ERM Bagging-K=2 Bagging-K=5 Twin-bootstrap λ=3

FreeSolv (N =642)

Lipo (N =4200)

id-MAE

churn

id-MAE

churn

id-MAE

churn

0.94 0.89 0.86 0.86

0.52 0.35 0.22 0.30

1.65 1.60 1.54 1.56

0.69 0.49 0.29 0.42

0.66 0.63 0.61 0.62

0.38 0.26 0.17 0.23

Paired ∆ churn vs. ERM (45 pairs of 10 retrainings, mean [ 95% CI ]): Bagging-K=2 -0.17 [-0.19, -0.15] (-32%) -0.20 [-0.23, -0.16] (-29%) Bagging-K=5 -0.30 [-0.32, -0.28] (-57%) -0.40 [-0.43, -0.36] (-58%) Twin-bootstrap λ=3 -0.22 [-0.23, -0.21] (-42%) -0.27 [-0.30, -0.23] (-39%)

-0.11 [-0.12, -0.11] (-30%) -0.21 [-0.21, -0.20] (-55%) -0.15 [-0.16, -0.15] (-40%)

At matched compute, twin-bootstrap cuts regression churn by 39–42% while bagging-K=2 cuts 29–32%; twin-bootstrap wins on every dataset by ∼10pp of relative reduction. This mirrors the classification main result (Section 6). Bagging-K=5 at 5× ERM compute is still strongest at 55–58%, consistent with the classification pattern, where 5×-compute bagging beats 2×-compute twin-bootstrap on 3/8 held-out datasets (DILI, Pgp, TADF). All four methods improve id-MAE by 3–9% over ERM; no accuracy regression on this task.

Compute footprint

Twin-bootstrap at K=2 trains two networks jointly: each batch is processed by both networks (four forwards), and a single backward pass on the joint loss updates both networks’ parameters. 27

Table 19: Twin-bootstrap costs ∼2× ERM wall-clock — matched compute against baggingK=2. Per-step compute footprint of each method, expressed in units of an ERM training step. Forward / backward counts are per training step on a single device; “Test models” is the number of trained networks queried per test prediction. Wall-clock is given for sequential training; ensemble methods can be parallelised over K devices to recover ERM wall-clock. Method ERM Deep ensemble K=5 Bagging K=2 Bagging K=5 Twin-bootstrap (K=2, joint)

Train fwd/step

Train bwd/step

Test models

Wall-clock vs. ERM

1 5 2 5 4

1 5 2 5 1 (joint)

1 5 2 5 2

1× 5× (sequential) 2× (sequential) 5× (sequential) ∼2×

In wall-clock terms this is roughly 2× ERM, matched to bagging-K=2. The “matched compute” comparison in Section 6 (Table 2) is at this 2×-ERM cost. The 5× Bagging-K=5 results in the same table are not matched compute; we include them as a stronger no-cost baseline.

Record · ID 180643 · SHA-256 aa2ee77c290dbe10
Retrieved via Conceptio — every document is proof-bundled with source, license, and retrieval metadata.