---
title: "Worked Examples Across Designs"
author: "Ken Kelley, Samantha F. Anderson, and Scott E. Maxwell"
output: rmarkdown::html_vignette
vignette: >
  %\VignetteIndexEntry{Worked Examples Across Designs}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---

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

Which planner fits a given prior study, and what does it need to know about
that study? This vignette walks through one example of each of the eighteen
`ss_buc_*` planners. The pattern is always the same: summarize the prior study
with a single test statistic and its design, and the function returns the
necessary planned-study sample size and the bias and uncertainty adjusted
noncentrality parameter (Anderson, Kelley, & Maxwell, 2017; Anderson & Kelley,
2024). Every example uses 80% power and 80% assurance unless noted. For the meaning of the `assurance` and
`alpha_prior` settings, see
`vignette("understanding-assurance", package = "BUCSS")`.

## Independent *t* Test

A prior two-group study reported *t* = 3.00 with 20 participants per group. The
returned sample size is per group.

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

`ss_buc_smd` is the same planner under an effect size name. It still takes
*t*, so convert a reported standardized mean difference *d* first:
`t = d / sqrt(1/n_1 + 1/n_2)`, with `n_1` and `n_2` the two group sizes.

## Dependent (Paired) *t* Test

A prior within-pairs study of 40 pairs reported *t* = 3.00. The returned sample
size is the number of pairs.

```{r}
ss_buc_paired_t(t_observed = 3, N = 40, alpha_prior = .05, alpha_planned = .05,
                assurance = .80, desired_power = .80)
```

The alias `ss_buc_smd_paired` serves the standardized mean difference framing of
the same comparison.

## One-Sample *t* Test

A prior study of 30 participants compared their mean with a fixed reference
value, such as a scale midpoint or a published norm, and reported *t* = 3.20.
The returned sample size is the total number of participants.

```{r}
ss_buc_one_sample_t(t_observed = 3.2, N = 30, alpha_prior = .05,
                    alpha_planned = .05, assurance = .80, desired_power = .80)
```

A paired *t* test is a one-sample *t* test on the difference scores, so this
planner is the paired planner under another name: for the same prior *t* and
*N* the two return the same plan.

## Between-Subjects ANOVA

A one-way fully between-subjects design takes the number of levels of the single
factor. Here a four-group design with total *N* = 120 reported *F* = 5.00 for
the omnibus effect. The returned sample size is per group.

```{r}
ss_buc_one_way_anova(F_observed = 5, N = 120, levels_A = 4, alpha_prior = .05,
                     alpha_planned = .05, assurance = .80, desired_power = .80)
```

A two-way fully between-subjects design takes the levels of both factors and the
effect of interest (`factor_A`, `factor_B`, or `interaction`). Here a 2 x 3
design with total *N* = 120 reported *F* = 5.00 for the main effect of Factor B.
The returned sample size is per cell.

```{r}
ss_buc_factorial_anova(F_observed = 5, N = 120, levels_A = 2, levels_B = 3,
                       effect = "factor_B", alpha_prior = .05,
                       alpha_planned = .05, assurance = .80, desired_power = .80)
```

When the effect of interest is not a standard main effect or interaction, or the
design has more than two factors, use the general form and supply the degrees of
freedom directly. The denominator degrees of freedom equal *N* minus the number
of cells when the cell means are the only estimated parameters, so supplying
`df_denominator = N - cells` (here 120 - 6 = 114) reproduces the two-way result
above. A smaller `df_denominator` accounts for covariates or other nuisance
parameters and calls for a larger sample.

```{r}
ss_buc_factorial_anova_general(F_observed = 5, N = 120, cells = 6,
                               df_numerator = 2, df_denominator = 114,
                               alpha_prior = .05, alpha_planned = .05,
                               assurance = .80, desired_power = .80)
```

## Within-Subjects ANOVA

