One-Way ANOVA

What one-way ANOVA is and what question it answers

One-way analysis of variance (ANOVA) answers a single omnibus question: "Is there any difference in the population mean of outcome Y across k ≥ 3 independent groups?" The test does not specify which groups differ or by how much — it only tells you whether the pattern of group means is more extreme than chance would produce if all groups were drawn from the same population. That omnibus structure is both its strength (controls the family-wise type I error for the overall comparison) and its limitation (a significant F requires follow-up post-hoc comparisons to be actionable).

Variance decomposition: between vs within

ANOVA rests on one algebraic identity: total variability in the data can be split exactly into two non-overlapping components.

• Between-group sum of squares (SS_between): measures how much the group means scatter around the grand mean. Formally: SS_between = Σ_i n_i (ȳ_i − ȳ_grand)², summing over k groups. If all groups have the same mean, SS_between = 0. Large SS_between signals that group membership explains a meaningful share of the outcome's variability.

• Within-group sum of squares (SS_within): measures how much individuals scatter around their own group mean. Formally: SS_within = Σ_i Σ_j (y_ij − ȳ_i)². This is the residual variability that remains after accounting for group membership — the "noise" against which the group signal is compared.

• Degrees of freedom: df_between = k − 1 (groups minus one); df_within = N − k (total observations minus number of groups). The mean squares are MS_between = SS_between / df_between and MS_within = SS_within / df_within.

• The F ratio: F = MS_between / MS_within. Under the null hypothesis that all group means are equal, the ratio MS_between / MS_within follows an F-distribution with (k−1, N−k) degrees of freedom. A large F means the between-group signal is large relative to the within-group noise; the p-value is the probability of observing an F this extreme or larger by chance.

• Additive identity check: SS_total = SS_between + SS_within, where SS_total = Σ_i Σ_j (y_ij − ȳ_grand)². This must hold exactly (within floating-point precision) in any correct implementation and serves as an internal audit check.

Assumptions and when they are most critical

ANOVA makes three classical assumptions:

1. observations are independent within and across groups;
2. within each group the outcome is approximately normally distributed;
3. all groups share the same population variance (homoscedasticity). In practice, the robustness of these assumptions varies by sample size and violation type:

• Independence: this is the most critical and least recoverable assumption. If patients are clustered within hospitals or practices, or if the same patient appears in multiple groups (repeated measures), standard ANOVA produces invalid inferences regardless of sample size. Mixed models or repeated-measures ANOVA are the correct alternatives.

• Normality: at large n (typically ≥ 20–30 per group), the central limit theorem ensures that the sampling distribution of the group mean is approximately normal even if individual observations are skewed, so the F-test remains valid. At small n with severely non-normal distributions, the Kruskal-Wallis H test is the nonparametric alternative.

• Homoscedasticity (equal variances): Welch ANOVA (also called the Welch F-test or oneway.test in R) relaxes this assumption by using group-specific variance estimates and adjusting degrees of freedom via the Welch-Satterthwaite equation, analogous to Welch's t-test for two groups. As a practical rule, check Levene's test as a diagnostic — not as a decision gate — and use Welch ANOVA when group variances differ substantially (ratio > 4:1) or when group sizes are unequal. The Brown-Forsythe test is an alternative diagnostic that is more robust to non-normality than Levene's test.

The post-hoc problem and multiplicity control

A significant omnibus F answers "yes, there is a difference somewhere" but provides no information about which pairs differ. Post-hoc pairwise comparisons are required, and they must account for the multiplicity created by testing all C(k,2) = k(k−1)/2 pairs simultaneously.

• Tukey HSD (Honestly Significant Difference): the standard choice for balanced designs (equal n per group). Controls the family-wise error rate (FWER) across all pairwise comparisons at exactly α. Produces simultaneous confidence intervals for all pairwise mean differences. Preferred when all pairwise comparisons are of interest.

• Bonferroni correction: divides the α level by the number of comparisons tested. Simple and applicable to any test family (not limited to pairwise means), but conservative — it understates the power of the multiple comparison procedure relative to Tukey HSD. Appropriate when only a pre-specified subset of comparisons is of interest or when mixed comparison types are needed.

