How can I guarantee a significant result?

Publish or perish, they say. Why not both?
model selection
glass-box modeling
post-selection inference
analysis
R
satire
Author

Perry Hackman, PhD

Published

December 22, 2025

Dear student,

Before you embark on your “career” as a statistician, you must purge yourself of a childish misconception: that our job is to seek truth. Truth is stubborn, unpredictable, and worst of all, often unpublishable. Scientists crave confidence, the journals crave significance, and we, if we are clever, can provide both without the nuisance of real rigor.

In this post series, I will instruct in the Statistical Dark Arts. Today, I’ll describe how to always ensure a publishable result with model selection and unadjusted post selection inferences (UPSIs).

Why UPSIs?

You will frequently come across scientists with ambitious ideas, saying things like:

I’ve conceived a brilliant new way to treat chronic pain effectively. I don’t want to waste time, efforts, and money on a result that ends up being insignificant… How can I guarantee a significant result??

No statistician wants to contribute to a null finding. Interpreting those is too hard! I’ll let you in on a statistical secret: there IS a way to guarantee statistical significance – and it’s actually extremely common in science.

It’s called unadjusted post-selection inference (UPSI).

With good model selection tools and UPSIs in your toolbelt, you can turn any negative study into a positive one – guaranteed.

Example: a study for treatment of chronic pain

Say the researcher who approached you is planning a cross-over study where patients were given one of two treatments for chronic pain. Patients will record their overall pain each day for 1 week while undergoing treatment A, and 1 week while undergoing treatment B.

Investigators are hoping to determine the difference between treatments in the average pain score. Our outcome \(Y_i\) is thus a continuous measure for pain reduction: for each subject, this value reflects the difference in that subject’s weekly average pain score between the two candidate treatments.

In addition to determining whether the treatment works on average in their population, investigators wish to determine specific subgroups that might see the most benefit. They are primarily interested in subgroups by age and sex, as well as patient self-reported data such as alcohol consumption and physical activity.

Note

Make sure they have an exhaustive list of candidate effect modifiers; I’d even suggest a few new ones like left/right handed-ness and coffee consumption.

Here’s a way to list all the possible subgroups of this set of variables:

subgroups <- expand.grid(
  age = c("18-35", "36-50", "51-65", "66-80", "80+"),
  sex = c("Male", "Female"),
  hand = c("Left", "Right"), 
  coffee = c("0", "1", "2+"),
  alcohol = c("0", "1", "2", "3+"),
  physical_activity = c("0", "1", "2", "3+")
)
nrow(subgroups)
[1] 960

This is a LOT of subgroups! And it’s way too many subgroups to possibly sift through, at least without tipping off the statistical reviewers about multiplicity.

CautionMultiplicity

Gone, alas, are the carefree days when we could run tons of tests, report the precious few that were significant, and quietly ignore the rest. Some meddlesome truth-seekers eventually caught on and began scolding us for our ingenuity, slapping the sinister label multiplicity onto our beloved golden-egg-laying goose. Now the reputable journals know to look for it. Our approach must become more subtle.

Here’s the key: turn to model selection. Write in the analysis plan:

We will use forward selection to determine which variables or interactions improve our model’s predictions, adding each candidate predictor in a stepwise fashion until the AIC indicates the model’s predictions can no longer be improved.

Sounds rigorous and honest, right?

This simple trick can guarantee significant results, and help you ensure your study will find a statistically significant result by the end, even when none exists!

How about that?! ZERO RISK! (well, to us anyway).

As a bonus, the more interactions we include that are truly null, the more unnecessarily opaque our model becomes – talk about a win-win.

Virtual Trial 1

Let me illustrate via a simulation. Let’s virtually collect \(n = 100\) patients in our trial.

library(tidyverse)

set.seed(21)
n <- 100        # Sample size  

# Simulate recruitment (y: outcome)
y <- rnorm(n)

# Simulate recruitment (assume each new person has random subgroup)
x_idx <- sample(1:nrow(subgroups), n, replace = TRUE) 
X <- subgroups[x_idx,]
simdata <- tibble(y = y, X)
simdata
# A tibble: 100 × 7
          y age   sex    hand  coffee alcohol physical_activity
      <dbl> <fct> <fct>  <fct> <fct>  <fct>   <fct>            
 1  0.793   36-50 Female Left  0      2       1                
 2  0.522   80+   Female Left  2+     1       3+               
 3  1.75    36-50 Female Left  2+     1       0                
 4 -1.27    66-80 Female Right 0      3+      0                
 5  2.20    51-65 Male   Right 2+     0       2                
 6  0.433   36-50 Male   Left  0      2       1                
 7 -1.57    66-80 Male   Right 0      2       0                
 8 -0.935   80+   Female Left  1      3+      3+               
 9  0.0635  36-50 Female Left  2+     0       2                
