Understanding POWERUP prediction uncertainty • powerup

In addition to point predictions, each predicted response is paired with a sample-specific prediction standard deviation, an approximate Gaussian prediction interval, and a “target-event” probability (i.e. the probability of falling on the response-of-interest side of a chosen cutoff).

Please note these quantities describe uncertainty in the predicted response, not uncertainty that the model is “correct.” For information about model performance interpretation, see Interpreting POWERUP results.

Generate prediction distributions

The bundled example can be used to fit a CTNNB1 model and predict two held-out samples.

expression_path <- system.file("extdata", "powerup_example_gene_expression.csv", package = "powerup")
dependency_path <- system.file("extdata", "powerup_example_gene_dependency.csv", package = "powerup")
user_path <- system.file("extdata", "powerup_example_user_matrix.csv", package = "powerup")

example_expression <- read.csv(expression_path, check.names = FALSE)
example_dependency <- read.csv(dependency_path, check.names = FALSE)
example_user <- read.csv(user_path, check.names = FALSE)

prepared <- prepare_powerup_data(
  gene_expression = example_expression,
  response = example_dependency,
  targets = "CTNNB1",
  user_matrix = example_user
)
#> [powerup] target resolution summary requestedUnique=1 resolvedUnique=1 metadataNameMatches=0 canonicalFallbackMatches=1 looseMetadataFallbackResolved=0 exactMetadataMultipleMatches=0 unresolvedSkipped=0
#> [powerup] top-variable feature selection source=reference_matrix reason=user_selected_reference_source reference_samples=300 n_ranked_reference_features=300 n_user_available_features=300

models <- fit_powerup_models(prepared, seed = 123L)
models <- add_powerup_predictions(models, prepared)

Use the long prediction summary to inspect the uncertainty outputs:

predictions <- summarize_predictions(models, format = "long")
predictions[, c(
  "sample",
  "perturbation",
  "pred_mean",
  "pred_sd",
  "pred_pi_lower_95",
  "pred_pi_upper_95",
  "prob_target_event"
)]
#> # A tibble: 2 × 7
#>   sample     perturbation pred_mean pred_sd pred_pi_lower_95 pred_pi_upper_95
#>   <chr>      <chr>            <dbl>   <dbl>            <dbl>            <dbl>
#> 1 ACH-000957 CTNNB1           0.983   0.134            0.720            1.25 
#> 2 ACH-002024 CTNNB1           0.193   0.167           -0.133            0.519
#> # ℹ 1 more variable: prob_target_event <dbl>

What the outputs mean

Output	Interpretation
`pred_mean`	POWERUP’s point estimate of the continuous response
`pred_sd`	Estimated scale of prediction error for that sample
`pred_pi_lower_95`, `pred_pi_upper_95`	Approximate Gaussian 95% prediction interval
`prob_below_cutoff`	Probability mass below the response cutoff
`prob_above_cutoff`	Probability mass above the response cutoff
`prob_target_event`	Probability mass on the response-of-interest side of the cutoff
`response_cutoff`	Cutoff used to define the response of interest
`decreasing`	Whether the response of interest is at or below the cutoff

pred_sd is not a standard error around an estimated mean. It describes how uncertain POWERUP expects the response prediction itself to be for that sample.

The 95% interval is therefore a prediction interval, not a confidence interval for the fitted mean.

Where the uncertainty comes from

POWERUP learns prediction uncertainty from out-of-fold errors rather than from the final model’s in-sample residuals. For an observed held-out response $y_i$ and its out-of-fold prediction $\hat{y}_i$ , the prediction error is

$e_i = y_i - \hat{y}_i.$

The uncertainty model is trained on the transformed squared error

$u_i = \log(e_i^2 + \epsilon),$

where $\epsilon = 10^{-8}$ prevents taking the logarithm of zero. Predicting log squared error keeps the modeled quantity unconstrained while ensuring that converting it back to the variance scale gives a positive value.

For each target, the procedure is:

repeated cross-validation generates out-of-fold response predictions;
POWERUP calculates the squared error for each out-of-fold prediction;
a second XGBoost model predicts log(squared error + epsilon) from the molecular features and the corresponding response prediction;
repeated out-of-fold rows from the same training sample are weighted so that each sample has the same total weight;
the uncertainty model is cross-fitted to obtain held-out estimates of its predicted error scale;
a global calibration factor rescales those error estimates;
the final uncertainty model is fit to all out-of-fold error rows and used to generate sample-specific prediction SDs.

The response prediction is included as an input to the uncertainty model because predictive error can depend both on molecular context and on where the sample lies on the modeled response scale. If the uncertainty model predicts $\hat{u}_i$ , the corresponding raw variance and SD are

$\tilde{\sigma}_i^2 = \exp(\hat{u}_i), \qquad \tilde{\sigma}_i = \sqrt{\exp(\hat{u}_i)}.$

POWERUP applies a small variance floor of $10^{-8}$ when converting these predictions back to the variance scale.

What calibration does

Before calibration, the uncertainty model produces a relative estimate of how large the prediction error should be for each sample. POWERUP first standardizes each cross-fitted absolute error by its raw predicted SD:

$r_i = \frac{|e_i|}{\tilde{\sigma}_i}.$

Repeated out-of-fold rows from the same training sample receive weights $w_i = 1/m_i$ , where $m_i$ is that sample’s number of out-of-fold rows. This gives each training sample the same total weight in calibration.

POWERUP then calculates the weighted 80th percentile $Q_{0.80}^{(w)}$ of these standardized errors. For a standard Gaussian, the symmetric central 80% interval extends to

$z_{0.80} = \Phi^{-1}(0.90) \approx 1.282.$

The global calibration factor is therefore

$c = \frac{Q_{0.80}^{(w)}(r_i)}{\Phi^{-1}(0.90)}.$

The final sample-specific uncertainty is obtained by multiplying the raw SD by this factor:

$\sigma_i = c\tilde{\sigma}_i, \qquad \sigma_i^2 = c^2\tilde{\sigma}_i^2.$

In practical terms, $c > 1$ widens the uncertainty distributions when the raw error model is too narrow, whereas $c < 1$ narrows them when it is too wide. Each target receives one global calibration factor across its cross-fitted training errors. Calibration aligns the overall 80% error scale in the training setting; it does not guarantee exact nominal coverage for every subgroup, individual sample, or out-of-distribution sample.

POWERUP currently exposes the resulting uncertainty estimates but does not report an interval-coverage diagnostic in summarize_models().

From prediction SD to a prediction interval

POWERUP treats the fitted predictive distribution for sample $i$ as Gaussian:

$Y_i \sim \mathcal{N}(\mu_i, \sigma_i^2),$

where $\mu_i$ is pred_mean and $\sigma_i$ is pred_sd. A two-sided prediction interval with nominal coverage $1-\alpha$ is

$\mu_i \pm \Phi^{-1}(1-\alpha/2)\sigma_i.$

For the reported 95% interval, $\Phi^{-1}(0.975) \approx 1.96$ , giving

$\text{pred\_mean} \pm 1.96 \times \text{pred\_sd}.$

A larger pred_sd therefore produces a wider interval. A smaller pred_sd produces a narrower interval.

Neither is automatically evidence of a better or worse response model. mean_pred_sd in summarize_models() is a descriptive summary of how wide the model’s predictive distributions tend to be, whereas r, R2, and RMSE evaluate held-out response prediction performance.

From the predictive distribution to an event probability

prob_target_event is calculated from the Gaussian predictive distribution and the response cutoff.

Let $t$ be the response cutoff and let $\Phi$ denote the standard normal cumulative distribution function. When decreasing = FALSE, the response of interest is $Y_i \ge t$ , so

$P(Y_i \ge t) = 1 - \Phi\left(\frac{t-\mu_i}{\sigma_i}\right).$

When decreasing = TRUE, the response of interest is $Y_i \le t$ , so

