Planning a Correlation Study

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. Doing so gives a slightly different, and slightly too small, answer:

t_equivalent <- .35 * sqrt(98 / (1 - .35^2))
c(correlation = ss_buc_correlation(t_observed = t_equivalent, N = 100)$value[1],
  regression  = ss_buc_reg_coef(t_observed = t_equivalent, N = 100, p = 1)$value[1])
#> correlation  regression 
#>         122         119

The gap is small here, a few participants, but it is systematic, and the reason is worth understanding because it decides which planner is correct for your design.

The regression planner 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. Planning as though the statistic were a plain noncentral F overstates the planned study’s power, by up to about .04, and so returns a sample size that is a few participants short.

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 is the planner for it.

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 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.

alphas <- c(.05, .20, 1)
sizes <- sapply(alphas, function(ap) {
  ss_buc_correlation(r_observed = .35, N = 100, alpha_prior = ap)$value[1]
})
data.frame(alpha_prior = alphas, necessary_N = sizes)
#>   alpha_prior necessary_N
#> 1        0.05         122
#> 2        0.20         108
#> 3        1.00         103

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"]
ss_buc_sensitivity(plan, true_ncp = ncp, replications = 2000, seed = 2017)
#>  term                     value 
#>  true_ncp                 6.71  
#>  ncp_adjusted             6.71  
#>  replications             2000  
#>  publication_rate         0.718 
#>  refusal_rate             0.324 
#>  attainment               0.708 
#>  attainment_mcse          0.0124
#>  attainment_with_refusals 0.803 
#>  size_at_true_effect      122   
#>  size_q25                 110   
#>  size_median              212   
#>  size_q75                 564   
#> 
#> Design: Pearson correlation
#> Sample size unit: total
#> Requested assurance: .80
#> Target power: .80

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 near the requested assurance of .80.
• attainment is the same fraction among only the studies that received a plan. It is lower, because it excludes the refusals, and it is what a researcher who did get a plan actually experiences.

The gap between the two is the refusal rate, reported alongside them. A large gap 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)
smaller$value[smaller$term %in% c("attainment", "attainment_with_refusals",
                                  "refusal_rate", "size_at_true_effect")]
#> [1]   0.457500   0.637788   0.803500 200.000000

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.

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

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.