---
title: "Getting Started With BUCSS"
author: "Ken Kelley, Samantha F. Anderson, and Scott E. Maxwell"
output: rmarkdown::html_vignette
vignette: >
  %\VignetteIndexEntry{Getting Started With BUCSS}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---

```{r setup, include = FALSE}
knitr::opts_chunk$set(collapse = TRUE, comment = "#>")
```

## Why Not Plan From the Observed Effect Size?

A natural way to plan a new study is to take the effect size from a prior study
and compute the sample size that would give, say, 80% power for an effect of
that size. The trouble is that the prior effect size is rarely a trustworthy
target, for two reasons.

* **Publication bias.** The published literature is filtered. Studies that
  reached statistical significance are far more likely to appear in print than
  those that did not, so the effects you can read about are systematically
  larger than the effects that were actually run. Planning from a published
  estimate inherits that upward bias.
* **Uncertainty.** Even without any filtering, a single study's effect size is
  one draw from a sampling distribution. Treating that one number as if it were
  the population value ignores how much it could have come out differently.

Both forces push the same way: they make the prior effect look bigger and more
certain than it is, so a study planned at face value tends to be underpowered.

The uncertainty piece is easy to underestimate. Even with no publication bias,
the sample size needed for a target power is a nonlinear function of the effect
size, so an overestimate reduces power (by planning too small a study) by more
than an equally sized underestimate raises it (by planning a larger one).
Symmetric error in the estimate therefore produces an asymmetric loss of power,
which is why correcting for uncertainty, not just bias, matters.

## What BUCSS Does

BUCSS ("Bias and Uncertainty Corrected Sample Size") corrects the prior study's
effect before planning. It uses the prior study's test statistic and design to
build the likelihood of a truncated noncentral distribution, truncated because
small, nonsignificant effects go unpublished (Taylor & Muller, 1996), then reads
off a corrected noncentrality parameter at a chosen level of *assurance*. That corrected
noncentrality parameter, rather than the face-value estimate, drives the sample
size calculation for the planned study. The method is described in Anderson,
Kelley, and Maxwell (2017), and Anderson and Kelley (2024) extend it to the
design of replication studies.

```{r}
library(BUCSS)
```

## A First Example

Suppose a prior study used an independent samples *t* test with 20 participants
per group and reported *t* = 3.00. We want to plan a new study with 80% power,
and we want 80% assurance, that is, we want the realized power to reach or
exceed 80% in 80% of replications.

```{r}
result <- ss_buc_independent_t(t_observed = 3, n = 20, alpha_prior = .05,
                         alpha_planned = .05, assurance = .80, desired_power = .80)
result
```

The correction calls for 130 participants per group, far more than a face-value
calculation from *t* = 3.00 would suggest.

## Reading the Result

Every `ss_buc_*` planner returns a tidy object shaped like the result tables
in the **DMAR** package: `term`/`value` rows holding the design results first
(the necessary sample size, named for its unit; the implied total sample size
where the unit is per group or per cell; the achieved power at the returned
size; the bias and uncertainty adjusted noncentrality parameter; and the
planning effect size, the effect size the plan is built on), followed by rows
echoing the planning inputs, so the assumptions travel with the result.
The print method shows the aligned table and closes with the design, for the
regression planners the predictor frame, the unit, and the largest assurance
the prior result can support. Because the object is
an ordinary `data.frame` underneath, you can pull any quantity out directly:

```{r}
result$value[result$term == "necessary_n_per_group"]
result$value[result$term == "ncp_adjusted"]
result$value[result$term == "smd_planning"]
```

The last of these, `smd_planning`, is the **planning effect size**, an
adjusted effect size: the value, in the design's usual effect size metric (here
a standardized mean difference), used to plan the sample size given the
situation at hand (the prior study's result, the publication bias assumed
through `alpha_prior`, and the uncertainty allowed for through `assurance`). It
is the adjusted noncentrality parameter restated in that metric, and it is what
the plan is built on, not an estimate of the population effect to report: above
an assurance of .5 it is deliberately conservative (a lower confidence limit),
and at .5 it is the median unbiased estimate. A study of 130 per group has
exactly the returned `actual_power` when the standardized mean difference is
`r format(signif(result$value[result$term == "smd_planning"], 3))`. The prior
study's own standardized mean difference was 3 * sqrt(2/20), about
`r format(signif(3 * sqrt(2 / 20), 3))`. Each planner names the row for its
metric (`f_planning` for Cohen's *f* in the ANOVA planners, `rho_planning` for
a correlation, and so on); `vignette("reading-output", package = "BUCSS")`
lists them.

The unit the sample size is counted in depends on the design (per group, per
cell, number of pairs, or total): the size row is named accordingly
(`necessary_n_per_group`, `necessary_n_per_cell`, or `necessary_N`), and the
unit also travels on the object so you never have to guess:

```{r}
attr(result, "sample_size_unit")
```

