Planning a Correlation Study
In BUCSS: Bias and Uncertainty Corrected Sample Size

knitr::opts_chunk$set(collapse = TRUE, comment = "#>")
library(BUCSS)
# Helpers for numbers quoted in the text: one row of a result by its term, and
# a proportion at the printed tables' precision, without the leading zero.
row_value <- function(x, term) x$value[x$term == term]
apa <- function(x, digits = 3) {
  sub("^0\\.", ".", format(signif(x, digits), scientific = FALSE))
}

A Pearson correlation is one of the most commonly reported effects in the social sciences, and one of the most commonly replicated. Suppose a published study of 100 participants reported a correlation of r = .35, and you want to plan a study powered to find that relationship again. How large should it be? As with every planner in this package, the reported correlation should not be taken at face value: it was selected for publication in part because it was large, so planning from it without correction tends to yield a study that is too small. ss_buc_correlation() applies the bias and uncertainty correction of Anderson, Kelley, and Maxwell (2017), which reads the prior result as a draw from a distribution truncated at the publication threshold (Taylor & Muller, 1996), and returns the total sample size the replication needs.

first_plan <- ss_buc_correlation(r_observed = .35, N = 100, alpha_prior = .05,
                                 alpha_planned = .05, assurance = .80,
                                 desired_power = .80)
first_plan

The returned necessary_N is the total sample size and ncp_adjusted the bias and uncertainty adjusted noncentrality parameter. rho_planning is the planning effect size, the adjusted correlation the plan is built on given the prior result, alpha_prior, and assurance (not an estimate of the population correlation to report), sqrt(ncp_adjusted / (ncp_adjusted + N)) with the prior study's N: here r apa(row_value(first_plan, "rho_planning")), against the published .35. A study of r row_value(first_plan, "necessary_N") participants has exactly the returned actual_power when the population correlation is r apa(row_value(first_plan, "rho_planning")). For the meaning of assurance and alpha_prior, and how each moves the recommended sample size, see vignette("understanding-assurance", package = "BUCSS"). Briefly, raising assurance buys protection against an unlucky prior estimate at the cost of a larger study:

assurances <- c(.5, .7, .8, .9)
sizes <- sapply(assurances, function(a) {
  plan <- ss_buc_correlation(r_observed = .35, N = 100, assurance = a)
  plan$value[plan$term == "necessary_N"]
})
data.frame(assurance = assurances, necessary_N = sizes)

A Correlation Is Not a Regression Slope

The test of whether a correlation differs from zero is arithmetically the test of the slope in a simple regression, so it is tempting to plan a correlation study with the single-coefficient regression planner, treating the predictor as fixed. Doing so gives a slightly different, and slightly too small, answer:

t_equivalent <- .35 * sqrt(98 / (1 - .35^2))
by_correlation <- ss_buc_correlation(t_observed = t_equivalent, N = 100)
by_regression  <- ss_buc_reg_coef(t_observed = t_equivalent, N = 100, p = 1,
                                  random_predictors = FALSE)
c(correlation = by_correlation$value[by_correlation$term == "necessary_N"],
  regression  = by_regression$value[by_regression$term == "necessary_N"])

The gap is small here, three participants, but it is systematic (only near the assurance ceiling does it reverse, and then by a few percent, which can be hundreds of participants because plans there are very large), and the reason is worth understanding because it decides which planner is correct for your design.

The regression planner, with random_predictors = FALSE, is exact when the predictor is fixed by design: the researcher chooses the values of the predictor, as in an experiment with assigned doses. A correlation study is different. It samples both variables, and the sampled spread of the predictor is itself random from study to study. Conditional on the spread that happens to be drawn, the test statistic is a noncentral F, but its noncentrality parameter rises and falls with that spread, which follows a chi square distribution on N - 1 degrees of freedom. The distribution that actually governs the statistic is therefore a mixture of noncentral F distributions, not a single one, and a mixture is more spread out than its center. The two frames give the same estimates, and the same distribution when there is no effect, but different power when there is one (Sampson, 1974; Gatsonis & Sampson, 1989; Kelley & Maxwell, 2003). Planning as though the statistic were a plain noncentral F overstates the planned study's power by roughly 1/n for a plan of n participants: about .01 at n = 100, .02 at 50, .05 at 20, and .10 at 10 (for the regression plan of 119 above, about .009). It therefore returns a sample size that is a few participants short; the correct plan is typically 2 to 5 participants larger, and 1 to 3 larger for the smallest plans.