• Holm-Bonferroni: sequentially-rejective version of Bonferroni, uniformly more powerful while still controlling FWER. A nearly free upgrade over plain Bonferroni in most software.

• Dunnett's test: designed specifically for comparing k−1 treatment groups to a single control group, with higher power than Tukey or Bonferroni for that specific structure.

The catalog's general principle on multiplicity applies directly: the appropriate correction depends on the inferential structure of the comparison (all pairwise vs control-vs-treatment vs pre-specified hypotheses) and the tolerable trade-off between FWER and power. Never apply post-hoc tests without a significant omnibus F; doing so inflates type I error.

Effect size and estimands beyond the p-value

The F-test p-value alone is insufficient for HEOR reporting. The appropriate effect estimates are:

• η² (eta-squared): SS_between / SS_total. The proportion of total variance explained by group membership. Easy to compute but biased upward in small samples.
• ω² (omega-squared): a bias-corrected alternative: ω² = (SS_between − df_between × MS_within) / (SS_total + MS_within). Preferred for reporting in confirmatory work.
• Partial η²: used in factorial ANOVA where multiple factors are present; equals SS_between / (SS_between + SS_within) for one-way ANOVA (= η²).
• Pairwise mean differences with CIs: the most interpretable quantities for decision-making. Tukey HSD produces simultaneous 95% CIs for all k(k−1)/2 pairwise differences. Each CI answers: "What is the plausible range of the true mean difference between these two groups?"

Why ANOVA is rarely the endpoint in RWE: confounding and the extension to regression

In a randomized trial, group membership is assigned by the study design, so ANOVA on the raw outcome estimates the causal effect of assignment. In observational RWE, groups typically correspond to self-selected treatment choices — first-line vs second-line vs third-line therapy, or different drug classes chosen by different types of patients. The group means differ for reasons that include both treatment effectiveness and the confounded characteristics of patients who receive each treatment.

An unadjusted one-way ANOVA comparing outcomes across treatment groups in a claims or EHR dataset is therefore descriptive, not causal. It describes the association between group membership and the outcome but cannot distinguish the treatment effect from confounding by indication, disease severity, prior treatment history, or provider practice patterns.

The natural extension is the general linear model (ANCOVA), which includes group indicators as categorical predictors alongside continuous and categorical covariates. In the language of regression: one-way ANOVA is the special case of OLS regression with only group dummy variables. Adding covariates moves the analysis toward adjusted mean differences that better account for confounding, though regression adjustment alone does not fully remove confounding from observational comparisons — propensity-score methods, instrumental variables, or target trial emulation are required when confounding by indication is substantial.

For HEOR specifically: comparing total costs, utilization rates, or quality-adjusted life years across 3+ lines of therapy requires at minimum covariate adjustment, and usually weighting or matching, before the group comparison is interpretable as anything close to a treatment effect. Document explicitly whether the ANOVA is unadjusted (descriptive, appropriate for Table 2 or supplementary exploratory analysis) or whether regression-adjusted means (LS means) from a GLM are being presented as a controlled comparison.

Welch ANOVA for heteroscedastic groups

When group variances differ — which is common in HEOR when comparing groups of patients with different disease severity or different lengths of follow-up — the standard one-way ANOVA F-test has inflated type I error if the group with the largest variance also has the smallest sample size (and deflated type I error in the reverse configuration). Welch ANOVA corrects for this by using a weighted combination of within-group variance estimates and adjusting the degrees of freedom. In R, oneway.test() implements Welch ANOVA by default (var.equal=FALSE); in Python, scipy's f_oneway() does NOT perform the Welch correction — use pingouin's welch_anova() or statsmodels for the Welch version. In SAS, PROC GLM produces the standard ANOVA; use the WELCH option in PROC GLM or PROC MIXED for the heteroscedastic case.

HEOR usage patterns: costs, utilization, and quality of life across groups

One-way ANOVA appears in HEOR for:

• Comparing mean total annual healthcare costs across 3+ lines of therapy or therapeutic areas (with the caveat that GLM with gamma/Tweedie link is the modern standard for cost inference).
• Comparing mean utilization counts (hospitalizations, ED visits, specialist visits) across geographic regions, payer types, or disease severity categories.
• Comparing mean quality-of-life scores (EQ-5D utility, PROMIS scores) across treatment groups in registry or survey studies.
• Comparing mean biomarker levels or laboratory values across categories.

A recurring HEOR challenge is that cost and utilization distributions are right-skewed with heavy tails, violating the normality assumption at small n. At large n (typical for claims), the CLT protects the ANOVA F-test for inference on means, but the GLM framework (gamma with log link for costs; negative binomial for counts) remains preferred because it:

1. respects the distributional shape,
2. naturally accommodates a log-scale mean ratio interpretation, and
3. extends seamlessly to covariate adjustment.

Clustered data and the independence violation in site-level studies

Many HEOR datasets have a hierarchical structure: patients are nested within physicians, who are nested within hospitals or health systems. When groups are defined by site-level factors (geographic region, hospital type, payer), patients within the same site are correlated — they share unmeasured site-level characteristics (practice patterns, formulary, local referral networks). This within-site correlation violates the independence assumption of one-way ANOVA and causes the standard F-test to underestimate the true variance, inflating the type I error rate.

Appropriate alternatives for clustered data:

• Cluster-robust standard errors (adjusting the variance estimator without changing the point estimate) when the number of clusters is large (≥ 30–50) and cluster sizes are balanced.
• Linear mixed models with a random intercept for cluster (site, practice, or physician) as the preferred full solution — this explicitly models the within-cluster correlation structure and produces valid inference even at moderate cluster counts.

Pros, cons, and trade-offs

Pros:

• The F-test controls the family-wise type I error for the omnibus multi-group comparison in a single test, avoiding the inflated error of k(k-1)/2 uncorrected pairwise t-tests.
• Directly interpretable effect estimates (mean differences, η²) are available.
• Well-understood power properties; sample size formulas are established for both the omnibus test and Tukey HSD post-hoc comparisons.
• Straightforward extension to factorial ANOVA (multiple grouping factors) and ANCOVA (covariate adjustment) within the general linear model framework.
• Computationally trivial; universally available in all statistical software.
• Robust to normality violations at large n (CLT protects the F-statistic).
• Welch ANOVA extends robustness to variance heterogeneity.

Cons:

• Omnibus F does not identify which pairs differ; post-hoc comparisons are always required for actionable conclusions, adding complexity and reducing power per comparison.
• Assumes independence; invalid for repeated measures, matched cohorts, or clustered data.
• Standard ANOVA assumes equal variances; Welch correction must be explicitly requested.
• For skewed outcomes (costs, utilization counts), a GLM is preferred for primary inference on the mean when n is small-to-moderate.
• In confounded observational data, group mean comparisons are descriptive, not causal.

When NOT to use

Do not use one-way ANOVA as the primary inferential method when:

• Repeated measures on the same patients: using standard ANOVA on measurements taken from the same patients at multiple time points violates independence and inflates type I error. Use repeated-measures ANOVA, mixed-effects models (MMRM), or GEE instead.
• Clustered data without correction: patients nested within sites, practices, or geographic units violate independence. Use cluster-robust standard errors or mixed models with random intercepts for the cluster level.
• Confounded group comparisons presented as causal claims: an unadjusted ANOVA across treatment groups in an observational dataset is not a valid estimate of the treatment effect. Route to propensity-score methods, g-methods, or regression adjustment for causal inference.
• Heavy skew with small samples: at n < 15–20 per group with a severely right-skewed or heavy- tailed distribution (common for costs in small pilot studies), the CLT has not fully kicked in and the F-test's normality assumption is load-bearing. Use Kruskal-Wallis H as the nonparametric alternative, or a GLM with an appropriate distributional family.
• Binary or count outcomes as the primary estimand: if the outcome is binary (response yes/no) or a count (number of events), logistic regression or Poisson/negative-binomial regression produce direct odds ratios or rate ratios with CIs. ANOVA on a binary 0/1 outcome is arithmetically equivalent to a linear probability model, which can produce predicted probabilities outside [0,1] and is not the preferred model for binary outcomes.
• When only two groups exist: one-way ANOVA with k=2 reduces to the pooled t-test (F = t² exactly). Use Welch's t-test directly for two groups; it is simpler and the Welch correction is the default in most software.
• When the research question is about spread, not means: ANOVA tests equality of means. If the question is whether variability differs across groups (e.g., consistency of a biomarker), Levene's test or Bartlett's test addresses variance heterogeneity directly.