For a fully within-subjects design, supply the factor levels and the effect.
Here a 2 x 3 within-subjects design with *N* = 60 reported *F* = 5.00 for the
main effect of Factor B. The returned sample size is the total number of
subjects.

```{r}
ss_buc_rm_anova(F_observed = 5, N = 60, levels_A = 2, levels_B = 3,
                effect = "factor_B", alpha_prior = .05, alpha_planned = .05,
                assurance = .80, desired_power = .80)
```

The general form takes the numerator degrees of freedom directly and works for
any number of within-subjects factors. This example uses `assurance = .50`,
which corrects for publication bias only.

```{r}
ss_buc_rm_anova_general(F_observed = 6.5, N = 80, df_numerator = 1,
                        alpha_prior = .05, alpha_planned = .05, assurance = .50,
                        desired_power = .80)
```

## Split-Plot (Mixed) ANOVA

A split-plot design crosses a between-subjects factor with a within-subjects
factor. Supply the levels of each and the effect (`between`, `within`, or
`interaction`). Here a design with 2 between-subjects levels and 3
within-subjects levels and *N* = 60 reported *F* = 5.00 for the within-subjects
effect. The returned sample size is per between-subjects cell.

```{r}
ss_buc_mixed_anova(F_observed = 5, N = 60, levels_between = 2, levels_within = 3,
                   effect = "within", alpha_prior = .05, alpha_planned = .05,
                   assurance = .80, desired_power = .80)
```

The general form takes the numerator degrees of freedom and the number of
between-subjects groups directly. The effect is one of `between_only`,
`within_only`, or `between_within`; supply `df_num_within`, the numerator
degrees of freedom of the within-subjects component, only when
`effect = "between_within"` (it is required there, and ignored with a warning
otherwise). The returned sample size is per group.

```{r}
ss_buc_mixed_anova_general(F_observed = 5, N = 90, df_numerator = 2,
                           num_groups = 3, effect = "between_only",
                           alpha_prior = .05, alpha_planned = .05,
                           assurance = .80, desired_power = .80)
```

## Departures From Sphericity

The within-subjects and split-plot planners assume sphericity by default, as
the uncorrected within-subjects test does, and echo that assumption as an
`epsilon` row of 1. When a paper reports a Greenhouse-Geisser (1959) or
Huynh-Feldt (1976) correction instead, it shows fractional degrees of freedom:
a three-level within-subjects factor measured on 30 subjects might be reported
as *F*(1.73, 50.2) = 6.50. The uncorrected degrees of freedom are 2 and 58, so
the correction is `epsilon = 1.73/2`, about .865. The denominator agrees
(50.2/58 is also about .865), and 50.2/1.73 = 29.0 is one less than the 30
subjects, which confirms how many subjects the analysis used.

```{r}
ss_buc_rm_anova(F_observed = 6.5, N = 30, levels_A = 3, epsilon = 1.73 / 2)
```

Ignoring the correction plans 92 subjects rather than 162:

```{r}
ss_buc_rm_anova(F_observed = 6.5, N = 30, levels_A = 3)
```

The difference has two sources. The same *F* on fewer degrees of freedom is
weaker evidence: its *p* value rises from .0028 to .0045, the largest
supportable assurance falls from .94 to .90, and the correction for
uncertainty at an `assurance` of .80 takes correspondingly more off. And the
planned study, analyzed with the same correction, has less power at any size.
Both follow the approximation of Muller and Barton (1989), in which both
degrees of freedom and the noncentrality parameter are multiplied by
`epsilon`; `ncp_adjusted` stays on the sphericity scale, so the two plans'
parameters are comparable. The plan is sensitive to `epsilon`, so it is worth
seeing across the values you find plausible:

```{r}
eps <- c(1, .95, .9, .85, .8, .75)
sapply(eps, function(e) {
  res <- ss_buc_rm_anova(F_observed = 6.5, N = 30, levels_A = 3, epsilon = e)
  res$value[res$term == "necessary_N"]
})
```