10 -0.00239 66-80 Male   Left  0      2       0                
# ℹ 90 more rows

We need to create interactions to select from:

# Create model matrix w/all interactions
X2 <- model.matrix(y ~ . * ., data = simdata)[,-1]

# Sanity check: make sure these columns represent main effects + interactions
sample(colnames(X2), 10) 
 [1] "alcohol2"                     "sexFemale:coffee1"           
 [3] "alcohol1:physical_activity3+" "age66-80:handRight"          
 [5] "alcohol1:physical_activity1"  "age51-65:sexFemale"          
 [7] "handRight:coffee1"            "age66-80:physical_activity2" 
 [9] "age51-65:coffee2+"            "age36-50:coffee1"            
# combine into data set for model fitting
simdata2 <- data.frame(y=y, X2)

OK! We have 93 candidate predictors for model selection. Don’t worry, unless our paper’s reviewers are extremely thorough, we don’t have to report this number. Our final model will have many fewer.

Model selection is a simple task. The following code uses forward step-wise selection with AIC to build an optimally-predicting model, per our analysis plan.

library(selectInferToolkit)
fit_aic <- select_stepwise_ic(y ~ ., data = simdata2, direction = "forward")
Caution

While other packages provide p-values after selection, most do so without being transparent, making it difficult to ensure they are truly UPSIs. If you aren’t careful with these less transparent tools, you might report p-values that are somehow diabolically adjusted towards insignificance. In contrast, the selectInferToolkit package, available on GitHub, helps you explicitly use UPSIs for inference:

infer_upsi(fit_aic, data = simdata2) %>% 
  tidy() %>% 
  filter(coef != 0)
term coef ci_low ci_high p_value
(Intercept) 0.07 −0.07 0.21 0.315
age36.50 0.17 −0.06 0.41 0.156
age66.80 −0.35 −0.55 −0.15 <0.001
handRight 0.23 −0.02 0.47 0.076
age80..sexFemale −0.29 −0.51 −0.07 0.013
age36.50.handRight −0.21 −0.43 0.01 0.071
age80..handRight −0.22 −0.44 0.00 0.053
age66.80.coffee1 0.45 0.19 0.71 <0.001
age36.50.coffee2. 0.19 −0.01 0.40 0.070
age51.65.coffee2. 0.33 0.16 0.49 <0.001
age36.50.alcohol1 −0.13 −0.31 0.06 0.178
age66.80.alcohol1 −0.26 −0.50 −0.01 0.042
age66.80.alcohol2 −0.19 −0.40 0.02 0.081
age51.65.physical_activity1 −0.33 −0.50 −0.17 <0.001
age51.65.physical_activity2 −0.18 −0.35 −0.01 0.039
sexFemale.alcohol1 0.33 0.15 0.50 <0.001
sexFemale.alcohol2 0.19 0.02 0.37 0.036
sexFemale.alcohol3. −0.17 −0.37 0.03 0.091
sexFemale.physical_activity2 −0.27 −0.44 −0.10 0.003
handRight.alcohol3. −0.16 −0.36 0.04 0.123
handRight.physical_activity3. −0.29 −0.48 −0.10 0.004
coffee2..alcohol3. 0.12 −0.05 0.29 0.172
coffee1.physical_activity1 0.17 0.01 0.33 0.040
alcohol3..physical_activity3. 0.25 0.08 0.42 0.006

And there you have it!

For many subgroups, treatment A worked great (increased pain scores, significant positive effects). Treatment A didn’t work so well for others (with negative significant effects), for whom it actually appears to decrease pain scores. These results indicate we could stand to maximize pain by giving certain patients A, and other patients B. I wouldn’t fret too much about how we included interactions without their constituent main effects; leave that subtlety to the media to figure out.

Caution

Sometimes, researchers want to actually decrease pain scores; you might want to check on this.

What’s with the high p-values?

I know what you’re thinking - won’t it be confusing to include all those results with high p-values?

An expert tip: if you want to make fewer discoveries and get even lower p-values, simply use BIC instead of AIC. BIC only lets the most significant results into the model, and you can say it’s asymptotically consistent, which sounds equally rigorous to AIC, which is asymptotically efficient.

