contrast_test: Tests One or More Contrasts of Group Means in a One-Way Design

Description

Given a fitted one-way aov or lm object and a set of contrast weights, computes for every contrast the estimate $\hat{\psi} = \sum_i c_i \bar{Y}_i$, its standard error, t-statistic, degrees of freedom, two-sided p-value, and confidence interval. Supports several common multiple-comparison adjustments and either equal-variance (pooled) or Welch-style unequal-variance inference.

Usage

contrast_test(
  object,
  contrasts = "pairwise",
  adjust = "none",
  conf_level = 0.95,
  var_equal = TRUE
)

Value

A data.frame with one row per contrast and columns

contrast, estimate, se, t, df,

p_value, p_adjusted, ci_lower, and

ci_upper. The adjustment, confidence level, and variance assumption are stored as

attr(*, "adjust"), attr(*, "conf_level"), and

attr(*, "var_equal"). The table prints through the

dmar_tbl display layer and works with

tidy and glance (see

dmar_tidiers).

Arguments

object

A fitted aov or lm object for a one-way design (a single grouping factor on the right-hand side of the formula).

contrasts

Specification of one or more contrasts. Any of:

"pairwise" (default): All pairwise comparisons among the group means.

a named list of numeric vectors

Each vector is one contrast and its name is used as the row label.

a numeric matrix

Each row is one contrast; rownames, if present, are used as labels.

a numeric vector

Treated as a single contrast.

Each contrast vector must have length equal to the number of groups, and the weights are typically chosen to sum to zero.

adjust

Multiple-comparison adjustment. One of "none" (default), "bonferroni", "scheffe", "tukey" (pairwise contrasts only), or any of the sequential methods supported by p.adjust ("holm", "hochberg", "BH", "BY").

conf_level

Confidence level for the interval (default 0.95).

var_equal

Logical. If TRUE (default), uses the pooled error variance $\mathit{MS}_{\text{error}}$ and the residual degrees of freedom from object. If FALSE, uses each group's own sample variance and a Welch-Satterthwaite approximate $df$ per contrast.

Author

Ken Kelley [email protected]

Details

Test statistic. For a contrast with weights $c_1, \ldots, c_k$ ($k$ = number of groups), the estimate is $\hat{\psi} = \sum_i c_i \bar{Y}_i$. Under equal variances, the standard error is $\sqrt{\mathit{MS}_{\text{error}} \sum_i c_i^2 / n_i}$ with $df = N - k$; under unequal variances, the standard error is $\sqrt{\sum_i c_i^2 s_i^2 / n_i}$ with the Welch-Satterthwaite df, $$df_{\text{Welch}} = \frac{\left(\sum_i c_i^2 s_i^2 / n_i\right)^2}{\sum_i (c_i^2 s_i^2 / n_i)^2 / (n_i - 1)}.$$ The unadjusted p-value is two-sided based on the t reference distribution.

Adjustments. The p_adjusted and confidence interval critical value are computed as follows.

"none": no adjustment; the CI uses $t_{1-\alpha/2,df}$.
"bonferroni": $p_{\text{adj}} = \min(1, m\, p)$ for $m$ contrasts, with CI based on $t_{1-\alpha/(2m),df}$.
"scheffe": appropriate for any contrast (or family of contrasts). $p_{\text{adj}}$ comes from the upper tail of an F reference distribution applied to $t^2 / (k-1)$, and the CI uses $\sqrt{(k-1)\, F_{1-\alpha,\,k-1,df}}$.
"tukey": requires every contrast to be pairwise. Uses the studentized range distribution (ptukey/qtukey) so that $p_{\text{adj}} = 1 - \mathrm{ptukey}(|t|\sqrt{2}; k, df)$ and the CI uses $q_{1-\alpha,\,k,df} / \sqrt{2}$.
"holm", "hochberg", "BH", "BY": p.adjust is applied to the unadjusted p-values; the CI uses the unadjusted $t$-critical value because these methods do not give simultaneous CIs in closed form.

Variance assumption with adjustments. The Tukey and Scheffé procedures assume equal variances; combining them with var_equal = FALSE is at the user's risk (the resulting Type I error rate is no longer guaranteed). For unequal variances, common alternatives are Games-Howell (Tukey-style) and Brown-Forsythe (Scheffé-style); these are not currently supported here.

Scope. Only one-way designs are supported in v1 (one outcome, one grouping factor). Multi-way designs throw an informative error.

References

Hsu, J. C. (1996). Multiple comparisons: Theory and methods. Chapman & Hall.

Maxwell, S. E., Delaney, H. D., & Kelley, K. (2027). Designing experiments and analyzing data: A model comparison perspective (4th ed.). Routledge.

Scheffe, H. (1953). A method for judging all contrasts in the analysis of variance. Biometrika, 40, 87--104.

Tukey, J. W. (1953). The problem of multiple comparisons. Unpublished manuscript, Princeton University.

Examples

Run this code

# All pairwise comparisons among the three arms of the depression_bdi
# treatment study.
fit <- aov(bdi_post ~ condition, data = depression_bdi)
contrast_test(fit, contrasts = "pairwise")

# Custom contrasts with a Tukey-protected family-wise error rate.
contrast_test(fit, contrasts = "pairwise", adjust = "tukey")

# A user-defined contrast: the SSRI arm vs. the average of the placebo
# and wait list arms. With levels ordered ssri, placebo, wait_list, the
# weights c(1, -0.5, -0.5) estimate that difference.
contrast_test(
  fit,
  contrasts = list("ssri vs non-drug arms" = c(1, -0.5, -0.5)),
  adjust = "scheffe"
)

# Welch-style inference: the wait list variance is about twice the
# SSRI variance, so the pooled error term is worth questioning.
contrast_test(fit, contrasts = "pairwise", var_equal = FALSE)

# Pairwise treatment comparisons in the Smith, Meyers, and Delaney
# (1998) drinking trial, on the normalizing log scale. Each row is
# one pairwise contrast of the three treatment means.
fit_drinks <- aov(log_drinks ~ treatment, data = drinks_trial)
contrast_test(fit_drinks, contrasts = "pairwise")

# An a priori contrast: the two active CRA arms (averaged) versus
# standard care. With levels ordered Standard, CRA, CRA + Disulfiram,
# the weights c(-1, 0.5, 0.5) compare the active arms against Standard.
contrast_test(
  fit_drinks,
  contrasts = list("CRA arms vs Standard" = c(-1, 0.5, 0.5))
)

Run the code above in your browser using DataLab