Interpreting the output

In the worked example, patients across three COPD severity groups (Mild, Moderate, Severe, four patients each) had mean specialist visit counts of 5, 9, and 13, respectively. The between-group mean square is 64 and the within-group mean square is approximately 6.667, yielding F = 9.6 on df = (2, 9) with p ≈ 0.006. Applying Tukey's HSD post-hoc test, the critical mean difference is approximately 5.10 visits. The Mild vs Severe comparison (8 visits) exceeds this threshold and is significant; the Mild vs Moderate and Moderate vs Severe comparisons (4 visits each) do not.

(1) Formal interpretation. The F ratio of 9.6 indicates that between-group variance is 9.6 times larger than pooled within-group variance. Under the null hypothesis that all three population means are equal, an F this extreme or larger on df = (2, 9) arises by chance in approximately 0.6% of samples (p ≈ 0.006). The omnibus test provides evidence against the null but does not itself identify which pairs differ. The Tukey HSD post-hoc comparisons — which simultaneously control family-wise type I error across all three pairwise tests — identify only the Mild vs Severe contrast as significantly different at the adjusted threshold.

(2) Practical interpretation. Patients with severe COPD averaged about 8 more specialist visits per year than patients with mild COPD — a gap large enough to survive the Tukey multiplicity correction. The adjacent pairs (Mild vs Moderate and Moderate vs Severe) show differences of 4 visits each that did not reach significance at n = 4 per group; the study is underpowered for those intermediate comparisons, though 4 visits per year is not necessarily clinically unimportant. Because this is an unadjusted observational comparison, differences in age, comorbidities, and care-seeking behavior across severity groups likely contribute to the observed gradient. The ANOVA result is descriptive; the gradient may partially or fully reflect confounding rather than a direct effect of COPD severity on specialist utilization.

Decision diagram

flowchart TD
  A[Outcome: continuous<br/>Groups: k ≥ 3 independent] --> B{Independence<br/>assumption met?}
  B -- No: repeated<br/>measures / clustering --> C[Mixed model / MMRM<br/>or cluster-robust SE]
  B -- Yes --> D{Equal variances<br/>across groups?}
  D -- Yes / unknown --> E["Fisher one-way ANOVA<br/>aov() / PROC GLM"]
  D -- No or unequal n --> F["Welch ANOVA<br/>oneway.test() / pingouin<br/>PROC GLM WELCH"]
  E --> G{Omnibus F<br/>significant?}
  F --> G
  G -- No --> H[Stop: no evidence of<br/>any group difference]
  G -- Yes --> I{Post-hoc<br/>structure?}
  I -- All pairwise<br/>equal variances --> J[Tukey HSD<br/>TukeyHSD / LSMEANS TUKEY]
  I -- All pairwise<br/>unequal variances --> K[Games-Howell<br/>rstatix / SAS macro]
  I -- Treat vs control<br/>only --> L["Dunnett's test<br/>PROC GLM DUNNETT"]
  I -- Pre-specified<br/>subset --> M["Bonferroni / Holm<br/>(p.adjust in R; PROC MULTTEST)"]
  J --> N[Report mean differences<br/>with 95% CI per pair<br/>and effect size omega-sq]
  K --> N
  L --> N
  M --> N

Worked example

Scenario