In this example, here’s the model selected via stepwise BIC:

fit_bic <- select_stepwise_ic(y ~ ., data = simdata2, 
                              direction = "forward", 
                              penalty = "BIC")
infer_upsi(fit_bic, data = simdata2) %>% 
  tidy() %>% 
  filter(coef != 0)
term coef ci_low ci_high p_value
(Intercept) 0.07 −0.09 0.24 0.393
age36.50 0.27 0.10 0.45 0.003
age66.80 −0.38 −0.57 −0.18 <0.001
age66.80.coffee1 0.23 0.04 0.42 0.019
age51.65.coffee2. 0.29 0.11 0.46 0.002
age51.65.physical_activity1 −0.24 −0.42 −0.07 0.007
sexFemale.physical_activity2 −0.32 −0.49 −0.15 <0.001

Voila – the headlines practically write themselves!

How UPSIs work so well

Let me let you in on a little secret. If you were paying attention, you’d have seen the outcomes in the preceding example, \(Y_i\), were generated completely randomly. There was, by design, no true relationship to discover at all. Yet, I was virtually guaranteed to find a significant result. This is due to multiplicity’s clever younger brother, selective inference.

Selective inference is a difficult problem to solve; in fact, some call it impossible. Others hold out hope. At any rate, UPSIs have inertia and incentives on their side.

Later posts to this blog will go into specific solutions to the selective inference problem and how to implement them in R. Fair warning though, these solutions are often less promising and will certainly be less significant.

For now, let me finish convincing you that UPSIs are indeed effective at resolving the “publish or perish” dilemma; with USPIs, we can thrive at both! It should suffice to repeat this simulation again, as though it were performed in a parallel universe.

Virtual Trial 2

simdata2$y <- rnorm(n)
fit_bic <- select_stepwise_ic(y ~ ., data = simdata2, 
                              direction = "forward", 
                              penalty = "BIC")
infer_upsi(fit_bic, data = simdata2) %>% 
  tidy() %>% 
  filter(coef != 0)
term coef ci_low ci_high p_value
(Intercept) 0.03 −0.14 0.21 0.715
age66.80.sexFemale −0.25 −0.43 −0.07 0.008
age80..coffee2. 0.30 0.11 0.50 0.003
age80..physical_activity1 −0.26 −0.46 −0.07 0.009

Indeed, we still discover “significant” relationships. And, these are completely different treatment effect modifiers. In this parallel universe, our scientist friends publish completely different effects and will no doubt pour their hard-fought resources into what we well know are wild goose chases.

This is not a fluke.

In case you are still skeptical, let’s repeat this 50 times, as though we ran the study in 50 parallel universes.

sim_results <- list()
for(s in 1:50) { 
  simdata2$y <- rnorm(n)
  fit_bic <- select_stepwise_ic(y ~ ., data = simdata2, 
                                direction = "forward", penalty = "BIC")
  sim_results[[s]] <- infer_upsi(fit_bic, data = simdata2) %>% 
    tidy() %>% 
    filter(coef != 0) 
}

all_selections <- bind_rows(sim_results, .id = "simulation")
main_results <- all_selections %>%
  filter(term != "(Intercept)") %>%
  group_by(simulation) %>%
  summarize(
    any_p_lt_05 = any(p_value < 0.05), 
    n_selections = n()
  )

Using BIC, we found at least one significant effect in 46 of 50 trials, and on average these trials found 2.6 significant effects. We achieved falsely significant results 92% of the time.

If you’re thinking “hmm, my study still might have a chance of failing”, well good news. AIC finds a significant result nearly 100% of the time. We also only considered pairwise interactions between 6 candidate predictors. To minimize the risk of a negative study, one could consider higher level interactions or a greater number of subgroups, and a positive result quickly becomes a sure thing.

As I’ve said, UPSIs have incentives and inertia on their side. The more broad we make this practice, the better we can do at ensuring that the publications keep flowing, regardless of the truth.

Despicably yours,

P. Hackman

Conclusion

Key Takeaways

  • Unadjusted post-selection inference nearly guarantees false positives
  • Searching for subgroup effects is a model selection problem
  • Model selection algorithms don’t make multiplicity issues go away; they make them more subtle and harder to adjust for.

Future Threads

  • What are some alternatives to UPSI-based inference, and how are they implemented in R?
  • What is the selectInferToolkit package?

