stat.APApr 13, 2026
SaveObserving the unobserved confounding through its effects: toward randomized trial-like estimates from real-world survival data
Organizations: Operations Research Center, Massachusetts Institute of Technology, Cambridge, MA, USA · Sloan School of Management, Massachusetts Institute of Technology, Cambridge, MA, USA · Department of Surgery, Memorial Sloan Kettering Cancer Center, New York, NY, USA · Charité-Universitätsmedizin Berlin, Corporate Member of Freie Universität Berlin and Humboldt-Universität zu Berlin, Berlin, Germany
Abstract
Background: Randomized controlled trials (RCTs) are costly, time-consuming, and often infeasible, while treatment-effect estimation from observational data is limited by unobserved confounding. Methods: We developed a three-step framework to address unobserved confounding in observational survival data. First, we infer a latent prognostic factor (U) from restricted mean survival time (RMST) discrepancies between patients with similar observed factors, the same treatment, and divergent outcomes, leveraging the idea that the aggregate effect of unmeasured factors can be inferred even if individual factors cannot. Second, we balance U with observed baseline covariates using prognostic matching, entropy balancing, or inverse probability of treatment weighting. Third, we apply multivariable survival analysis to estimate hazard ratios (HRs). We evaluated the framework in three observational cohorts with RCT benchmarks, two RCT cohorts, and six multicenter observational cohorts. Results: In three observational cohorts (nine comparisons), balancing U improved agreement with trial HRs in all cases; in the strongest settings, it reduced absolute log-HR error by approximately ten-fold versus using observed covariates alone (mean reduction 0.344; p=0.001). In two RCT cohorts, U was balanced across arms (most SMDs <0.1) and adjustment had minimal impact on log-HRs (mean absolute change 0.08). Across six multicenter cohorts, balancing U within centers reduced cross-center dispersion in chemotherapy log-HR estimates (mean reduction 0.147; p=0.016); when populations were directly balanced across centers to account for case-mix differences, cross-center survival differences were narrowed in 75%-100% of comparisons. Conclusions: Inferring and balancing a latent prognostic signal may reduce unobserved confounding and improve treatment-effect estimation from real-world data.
Figures & tables
Figure 1: Full pipeline of the proposed framework.
| Dataset | N | Surgery + adjuvant therapy | Surgery alone | Events (Surgery + adjuvant therapy) | Events (Surgery alone) |
|---|---|---|---|---|---|
| GIST | 480 | 105 (21.9%) | 375 (78.1%) | 34 (32.4%) | 52 (13.9%) |
| RPS | 201 | 29 (14.4%) | 172 (85.6%) | 3 (10.3%) | 49 (28.5%) |
| CRLM | 778 | 140 (18.0%) | 638 (82.0%) | 83 (59.3%) | 445 (69.8%) |
Table 1: Cohort characteristics of the observational datasets used in the study.
Figure 2: Comparison of treatment-effect estimation accuracy across observational datasets. Bars represent the mean absolute log-HR error relative to the benchmark RCT hazard ratio, averaged across hyperparameter configurations. Results are shown for matching, entropy balancing, and inverse propensity weighting using either baseline covariates or the augmented representation . The dashed line denotes multivariable Cox adjustment using baseline covariates only. Error bars represent the standard error across configurations.
Figure 3: Comparison of treatment-effect estimates across observational datasets. Bars represent the mean hazard ratio (HR), averaged across hyperparameter configurations. Results are shown for matching, entropy balancing, and inverse propensity weighting using either baseline covariates or the augmented representation . The dotted line denotes multivariable Cox adjustment using baseline covariates only, and the dashed green line indicates the benchmark RCT hazard ratio. Error bars represent the standard error across configurations.
| Statistic | Value |
|---|---|
| Directional improvements ( ) | 9 / 9 |
| Sign test -value (one-sided) | 0.00195 |
| Wilcoxon signed-rank -value (one-sided) | 0.00195 |
| Median cell-level | 0.2782 |
| Mean cell-level | 0.3443 |
Table 2: Validation approach 1: directional improvement in closeness to the benchmark RCT hazard ratio after adding .
| Treatment | Events | ||||
|---|---|---|---|---|---|
| Dataset | N | Surgery + periop. therapy | Surgery alone | Surgery + periop. therapy | Surgery alone |
| BRT | 162 | 75 (46.3%) | 87 (53.7%) | 12 (16.0%) | 26 (29.9%) |
| STRASS | 265 | 132 (49.8%) | 133 (50.2%) | 59 (44.7%) | 61 (45.9%) |
Table 3: Cohort characteristics of the RCT datasets.
Figure 4: Comparison of treatment-effect estimation accuracy across RCT datasets. Bars represent the mean absolute log-HR error relative to the trial-reported hazard ratio, averaged across hyperparameter configurations. Results are shown for matching, entropy balancing, and inverse propensity weighting using either baseline covariates or the augmented representation . The dotted line denotes multivariable Cox adjustment using baseline covariates only. Error bars denote the standard error across configurations.
Figure 5: Comparison of treatment-effect estimates across RCT datasets. Bars represent the mean hazard ratio (HR), averaged across hyperparameter configurations, for matching, entropy balancing, and inverse propensity weighting using either baseline covariates or the augmented representation . Error bars denote the standard error across configurations. The dotted line indicates the multivariable estimate, and the dashed green line indicates the trial-reported ground-truth HR.
| Metric | Value |
|---|---|
| Equivalence wins | 4 / 6 |
| Mean | 0.081 |
Table 4: Validation approach 2: equivalence testing for RCT datasets. A “win” indicates that the augmentation with passes the equivalence test, meaning the 95% confidence interval of the shift lies entirely within the pre-specified equivalence margin ( ).
Figure 6: Standardized mean differences (SMDs) between treatment arms for observed covariates and estimated latent factors in the RCT datasets. The SMD values for are generally below or close to the commonly used balance threshold of 0.1 and are comparable to those of the observed covariates.
| Dataset | Event | N | Surgery + chemotherapy | Surgery alone | Events (Surgery + chemotherapy) | Events (Surgery alone) |
|---|---|---|---|---|---|---|
| Graz and Vienna Universities | Death | 128 | 88 (68.8%) | 40 (31.2%) | 47 (53.4%) | 12 (30.0%) |
| Graz and Vienna Universities | Recurrence | 128 | 88 (68.8%) | 40 (31.2%) | 68 (77.3%) | 27 (67.5%) |
| Johns Hopkins | Death | 461 | 296 (64.2%) | 165 (35.8%) | 160 (54.0%) | 116 (70.3%) |
| Johns Hopkins | Recurrence | 461 | 296 (64.2%) | 165 (35.8%) | 183 (61.8%) | 78 (47.3%) |
| Kumamoto University | Death | 90 | 46 (51.1%) | 44 (48.9%) | 23 (50.0%) | 26 (59.1%) |
| Kumamoto University | Recurrence | 90 | 46 (51.1%) | 44 (48.9%) | 39 (84.8%) | 38 (86.4%) |
Table 5: Cohort characteristics of CRLM multi-center datasets.
Figure 7: Between-center heterogeneity in treatment-effect estimates for the CRLM multi-center cohort. Bars represent the mean pairwise absolute difference in log-hazard ratios (log-HR) across centers, averaged over all center pairs within each method and configuration. Results are shown for matching, entropy balancing, and inverse propensity weighting using either baseline covariates or the augmented representation . Error bars denote the standard error across pairwise comparisons.
| Statistic | Value |
|---|---|
| Cells with improvement ( ) | 6 / 6 |
| Sign test -value (one-sided) | 0.0156 |
| Wilcoxon signed-rank -value (one-sided) | 0.0156 |
| Median cell-level | 0.1343 |
| Mean cell-level | 0.1467 |
Table 6: Validation approach 3, first component: reduction in mean pairwise absolute deviation after adjusting for for 6 CRLM centers.
| Subset | Event | Method | Pairs Evaluated | Pairs Widened by | Fixed (%) | Binomial |
|---|---|---|---|---|---|---|
| Full cohort | Mortality | Matching | 30 | 19 (63.3%) | 16/19 (84.2%) | 0.0022 |
| Full cohort | Mortality | IPTW | 30 | 15 (50.0%) | 13/15 (86.7%) | 0.0037 |
| Full cohort | Mortality | Entropy Balancing | 16 | 10 (62.5%) | 9/10 (90.0%) | 0.0107 |
| Full cohort | Recurrence | Matching | 30 | 10 (33.3%) | 8/10 (80.0%) | 0.0547 |
| Full cohort | Recurrence | IPTW | 30 | 18 (60.0%) | 15/18 (83.3%) | 0.0038 |
| Full cohort | Recurrence | Entropy Balancing | 16 | 9 (56.2%) | 7/9 (77.8%) | 0.0898 |
Table 7: Validation 3, second component: Cross-center validation of the latent prognostic factor . Among center pairs where adjustment for observed covariates increased disagreement in adjusted 5-year survival, we evaluated whether adding reduced the gap. Results are shown for the full cohort and for the surgery alone subset. The adjustment happens in the form of weights calculated to balance the pairs of centers with different balancing methods. These weights are then used in the Kaplan-Meier estimation.
Appendix figures & tables18 assets
Supplementary material from the paper’s appendix.
Appendix
Figure 1.1 : Sensitivity of treatment-effect estimation to the number of bins used in prognostic matching for the observational and RCT cohorts. The plotted values represent the mean absolute log-HR error relative to the benchmark hazard ratio, averaged across hyperparameter configurations. Results are shown with and without the latent factor . Error bars denote the standard error across configurations.
Figure 1.2 : Sensitivity of treatment-effect estimation to the number of bins used in prognostic matching for the 6 CRLM centers. The plotted values represent the mean pairwise absolute difference of log-HR between centers, averaged across hyperparameter configurations. Results are shown with and without the latent factor . Error bars denote the standard error across configurations.
Figure 1.3 : Sensitivity of inverse propensity weighting (IPTW) to the clipping threshold applied to the estimated propensity scores for the observational and RCT cohorts. The plotted values represent the mean absolute log-HR error relative to the benchmark hazard ratio across hyperparameter configurations. Results are shown for models using baseline covariates and the augmented representation . Error bars denote the standard error across configurations.
Figure 1.4 : Sensitivity of inverse propensity weighting (IPTW) to the clipping threshold applied to the estimated propensity scores for the 6 CRLM centers. The plotted values represent the mean pairwise absolute difference of log-HR between centers, averaged across hyperparameter configurations. Results are shown for models using baseline covariates and the augmented representation . Error bars denote the standard error across configurations.
Figure 1.5 : Sensitivity of entropy balancing to the number of balanced moments for the observational and RCT datasets. Bars represent the mean absolute log-HR error relative to the benchmark hazard ratio, averaged across hyperparameter configurations. Results compare models using baseline covariates with those augmented by the latent factor . Error bars denote the standard error across configurations.
Figure 1.6 : Sensitivity of entropy balancing to the number of balanced moments for the 6 CRLM centers. Bars represent the mean pairwise absolute difference of log-HR between centers, averaged across hyperparameter configurations. Results compare models using baseline covariates with those augmented by the latent factor . Error bars denote the standard error across configurations.
Figure 1.7 : Sensitivity of the framework to the number of nearest neighbors ( ) used when constructing the latent prognostic factor for the observational and RCT datasets. Lines represent the mean absolute log-HR error relative to the benchmark hazard ratio across configurations for each balancing method. Error bars denote the standard error across configurations.
Figure 1.8 : Sensitivity of the framework to the number of nearest neighbors ( ) used when constructing the latent prognostic factor for the 6 CRLM centers. Lines represent the the mean pairwise absolute difference of log-HR between centers, averaged across hyperparameter configurations for each balancing method. Error bars denote the standard error across configurations.
Figure 1.9 : Impact of incorporating the prognostic score in the distance metric used to construct for the observational and RCT datasets. Bars represent the mean absolute log-HR error relative to the benchmark hazard ratio across hyperparameter configurations for each balancing method. Error bars denote the standard error across configurations.
Figure 1.10 : Impact of incorporating the prognostic score in the distance metric used to construct for the 6 CRLM centers. Bars represent the mean pairwise absolute difference of log-HR between centers, averaged across hyperparameter configurations for each balancing method. Error bars denote the standard error across configurations.
Figure 1.11 : Sensitivity of results to the choice of prognostic score model used to estimate the baseline risk component for the observational and RCT datasets. Bars represent the mean absolute log-HR error relative to the benchmark hazard ratio across hyperparameter configurations for each balancing method. Error bars denote the standard error across configurations.
Figure 1.12 : Sensitivity of results to the choice of prognostic score model used to estimate the baseline risk component for the 6 CRLM centers. Bars represent the mean pairwise absolute difference of log-HR between centers, averaged across hyperparameter configurations for each balancing method. Error bars denote the standard error across configurations.
Figure 2.13 : Permutation and ablation test for the observational datasets. The latent factor is randomly permuted across patients before balancing. Bars represent the mean absolute log-HR error relative to the benchmark hazard ratio. Performance degradation relative to the non-permuted case indicates that the predictive signal of contributes to improved treatment-effect estimation. Error bars denote the standard error across configurations.
Figure 2.14 : Permutation and ablation test for the RCT datasets. The latent factor is randomly permuted across patients before balancing. Because treatment assignment is randomized, permuting does not materially change treatment-effect estimates, providing a negative control for the method. Error bars denote the standard error across configurations.
Figure 2.15 : Permutation test for the CRLM multi-center analysis. The latent factor is randomly permuted across patients prior to balancing. Bars represent the mean absolute pairwise HR difference across centers. The increase in cross-center variability after permutation indicates that the estimated captures meaningful latent prognostic structure shared across cohorts. Error bars denote the standard error across center pairs.
| Dataset | Method | Mean | Std. Error |
|---|---|---|---|
| GIST | Matching | 0.334 | 0.012 |
| Entropy Balancing | 0.256 | 0.022 | |
| IPTW | 0.326 | 0.026 | |
| RPS | Matching | 0.084 | 0.008 |
| Entropy Balancing | 0.592 | 0.031 | |
| IPTW | 0.903 | 0.007 |
Table 2.1 : Comparison of balancing methods across datasets. represents directional improvement in closeness to the benchmark RCT hazard ratio after adding . Mean difference is reported on the log hazard ratio scale.
| Dataset | Method | Mean Diff. (log HR) | Std. Error | 95% CI |
|---|---|---|---|---|
| BRT | Matching | 0.075 | 0.0085 | [0.058, 0.092] |
| Entropy Balancing | 0.183 | 0.0221 | [0.139, 0.227] | |
| IPTW | 0.035 | 0.0163 | [0.002, 0.067] | |
| STRASS | Matching | -0.091 | 0.0058 | [-0.1025, -0.796] |
| Entropy Balancing | -0.038 | 0.0095 | [-0.057, -0.019] | |
| IPTW | -0.066 | 0.0079 | [-0.081, -0.050] |
Table 2.2 : Equivalence testing of treatment effect differences across datasets. Mean difference is calculated as .
| Event | Method | Mean | Std. Error |
|---|---|---|---|
| Mortality | Matching | 0.129 | 0.0087 |
| Entropy Balancing | 0.098 | 0.0327 | |
| IPTW | 0.224 | 0.0029 | |
| Recurrence | Matching | 0.139 | 0.0047 |
| Entropy Balancing | 0.169 | 0.0360 | |
| IPTW | 0.121 | 0.0035 |
Table 2.3 : Comparison of balancing methods across datasets. represents the reduction in mean pairwise absolute deviation after adjusting for for 6 CRLM centers. Mean difference is reported on the log hazard ratio scale.