A health outcomes analyst at a regional payer is examining whether the number of specialist visits in a 6-month period differs across patients with mild, moderate, and severe chronic obstructive pulmonary disease (COPD) enrolled in the same health plan. Four patients from each severity tier are randomly selected from the claims database. The analyst wants to know whether any severity tier has a meaningfully different mean visit count before proceeding to a fully adjusted regression model. She runs a one-way ANOVA and manually walks through the variance decomposition to understand the F ratio.

Dataset

Number of specialist visits in 6 months by COPD severity tier (n=4 per group). All values are counts extracted directly from outpatient claims. Grand mean across all 12 patients = 9.

Steps

Result

Group means: Mild=20/4=5, Moderate=36/4=9, Severe=52/4=13. Grand mean=108/12=9. SS_between=128, SS_within=60, SS_total=128+60=188. df_between=2, df_within=9. MS_between=128/2=64, MS_within=60/9=6.667. F=64/(60/9)=576/60=9.6 > F_crit(2,9)=4.26; reject the null. Tukey HSD post-hoc: only Mild vs Severe is significant (difference=8 > HSD=5.10); Mild vs Moderate and Moderate vs Severe do not reach significance at alpha=0.05. Interpretation: COPD severity is associated with specialist visit frequency, but the Mild vs Moderate and Moderate vs Severe contrasts are not individually detectable at this sample size. A larger sample or an adjusted regression model is needed before drawing clinical conclusions — and this unadjusted comparison may reflect patient case-mix differences as much as any direct effect of severity classification on care-seeking behavior.

Trade-offs

Runnable example

import numpy as np
from scipy import stats
from statsmodels.stats.multicomp import pairwise_tukeyhsd

# ── Dataset: specialist visits by COPD severity (n=4 per group) ──
mild     = [2, 4, 6, 8]    # group mean = 5
moderate = [6, 8, 10, 12]  # group mean = 9
severe   = [10, 12, 14, 16] # group mean = 13

# ── 1. One-way ANOVA (scipy: Fisher / equal-variance version) ──
f_stat, p_val = stats.f_oneway(mild, moderate, severe)
print(f"One-way ANOVA: F={f_stat:.4f}, p={p_val:.4f}")
# Expected: F=9.6, p=0.0057

# ── 2. Manual SS decomposition (matches worked example) ──
all_vals = mild + moderate + severe
grand_mean = np.mean(all_vals)
groups = [mild, moderate, severe]
group_means = [np.mean(g) for g in groups]
n_per = [len(g) for g in groups]

ss_between = sum(n * (gm - grand_mean)**2 for n, gm in zip(n_per, group_means))
ss_within  = sum(sum((y - gm)**2 for y in g) for g, gm in zip(groups, group_means))
ss_total   = sum((y - grand_mean)**2 for y in all_vals)
print(f"\nSS_between={ss_between:.1f}, SS_within={ss_within:.1f}, SS_total={ss_total:.1f}")
print(f"Decomposition check: {ss_between:.1f} + {ss_within:.1f} = {ss_between+ss_within:.1f}")
# Expected: 128.0 + 60.0 = 188.0

k = len(groups); N = len(all_vals)
df_between = k - 1; df_within = N - k
ms_between = ss_between / df_between
ms_within  = ss_within  / df_within
f_manual   = ms_between / ms_within
print(f"MS_between={ms_between:.3f}, MS_within={ms_within:.3f}, F={f_manual:.4f}")
# Expected: MS_between=64.000, MS_within=6.667, F=9.6000

# ── 3. Effect size: eta-squared and omega-squared ──
eta_sq  = ss_between / ss_total
omega_sq = (ss_between - df_between * ms_within) / (ss_total + ms_within)
print(f"\nEta-squared = {eta_sq:.3f}, Omega-squared = {omega_sq:.3f}")

# ── 4. Tukey HSD post-hoc (pairwise comparisons with FWER control) ──
import pandas as pd
data_flat  = np.array(all_vals)
groups_lbl = (["Mild"]*4 + ["Moderate"]*4 + ["Severe"]*4)
tukey = pairwise_tukeyhsd(data_flat, groups_lbl, alpha=0.05)
print("\nTukey HSD post-hoc:")
print(tukey.summary())
# Mild-Severe significant; Mild-Moderate and Moderate-Severe not significant