For programmatic use the broom verbs give the convenient views. `tidy()`
returns the compact estimate view (`term`, `estimate`, and `power`: the sample
size and its achieved power), as the **DMAR** package's sample size planners
do, and `glance()` returns the one-row wide view, which begins with those same
three columns and then keeps every row of the result as a column.

```{r}
tidy(result)
tidy(result)$estimate
glance(result)
```

For a manuscript's method section, `planning_sentence()` writes the sentence:

```{r}
planning_sentence(result)
```

## The Replication Example

When the prior study has unequal group sizes, you can pass both. Here a prior
study reported *t* = 3.00 with group sizes of 50 and 55, and we ask for 90%
assurance:

```{r}
replication <- ss_buc_independent_t(t_observed = 3, n = c(50, 55),
                                    alpha_prior = .05, alpha_planned = .05,
                                    desired_power = .80, assurance = .90)
replication
```

The corrected plan calls for 1485 per group, a vivid illustration of how much
publication bias and uncertainty can inflate the sample size needed for a
dependable replication. The published standardized mean difference was
3 * sqrt(1/50 + 1/55), about `r format(signif(3 * sqrt(1 / 50 + 1 / 55), 3))`,
while the plan is built on a planning effect size of
`r format(signif(replication$value[replication$term == "smd_planning"], 3))`:
at 90% assurance, after allowing for publication at the .05 level, that is the
effect the replication must be powered to detect. BUCSS 1.x reported 1482 for this example under its grid
search; as of BUCSS 2.0.0 the adjusted noncentrality parameter is found by
direct root finding, which gives 1485.

## When the Prior Result Is Not Enough

The correction can also decline to plan. The prior study of the first example
supports assurance up to about .90, the value its printed footer reports; ask
for .95 and the function stops rather than returning a meaningless number, and
its message says how far to step back:

```{r, error = TRUE}
ss_buc_independent_t(t_observed = 3, n = 20, assurance = .95)
```

`vignette("understanding-assurance", package = "BUCSS")` explains where that
ceiling comes from and draws it.

## The Function Family

BUCSS covers eighteen designs, one function per analysis. Each summarizes the
prior study with a single test statistic plus the design.

| Function | Analysis | Prior statistic |
|----------|----------|-----------------|
| `ss_buc_one_sample_t` | one-sample *t* test | `t_observed` |
| `ss_buc_independent_t` (alias `ss_buc_smd`) | independent *t* test / standardized mean difference | `t_observed` |
| `ss_buc_paired_t` (alias `ss_buc_smd_paired`) | dependent (paired) *t* test | `t_observed` |
| `ss_buc_welch_t` | Welch (unequal variance) *t* test | `t_observed` |
| `ss_buc_one_way_anova` | one-way between-subjects ANOVA | `F_observed` |
| `ss_buc_factorial_anova` | two-way between-subjects ANOVA | `F_observed` |
| `ss_buc_factorial_anova_general` | between-subjects ANOVA, any factors or contrast | `F_observed` |
| `ss_buc_ancova` | analysis of covariance | `F_observed` |
| `ss_buc_rm_anova` | within-subjects ANOVA (one-way or two-way) | `F_observed` |
| `ss_buc_rm_anova_general` | within-subjects ANOVA, any factors | `F_observed` |
| `ss_buc_mixed_anova` | split-plot (mixed) ANOVA | `F_observed` |
| `ss_buc_mixed_anova_general` | split-plot ANOVA, any factors | `F_observed` |
| `ss_buc_correlation` | Pearson correlation | `r_observed` (or `t_observed`) |
| `ss_buc_reg_coef` | single regression coefficient | `t_observed` |
| `ss_buc_R2` | omnibus model R-squared | `F_observed` |
| `ss_buc_reg_joint` | joint test of a subset of predictors | `F_observed` |
| `ss_buc_manova` | two-group multivariate comparison (Hotelling *T* squared) | `F_observed` (or `T2_observed`) |
| `ss_buc_chisq_diff` | nested model chi square difference test | `chisq_observed` |

Two further functions take a finished plan from any of these planners:

| Function | What it does |
|----------|--------------|
| `ss_buc_sensitivity()` | audits the plan by simulation: from a true effect you name, how often a plan built this way reaches the target power |
| `plot()` | draws the plan's assurance ceiling: the corrected noncentrality parameter and the sample size as assurance rises toward the largest value the prior result supports |

See `vignette("worked-examples", package = "BUCSS")` for one worked example of
each planner, `vignette("understanding-assurance", package = "BUCSS")` for how
the `assurance` and `alpha_prior` settings shape the correction,
`vignette("reading-output", package = "BUCSS")` for the print, `tidy()`, and
`glance()` views of the result,
`vignette("correlation", package = "BUCSS")` for why a correlation study needs
its own planner and how `ss_buc_sensitivity()` audits a plan, and
`vignette("planning-replications", package = "BUCSS")` for a worked
demonstration of why a literature of underpowered studies misleads and how the
correction repairs the replication sample size.

## 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.
https://doi.org/10.1080/03610929608831787