ss_buc_correlation() uses the mixture. It has an exact closed form (the Poisson weights that build the noncentral F become negative binomial weights when mixed over the chi square), so the returned sample size is exact for bivariate normal data rather than an approximation. The rule of thumb is simple:

If both variables are sampled, as in an observational correlation study, use ss_buc_correlation().
If the predictor is fixed by the design, the effect is a regression slope, and ss_buc_reg_coef() with p = 1 and random_predictors = FALSE is the planner for it. With random_predictors = TRUE the same call returns exactly what ss_buc_correlation() returns, and the multiple regression planners use the same mixture for several sampled predictors.

You may supply the prior result either as the correlation r_observed or as its t statistic t_observed; the two are equivalent.

Assuming Less Publication Bias

alpha_prior is the significance level a prior result had to reach to be published, so it encodes how strongly publication bias is assumed to have shaped the prior literature. The default of .05 assumes the field publishes results significant at the .05 level; a larger value assumes a more permissive literature and therefore a smaller correction. A prior study that escaped the file drawer entirely, such as a registered report or a pilot, is described by alpha_prior = 1: no p value exceeds 1, so that means no publication filter at all, and only the correction for uncertainty, set by assurance, remains.

alphas <- c(.05, .20, 1)
sizes <- sapply(alphas, function(ap) {
  plan <- ss_buc_correlation(r_observed = .35, N = 100, alpha_prior = ap)
  plan$value[plan$term == "necessary_N"]
})
data.frame(alpha_prior = alphas, necessary_N = sizes)

Auditing the Plan by Simulation

The plan above rests on a promise about the long run: at assurance = .80, about 80% of studies planned this way should reach the target power. That promise is easy to state and hard to picture. ss_buc_sensitivity() makes it concrete by simulating the very process the correction models. It draws many prior studies from a true effect you name, keeps only the ones a literature would have published, re-plans from each survivor, and reports how often those plans actually attain the target.

You must name the true effect; there is no default, because letting the method assume the corrected estimate is true would be asking it to grade its own work. A natural first choice is the plan's own corrected parameter, which asks "if the correction happened to be exactly right, how would plans built this way behave?"

plan <- ss_buc_correlation(r_observed = .35, N = 100, assurance = .80)
ncp <- plan$value[plan$term == "ncp_adjusted"]
audit <- ss_buc_sensitivity(plan, true_ncp = ncp, replications = 2000,
                            seed = 2017)
audit

Two attainment rates are reported, and the difference between them is the point of the exercise.

attainment_with_refusals is the fraction of published prior studies for which a plan built this way reaches the target, counting a refusal (a prior result too weak to plan from, which the function declines rather than guesses at) as a success, because declining is the most cautious plan of all. This is the quantity the assurance guarantee is about, and it lands on the requested assurance of .80 within simulation error: here it is r apa(row_value(audit, "attainment_with_refusals")), with a Monte Carlo standard error of r apa(row_value(audit, "attainment_with_refusals_mcse"), 2).
attainment is the same fraction among only the studies that received a plan. It is lower, r apa(row_value(audit, "attainment")) here (Monte Carlo standard error r apa(row_value(audit, "attainment_mcse"), 2)), because it excludes the refusals, and it is what a researcher who did get a plan actually experiences.

The two are linked by the refusal rate reported alongside them:

attainment_with_refusals = refusal_rate + (1 - refusal_rate) * attainment