# ── 5. Welch ANOVA (robust to unequal variances) ──
# scipy does not implement Welch ANOVA natively; use pingouin if available
try:
    import pingouin as pg
    df_long = pd.DataFrame({
        "visits": all_vals,
        "group":  groups_lbl
    })
    welch_res = pg.welch_anova(data=df_long, dv="visits", between="group")
    print("\nWelch ANOVA (pingouin):")
    print(welch_res[["Source", "ddof1", "ddof2", "F", "p-unc", "np2"]])
except ImportError:
    print("\npingouin not installed — install with: pip install pingouin")

# ── 6. Levene's test for equal variances (diagnostic, not a decision gate) ──
levene_stat, levene_p = stats.levene(*groups)
print(f"\nLevene's test for equal variances: W={levene_stat:.3f}, p={levene_p:.4f}")
print("(Diagnostic only — do not use to choose between Fisher and Welch ANOVA.)")

# ── Dataset: specialist visits by COPD severity (n=4 per group) ──
visits <- c(2, 4, 6, 8,      # Mild
            6, 8, 10, 12,    # Moderate
            10, 12, 14, 16)  # Severe
group  <- factor(rep(c("Mild", "Moderate", "Severe"), each = 4),
                 levels = c("Mild", "Moderate", "Severe"))

# ── 1. Fisher one-way ANOVA (assumes equal variances) ──
anova_fit <- aov(visits ~ group)
cat("Fisher one-way ANOVA (aov):\n")
print(summary(anova_fit))
# Expected: F=9.6, p=0.0057

# ── 2. Welch one-way ANOVA (DEFAULT in R's oneway.test; does NOT assume equal variances) ──
cat("\nWelch one-way ANOVA (oneway.test, var.equal=FALSE — the R default):\n")
print(oneway.test(visits ~ group))
# var.equal=FALSE is the default; note the adjusted (non-integer) df for the denominator

cat("\nFisher version via oneway.test (var.equal=TRUE, equivalent to aov above):\n")
print(oneway.test(visits ~ group, var.equal = TRUE))

# ── 3. Manual SS decomposition (matches worked example) ──
grand_mean  <- mean(visits)
group_means <- tapply(visits, group, mean)
n_per       <- tapply(visits, group, length)

ss_between <- sum(n_per * (group_means - grand_mean)^2)
ss_within  <- sum(tapply(seq_along(visits), group, function(idx) {
  gm <- mean(visits[idx]); sum((visits[idx] - gm)^2)
}))
ss_total   <- sum((visits - grand_mean)^2)
df_between <- length(levels(group)) - 1
df_within  <- length(visits) - length(levels(group))
ms_between <- ss_between / df_between
ms_within  <- ss_within  / df_within
f_manual   <- ms_between / ms_within

cat(sprintf("\nManual: SS_between=%.1f, SS_within=%.1f, SS_total=%.1f\n",
            ss_between, ss_within, ss_total))
cat(sprintf("MS_between=%.3f, MS_within=%.3f, F=%.4f\n",
            ms_between, ms_within, f_manual))

# ── 4. Effect sizes ──
eta_sq  <- ss_between / ss_total
omega_sq <- (ss_between - df_between * ms_within) / (ss_total + ms_within)
cat(sprintf("Eta-squared=%.3f, Omega-squared=%.3f\n", eta_sq, omega_sq))

# ── 5. Tukey HSD post-hoc (Fisher ANOVA; equal variances assumed) ──
cat("\nTukey HSD post-hoc:\n")
print(TukeyHSD(anova_fit))
# Mild-Severe should be significant; Mild-Moderate and Moderate-Severe not significant

# ── 6. Games-Howell post-hoc (for Welch ANOVA; does NOT assume equal variances) ──
# Requires rstatix package
if (requireNamespace("rstatix", quietly = TRUE)) {
  library(rstatix)
  df_data <- data.frame(visits = visits, group = group)
  cat("\nGames-Howell post-hoc (use with Welch ANOVA):\n")
  print(games_howell_test(df_data, visits ~ group))
} else {
  cat("\nInstall rstatix for Games-Howell: install.packages('rstatix')\n")
}