One `epsilon` describes both studies: it is a property of the covariance of
the repeated measures, which the planned study shares with the prior one.
When a paper reports both corrections, use the Greenhouse-Geisser value. The
correction is an approximation, and the help page for `ss_buc_rm_anova()`
reports how well it holds up in simulation. In a split-plot design it applies
to the within-subjects effect and the interaction, not to the between-subjects
effect, whose test is on subject means.

## Multiple Regression

Anderson (2021) is a step-by-step tutorial for planning regression studies with
BUCSS, and is the best starting point for the three regression planners below.

For a single regression coefficient, supply the *t* statistic for that
coefficient, the prior sample size, and the number of predictors `p`. A prior
study with *N* = 150 and 3 predictors reported *t* = 3.00 for one coefficient.
The returned sample size is the total.

```{r}
ss_buc_reg_coef(t_observed = 3, N = 150, p = 3, alpha_prior = .05,
                alpha_planned = .05, assurance = .80, desired_power = .80)
```

For the omnibus test of the model R-squared, supply the model *F* statistic.
Here a prior study with *N* = 150 and 4 predictors reported *F* = 5.00.

```{r}
ss_buc_R2(F_observed = 5, N = 150, p = 4, alpha_prior = .05,
          alpha_planned = .05, assurance = .80, desired_power = .80)
```

For a joint test of a subset of predictors, also supply `p_joint`, the number of
predictors tested together. Here 2 of the 4 predictors are tested jointly.

```{r}
ss_buc_reg_joint(F_observed = 5, N = 150, p = 4, p_joint = 2,
                 alpha_prior = .05, alpha_planned = .05, assurance = .80,
                 desired_power = .80)
```

All three regression planners take `random_predictors`. By default
(`random_predictors = FALSE`) they treat the predictors as fixed by design, as
in an experiment in which the researcher sets their values, the frame of BUCSS
1.x. When the predictors are sampled along with the outcome, as in most
observational research (Kelley & Maxwell, 2003), pass
`random_predictors = TRUE`: the test statistic's noncentrality parameter then
rises and falls with the sampled spread of the tested predictors, and the
planners use the exact distribution that follows (Sampson, 1974; Gatsonis &
Sampson, 1989).

```{r}
ss_buc_R2(F_observed = 5, N = 150, p = 4, random_predictors = TRUE)
```