Appendix

R Session Info 
sessioninfo::session_info()
─ Session info ───────────────────────────────────────────────────────────────
 setting  value
 version  R version 4.5.2 (2025-10-31)
 os       macOS Sequoia 15.7.2
 system   aarch64, darwin20
 ui       X11
 language (EN)
 collate  en_US.UTF-8
 ctype    en_US.UTF-8
 tz       America/Chicago
 date     2025-12-22
 pandoc   3.6.3 @ /Applications/RStudio.app/Contents/Resources/app/quarto/bin/tools/aarch64/ (via rmarkdown)
 quarto   1.7.32 @ /Applications/RStudio.app/Contents/Resources/app/quarto/bin/quarto

─ Packages ───────────────────────────────────────────────────────────────────
 package            * version     date (UTC) lib source
 adaptMCMC            1.5         2024-01-29 [1] CRAN (R 4.5.0)
 backports            1.5.0       2024-05-23 [1] CRAN (R 4.5.0)
 broom              * 1.0.11      2025-12-04 [1] CRAN (R 4.5.2)
 class                7.3-23      2025-01-01 [2] CRAN (R 4.5.2)
 cli                  3.6.5       2025-04-23 [1] CRAN (R 4.5.0)
 coda                 0.19-4.1    2024-01-31 [1] CRAN (R 4.5.0)
 codetools            0.2-20      2024-03-31 [2] CRAN (R 4.5.2)
 data.table           1.17.8      2025-07-10 [1] CRAN (R 4.5.0)
 digest               0.6.39      2025-11-19 [1] CRAN (R 4.5.2)
 dplyr              * 1.1.4       2023-11-17 [1] CRAN (R 4.5.0)
 evaluate             1.0.5       2025-08-27 [1] CRAN (R 4.5.0)
 farver               2.1.2       2024-05-13 [1] CRAN (R 4.5.0)
 fastmap              1.2.0       2024-05-15 [1] CRAN (R 4.5.0)
 forcats            * 1.0.1       2025-09-25 [1] CRAN (R 4.5.0)
 foreach              1.5.2       2022-02-02 [1] CRAN (R 4.5.0)
 fs                   1.6.6       2025-04-12 [1] CRAN (R 4.5.0)
 future               1.68.0      2025-11-17 [1] CRAN (R 4.5.2)
 future.apply         1.20.0      2025-06-06 [1] CRAN (R 4.5.0)
 generics             0.1.4       2025-05-09 [1] CRAN (R 4.5.0)
 ggplot2            * 4.0.1       2025-11-14 [1] CRAN (R 4.5.2)
 glmnet               4.1-10      2025-07-17 [1] CRAN (R 4.5.0)
 globals              0.18.0      2025-05-08 [1] CRAN (R 4.5.0)
 glue                 1.8.0       2024-09-30 [1] CRAN (R 4.5.0)
 gower                1.0.2       2024-12-17 [1] CRAN (R 4.5.0)
 gt                 * 1.1.0       2025-09-23 [1] CRAN (R 4.5.0)
 gtable               0.3.6       2024-10-25 [1] CRAN (R 4.5.0)
 hardhat              1.4.2       2025-08-20 [1] CRAN (R 4.5.0)
 hms                  1.1.4       2025-10-17 [1] CRAN (R 4.5.0)
 htmltools            0.5.9       2025-12-04 [1] CRAN (R 4.5.2)
 htmlwidgets          1.6.4       2023-12-06 [1] CRAN (R 4.5.0)
 intervals            0.15.5      2024-08-23 [1] CRAN (R 4.5.0)
 ipred                0.9-15      2024-07-18 [1] CRAN (R 4.5.0)
 iterators            1.0.14      2022-02-05 [1] CRAN (R 4.5.0)
 jsonlite             2.0.0       2025-03-27 [1] CRAN (R 4.5.0)
 knitr                1.50        2025-03-16 [1] CRAN (R 4.5.0)
 lattice              0.22-7      2025-04-02 [2] CRAN (R 4.5.2)
 lava                 1.8.2       2025-10-30 [1] CRAN (R 4.5.0)
 lifecycle            1.0.4       2023-11-07 [1] CRAN (R 4.5.0)
 listenv              0.10.0      2025-11-02 [1] CRAN (R 4.5.0)
 lubridate          * 1.9.4       2024-12-08 [1] CRAN (R 4.5.0)
 magrittr             2.0.4       2025-09-12 [1] CRAN (R 4.5.0)
 MASS                 7.3-65      2025-02-28 [2] CRAN (R 4.5.2)
 Matrix               1.7-4       2025-08-28 [2] CRAN (R 4.5.2)
 ncvreg               3.16.0      2025-10-09 [1] Github (pbreheny/ncvreg@5fecc8c)
 nnet                 7.3-20      2025-01-01 [1] CRAN (R 4.5.0)
 parallelly           1.45.1      2025-07-24 [1] CRAN (R 4.5.0)
 pbapply              1.7-4       2025-07-20 [1] CRAN (R 4.5.0)
 pillar               1.11.1      2025-09-17 [1] CRAN (R 4.5.0)
 pkgconfig            2.0.3       2019-09-22 [1] CRAN (R 4.5.0)
 prodlim              2025.04.28  2025-04-28 [1] CRAN (R 4.5.0)
 purrr              * 1.2.0       2025-11-04 [1] CRAN (R 4.5.0)
 R6                   2.6.1       2025-02-15 [1] CRAN (R 4.5.0)
 RColorBrewer         1.1-3       2022-04-03 [1] CRAN (R 4.5.0)
 Rcpp                 1.1.0       2025-07-02 [1] CRAN (R 4.5.0)
 readr              * 2.1.6       2025-11-14 [1] CRAN (R 4.5.2)
 recipes              1.3.1       2025-05-21 [1] CRAN (R 4.5.0)
 rlang                1.1.6       2025-04-11 [1] CRAN (R 4.5.0)
 rmarkdown            2.30        2025-09-28 [1] CRAN (R 4.5.0)
 rpart                4.1.24      2025-01-07 [2] CRAN (R 4.5.2)
 rstudioapi           0.17.1      2024-10-22 [1] CRAN (R 4.5.0)
 S7                   0.2.1       2025-11-14 [1] CRAN (R 4.5.2)
 sass                 0.4.10      2025-04-11 [1] CRAN (R 4.5.0)
 scales               1.4.0       2025-04-24 [1] CRAN (R 4.5.0)
 selectInferToolkit * 0.4.2       2025-12-15 [1] local
 selectiveInference   1.2.5       2019-09-07 [1] CRAN (R 4.5.0)
 sessioninfo          1.2.3       2025-02-05 [1] CRAN (R 4.5.0)
 shape                1.4.6.1     2024-02-23 [1] CRAN (R 4.5.0)
 sparsevctrs          0.3.4       2025-05-25 [1] CRAN (R 4.5.0)
 stringi              1.8.7       2025-03-27 [1] CRAN (R 4.5.0)
 stringr            * 1.6.0       2025-11-04 [1] CRAN (R 4.5.0)
 survival             3.8-3       2024-12-17 [2] CRAN (R 4.5.2)
 tibble             * 3.3.0       2025-06-08 [1] CRAN (R 4.5.0)
 tidyr              * 1.3.1       2024-01-24 [1] CRAN (R 4.5.0)
 tidyselect           1.2.1       2024-03-11 [1] CRAN (R 4.5.0)
 tidyverse          * 2.0.0       2023-02-22 [1] CRAN (R 4.5.0)
 timechange           0.3.0       2024-01-18 [1] CRAN (R 4.5.0)
 timeDate             4051.111    2025-10-17 [1] CRAN (R 4.5.0)
 tzdb                 0.5.0       2025-03-15 [1] CRAN (R 4.5.0)
 utf8                 1.2.6       2025-06-08 [1] CRAN (R 4.5.0)
 vctrs                0.6.5       2023-12-01 [1] CRAN (R 4.5.0)
 withr                3.0.2       2024-10-28 [1] CRAN (R 4.5.0)
 xfun                 0.54        2025-10-30 [1] CRAN (R 4.5.0)
 xml2                 1.5.1       2025-12-01 [1] CRAN (R 4.5.2)
 yaml                 2.3.11.9000 2025-12-10 [1] Github (r-lib/r-yaml@6dc4582)

 [1] /Users/rpterson/Library/R/arm64/4.5/library
 [2] /Library/Frameworks/R.framework/Versions/4.5-arm64/Resources/library
 * ── Packages attached to the search path.

──────────────────────────────────────────────────────────────────────────────

Thanks for reading

© 2025 Glassbox Modeling Working Group

Creative Commons License
Content on this site is licensed under the Creative Commons Attribution 4.0 International License (CC BY 4.0) .