$P(Y_i \le t) = \Phi\left(\frac{t-\mu_i}{\sigma_i}\right).$

This means that the predicted mean and target-event probability answer different questions:

pred_mean: What response value does POWERUP expect?
prob_target_event: Given the predicted response and its uncertainty, how much of the predictive distribution lies on the response-of-interest side of the cutoff?

Two samples can have similar predicted means but different event probabilities when their prediction SDs differ, particularly when the means lie near the cutoff.

What prediction uncertainty does not include

The reported SD and interval describe uncertainty learned from prediction errors in the reference training setting. They should not be interpreted as capturing every possible source of uncertainty.

In particular, they do not automatically account for:

dataset shift between the reference models and a new tumor cohort;
all technical uncertainty in the user’s expression measurements;
unmeasured biological variables absent from the feature matrix;
uncertainty in the SHAP explanation itself;
uncertainty in the chosen biological response cutoff.

An unusually narrow interval should therefore not be treated as a guarantee that a prediction is correct, especially for samples that differ substantially from the training distribution.

A practical reading order

For an individual prediction:

inspect pred_mean for the expected response;
inspect pred_sd or the 95% interval for the estimated prediction uncertainty;
inspect prob_target_event when the biological question is whether the response crosses a defined cutoff;
interpret all three in the context of held-out model performance.

If experimental evidence becomes available, these predictive distributions can be used as priors for experimental observations and posterior updating. Posterior updating can optionally transform the predictive mean and SD to an outcome-appropriate latent Gaussian scale, such as logit for responses constrained to (0, 1) or log/log10 for strictly positive responses.

Session information

sessionInfo()
#> R version 4.4.2 (2024-10-31)
#> Platform: aarch64-apple-darwin20
#> Running under: macOS Sequoia 15.7.3
#> 
#> Matrix products: default
#> BLAS:   /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/lib/libRblas.0.dylib 
#> LAPACK: /Library/Frameworks/R.framework/Versions/4.4-arm64/Resources/lib/libRlapack.dylib;  LAPACK version 3.12.0
#> 
#> locale:
#> [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
#> 
#> time zone: America/New_York
#> tzcode source: internal
#> 
#> attached base packages:
#> [1] stats     graphics  grDevices utils     datasets  methods   base     
#> 
#> other attached packages:
#> [1] powerup_1.0.99
#> 
#> loaded via a namespace (and not attached):
#>  [1] Matrix_1.7-3      jsonlite_2.0.0    dplyr_1.1.4       compiler_4.4.2   
#>  [5] tidyselect_1.2.1  parallel_4.4.2    tidyr_1.3.1       jquerylib_0.1.4  
#>  [9] globals_0.16.3    systemfonts_1.2.2 textshaping_1.0.0 yaml_2.3.10      
#> [13] fastmap_1.2.0     lattice_0.22-6    R6_2.6.1          generics_0.1.3   
#> [17] knitr_1.50        htmlwidgets_1.6.4 future_1.34.0     tibble_3.3.0     
#> [21] desc_1.4.3        bslib_0.9.0       pillar_1.10.1     rlang_1.1.6      
#> [25] utf8_1.2.4        cachem_1.1.0      xfun_0.51         fs_1.6.5         
#> [29] sass_0.4.9        cli_3.6.5         pkgdown_2.2.0     withr_3.0.2      
#> [33] magrittr_2.0.3    grid_4.4.2        digest_0.6.37     xgboost_3.2.1.1  
#> [37] lifecycle_1.0.4   vctrs_0.6.5       data.table_1.17.0 evaluate_1.0.3   
#> [41] glue_1.8.0        rsample_1.2.1     listenv_0.9.1     furrr_0.3.1      
#> [45] codetools_0.2-20  ragg_1.5.1        parallelly_1.42.0 rmarkdown_2.29   
#> [49] purrr_1.1.0       tools_4.4.2       pkgconfig_2.0.3   htmltools_0.5.8.1