Here that is r apa(row_value(audit, "refusal_rate")) + (1 - r apa(row_value(audit, "refusal_rate"))) times r apa(row_value(audit, "attainment")), or r apa(row_value(audit, "attainment_with_refusals")); the refusal rate's own Monte Carlo standard error is r apa(sqrt(row_value(audit, "refusal_rate") * (1 - row_value(audit, "refusal_rate")) / row_value(audit, "replications")), 2). The two rates coincide when nothing is refused and separate as refusals rise. That separation is not a flaw in the method; it is the visible price of a method that refuses to plan from a prior result too weak to support the request, rather than returning a confident but hollow number. size_at_true_effect is the sample size the true effect really requires, the benchmark each simulated plan is measured against, and the quartiles show how widely the plans scatter around it.

The most instructive use of the function is to vary the assumed truth. Ask what happens if the real correlation is smaller than the corrected estimate, which is the case a cautious planner most wants to guard against:

smaller <- ss_buc_sensitivity(plan, true_ncp = ncp * 0.6, replications = 2000,
                              seed = 2017)
rows <- c("attainment_with_refusals", "attainment", "refusal_rate",
          "size_at_true_effect")
data.frame(term = rows, value = smaller$value[match(rows, smaller$term)])

A smaller truth needs a larger study, r row_value(smaller, "size_at_true_effect") participants rather than r row_value(plan, "necessary_N"), yet attainment_with_refusals holds at r apa(row_value(smaller, "attainment_with_refusals")) (Monte Carlo standard error r apa(row_value(smaller, "attainment_with_refusals_mcse"), 2)): the guarantee does not depend on the corrected estimate being right. What changes is how it is kept. More of the published prior results are now too weak to plan from, so the refusal rate rises to r apa(row_value(smaller, "refusal_rate")) (Monte Carlo standard error r apa(sqrt(row_value(smaller, "refusal_rate") * (1 - row_value(smaller, "refusal_rate")) / row_value(smaller, "replications")), 2)), and among the plans actually issued only r apa(row_value(smaller, "attainment")) reach the target (Monte Carlo standard error r apa(row_value(smaller, "attainment_mcse"), 2)).

When the Prior Result Cannot Support the Request

Every prior result has a largest assurance it can support, the closed-form value 1 - p/alpha_prior, where p is the prior study's p value. Ask for more and the correction drives the noncentrality parameter to zero: the prior effect cannot be distinguished from nothing at that level of confidence, so no sample size attains the target and the function stops rather than returning a meaningless answer. The strong prior used throughout this vignette supports assurance up to about .99, so this does not arise for it, but a weaker prior reaches its ceiling quickly. plot() on any plan draws the ceiling and the way the required sample size climbs toward it; see vignette("understanding-assurance", package = "BUCSS") for that figure and a fuller discussion. A refusal is informative in its own right: it says the prior study cannot anchor a replication aimed at detecting the effect, and Anderson and Kelley (2024) describe the other goals a replication can have, each with its own sample size logic.

References

Anderson, S. F., & Kelley, K. (2024). Sample size planning for replication studies: The devil is in the design. Psychological Methods, 29, 844--867. https://doi.org/10.1037/met0000520

Anderson, S. F., Kelley, K., & Maxwell, S. E. (2017). Sample-size planning for more accurate statistical power: A method adjusting sample effect sizes for publication bias and uncertainty. Psychological Science, 28, 1547--1562. https://doi.org/10.1177/0956797617723724

Gatsonis, C., & Sampson, A. R. (1989). Multiple correlation: Exact power and sample size calculations. Psychological Bulletin, 106, 516--524. https://doi.org/10.1037/0033-2909.106.3.516

Kelley, K., & Maxwell, S. E. (2003). Sample size for multiple regression: Obtaining regression coefficients that are accurate, not simply significant. Psychological Methods, 8, 305--321. https://doi.org/10.1037/1082-989X.8.3.305

Sampson, A. R. (1974). A tale of two regressions. Journal of the American Statistical Association, 69, 682--689. https://doi.org/10.1080/01621459.1974.10480189

Taylor, D. J., & Muller, K. E. (1996). Bias in linear model power and sample size calculation due to estimating noncentrality. Communications in Statistics: Theory and Methods, 25, 1595--1610. https://doi.org/10.1080/03610929608831787