# ── 7. Levene's test (diagnostic for variance heterogeneity) ──
if (requireNamespace("car", quietly = TRUE)) {
  cat("\nLevene's test (diagnostic only):\n")
  print(car::leveneTest(visits ~ group))
}

/* ── Create dataset: specialist visits by COPD severity (n=4 per group) ── */
data work.copd;
  input severity $ visits;
  datalines;
Mild     2
Mild     4
Mild     6
Mild     8
Moderate  6
Moderate  8
Moderate 10
Moderate 12
Severe   10
Severe   12
Severe   14
Severe   16
;
run;

/* ── 1. PROC ANOVA (balanced, equal n per group only) ── */
proc anova data=work.copd;
  class severity;
  model visits = severity;
  /* ANOVA (not GLM) is appropriate here only because n=4 per group (balanced design) */
  means severity / tukey;   /* Tukey HSD post-hoc pairwise comparisons */
run;
quit;

/* ── 2. PROC GLM: one-way ANOVA, Tukey HSD, LS means, and WELCH option ── */
proc glm data=work.copd;
  class severity;
  model visits = severity;
  /* Fisher one-way ANOVA F-test                                                       */
  means severity / hovtest=levene;  /* Levene's test for equal variances (diagnostic)  */
  means severity / tukey alpha=0.05; /* Tukey HSD pairwise comparisons (equal-var)    */
  lsmeans severity / pdiff=all cl adjust=tukey; /* LS means with Tukey-adjusted CIs   */
  /* WELCH option: Welch F-test robust to unequal variances                            */
  means severity / welch;
run;
quit;

/* ── 3. Manual SS verification (PROC MEANS + DATA step) ── */
proc means data=work.copd mean var n noprint;
  class severity;
  var visits;
  output out=work.gstats mean=gm n=n var=gv;
run;

proc means data=work.copd mean noprint;
  var visits;
  output out=work.grandmean mean=grand_mean;
run;

/* SS_between = sum(n_i * (mean_i - grand_mean)^2) */
/* SS_within  = sum((n_i - 1) * var_i)             */
data work.anova_manual;
  set work.gstats(where=(_type_=1));
  if _n_=1 then set work.grandmean(keep=grand_mean);
  ss_between_i = n * (gm - grand_mean)**2;
  ss_within_i  = (n - 1) * gv;
run;

proc means data=work.anova_manual sum noprint;
  var ss_between_i ss_within_i;
  output out=work.anova_totals sum=ss_between ss_within;
run;

data work.anova_results;
  set work.anova_totals;
  k = 3; N = 12;
  df_between = k - 1;
  df_within  = N - k;
  ms_between = ss_between / df_between;
  ms_within  = ss_within  / df_within;
  f_stat     = ms_between / ms_within;
  ss_total   = ss_between + ss_within;
  eta_sq     = ss_between / ss_total;
  omega_sq   = (ss_between - df_between * ms_within) / (ss_total + ms_within);
  put "SS_between=" ss_between "SS_within=" ss_within "SS_total=" ss_total;
  put "MS_between=" ms_between "MS_within=" ms_within "F=" f_stat;
  put "Eta-squared=" eta_sq "Omega-squared=" omega_sq;
run;

Citations

1. [1]Altman DG, Bland JM. Statistics Notes: Comparing several groups using analysis of variance. BMJ. 1996;312(7044):1472-1473.
2. [2]du Prel JB, Röhrig B, Hommel G, Blettner M. Choosing Statistical Tests. Deutsches Ärzteblatt International. 2010;107(19):343-348.

1. [3]Delacre M, Lakens D, Leys C. Why psychologists should by default use Welch's t-test instead of Student's t-test. International Review of Social Psychology. 2017;30(1):92-101.

1. [4]Bland JM, Altman DG. Multiple significance tests: the Bonferroni method. BMJ. 1995;310(6973):170.