Planned as though they were fixed, sampled predictors make the reported
`actual_power` overstate the planned study's power by an amount that shrinks
as 1/*n* for a plan of *n* participants and grows with the number of
predictors left untested (roughly 1/*n* with one predictor and 4/*n* for a
coefficient among eight: about .04 at *n* = 100, and less than .01 for each of
the three plans above). From prior studies as large as these the two frames
differ little: with `random_predictors = TRUE` the three plans above are 623,
239, and 3,898 rather than 624, 234, and 3,960, the joint test's differing
most because it is requested close to its assurance ceiling of .84. From a
small prior study, with predictors left untested and a plan several times the
prior *N*, the random-predictor plan can be smaller by up to about the number
of untested predictors plus one, divided by the prior *N*: as much as 40
percent for a coefficient among eight predictors from a prior study of 20
(the help pages of `ss_buc_reg_coef()` and `ss_buc_reg_joint()` describe the
pattern). With one predictor,
`ss_buc_reg_coef(p = 1, random_predictors = TRUE)` returns exactly what
`ss_buc_correlation()`, described in the next section, returns.

## Correlation

A prior study reporting a Pearson correlation is planned from the correlation
itself. A study with *N* = 100 reported *r* = .35. The returned sample size is
the total.

```{r}
ss_buc_correlation(r_observed = .35, N = 100, alpha_prior = .05,
                   alpha_planned = .05, assurance = .80, desired_power = .80)
```

The test of a correlation is arithmetically the test of a simple regression
slope, so it is tempting to plan it as one. That is exact only when the
predictor is fixed by design. A correlation study samples both variables, and
then the statistic's noncentrality parameter is itself random, proportional to
the sampled spread of the predictor. The resulting distribution is more
dispersed, and planning as though it were not overstates the power. The
planner uses the exact distribution for sampled variables, so it asks for a
few more participants than the regression planner does when that planner
treats the predictor as fixed:

```{r}
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"])
```

If your predictor really is fixed by design, the analysis is a regression
slope, and `ss_buc_reg_coef()` with `p = 1` and `random_predictors = FALSE` is
the planner for it.
`vignette("correlation", package = "BUCSS")` develops the comparison.

## Analysis of Covariance

Covariates cost error degrees of freedom in both the prior study and the
planned one, so the planner takes their number directly. A prior study with
*N* = 120, three groups, and two covariates reported *F* = 6.00 for the group
effect. The returned sample size is per cell.

```{r}
ss_buc_ancova(F_observed = 6, N = 120, cells = 3, n_covariates = 2,
              alpha_prior = .05, alpha_planned = .05, assurance = .80,
              desired_power = .80)
```

With `n_covariates = 0` this reproduces `ss_buc_one_way_anova()` exactly. For a
single-degree-of-freedom contrast among the adjusted means, supply that
contrast's *F* and set `df_numerator = 1`.

## Multivariate Comparison of Two Groups

When two groups are compared on several outcome variables at once, the prior
study's Hotelling *T* squared (Hotelling, 1931), or the *F* it maps to, plans
the replication.
A study with *N* = 80 comparing two groups on 3 outcomes reported *F* = 4.50.
The returned sample size is per group.

```{r}
ss_buc_manova(F_observed = 4.5, N = 80, p_variables = 3, alpha_prior = .05,
              alpha_planned = .05, assurance = .80, desired_power = .80)
```

This covers the multivariate hypothesis whose rank is one, which is what a
two-group comparison gives. A hypothesis of higher rank has a set of
noncentrality parameters rather than a single one, and the correction does not
apply to it.

## Nested Model Comparison

An effect tested by the chi square difference between two nested models, such
as a constrained path in a structural equation model, is planned from that
difference. Supply both models' degrees of freedom rather than only their
difference. A prior study with *N* = 250 fit a model with 42 degrees of
freedom and a constrained version with 43, a difference of 9.50 on 1 degree of
freedom.

```{r}
ss_buc_chisq_diff(chisq_observed = 9.5, N = 250, df_full = 42, df_restricted = 43,
                  alpha_prior = .05, alpha_planned = .05, assurance = .80,
                  desired_power = .80)
```

This is for the *difference* test, not for an omnibus model fit chi square: a
publishable fit statistic is a small one, which inverts the selection the
correction assumes. Requiring both models' degrees of freedom is what lets the
function enforce that. An omnibus fit chi square compares a model against the
*saturated* model, which has zero degrees of freedom, so `df_full = 0` is its
signature and the function refuses it:

```{r error = TRUE}
ss_buc_chisq_diff(chisq_observed = 9.5, N = 250, df_full = 0, df_restricted = 43)
```

A just-identified full model also has zero degrees of freedom, and there the
difference test is a genuine one; the message above, and the function's help
page, say how to plan that case.

## Unequal Variances

When the two groups' variances differ, the Welch test is the appropriate
comparison, and planning it needs the groups' relative spread in addition to
their sizes. A prior study with 40 and 55 participants, whose second group was
half again as variable as the first, reported *t* = 3.00.

```{r}
ss_buc_welch_t(t_observed = 3, n_1 = 40, n_2 = 55, sd_ratio = 1.5,
               alpha_prior = .05, alpha_planned = .05, assurance = .80,
               desired_power = .80)
```

This planner treats the Welch statistic as noncentral *t* on the
Welch-Satterthwaite degrees of freedom (Satterthwaite, 1946; Welch, 1947),
the usual approximation in power analysis for this test. A minimum planned
sample size of 10 per group is imposed because of estimation accuracy issues
at smaller sample sizes: with fewer per group the approximation overstates
the real test's power, while from 10 per group it is within .006 of it at an
`alpha_planned` of .05. At a
stricter `alpha_planned` the power of the smallest plans is overstated by
more (at 10 per group by up to about .012 at .01 and .025 at .001), and the
.006 holds again from 14 and 20 per group, respectively. When the target
power is reached with fewer than 10 per group, the planner returns 10.

Only the ratio enters the computation, but a paper usually reports the group
standard deviations or variances rather than their ratio, so those may be
given directly instead. All three calls describe the same study:

```{r}
by_ratio <- ss_buc_welch_t(t_observed = 3, n_1 = 40, n_2 = 55, sd_ratio = 1.5)
by_sds <- ss_buc_welch_t(t_observed = 3, n_1 = 40, n_2 = 55, sd_1 = 8.4,
                         sd_2 = 12.6)
by_vars <- ss_buc_welch_t(t_observed = 3, n_1 = 40, n_2 = 55, var_1 = 70.56,
                          var_2 = 158.76)
size_term <- "necessary_n_per_group"
c(ratio = by_ratio$value[by_ratio$term == size_term],
  sds   = by_sds$value[by_sds$term == size_term],
  vars  = by_vars$value[by_vars$term == size_term])
```

The assumed ratio matters more than any other input here, so it is worth
seeing how the plan moves across the values you find plausible:

```{r}
ratios <- c(1, 1.5, 2, 2.5, 3)
sapply(ratios, function(r) {
  res <- ss_buc_welch_t(t_observed = 3, n_1 = 40, n_2 = 55, sd_ratio = r)
  res$value[res$term == "necessary_n_per_group"]
})
```

One `sd_ratio` describes both studies: the prior study, whose statistic is
read under it, and the planned study. Getting it wrong can be costly. The plan
built assuming equal variances (213 per group) delivers about .76 power rather
than .80 if the ratio is in fact 2 in both studies. It delivers about .43 if
the prior study's groups were equally variable but the planned study's ratio
is 2, with the mean difference held fixed in units of the first group's
standard deviation. Re-running across the ratios you find plausible, as above,
and planning for the least favorable one you find credible guards against the
first kind of error. The second, a planned study more unequal than the prior
one, is not something the prior result can speak to, and needs its own
allowance.

## References

Anderson, S. F. (2021). Using prior information to plan appropriately powered
regression studies: A tutorial using BUCSS. *Psychological Methods, 26*,
513--526. https://doi.org/10.1037/met0000366

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

Greenhouse, S. W., & Geisser, S. (1959). On methods in the analysis of profile
data. *Psychometrika, 24*, 95--112. https://doi.org/10.1007/BF02289823

Hotelling, H. (1931). The generalization of Student's ratio. *The Annals of
Mathematical Statistics, 2*, 360--378. https://doi.org/10.1214/aoms/1177732979

Huynh, H., & Feldt, L. S. (1976). Estimation of the Box correction for degrees
of freedom from sample data in randomized block and split-plot designs.
*Journal of Educational Statistics, 1*, 69--82.
https://doi.org/10.3102/10769986001001069

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

Muller, K. E., & Barton, C. N. (1989). Approximate power for repeated-measures
ANOVA lacking sphericity. *Journal of the American Statistical Association,
84*, 549--555. https://doi.org/10.1080/01621459.1989.10478802

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

Satterthwaite, F. E. (1946). An approximate distribution of estimates of
variance components. *Biometrics Bulletin, 2*, 110--114.
https://doi.org/10.2307/3002019

Welch, B. L. (1947). The generalization of 'Student's' problem when several
different population variances are involved. *Biometrika, 34*, 28--35.
https://doi.org/10.1093/biomet/34.1-2.28