26  Testing Data Analysis Workflows

Program testing can be used to show the presence of bugs, but never to show their absence.

Edsger W. Dijkstra, Notes on Structured Programming (1970)

NoteSources

Adapted from author’s lecture notes and supporting materials for a graduate practicum in biostatistics.

26.1 Prerequisites

Answer the following questions to see if you can bypass this chapter. You can find the answers at the end of the chapter in Section 26.17.

  1. Why would you test a one-off data analysis script that will only ever run on a single dataset?
  2. What is the difference between a unit test, an integration test, and an end-to-end test in the context of a data analysis pipeline?
  3. What does testthat::expect_snapshot() capture, and when is a snapshot test preferable to a direct value-comparison test?

26.2 Learning objectives

By the end of this chapter you should be able to:

  • Write unit tests for analytic helper functions with testthat.
  • Add an integration test that exercises the full data pipeline on a small synthetic dataset.
  • Use expect_snapshot() to capture complex outputs (printed objects, plots, summary tables).
  • Run the test suite locally via testthat::test_local() and in CI via GitHub Actions.
  • Identify test smells (tests that pass for the wrong reason, tests that hide bugs, flaky tests).
  • Apply visual regression testing to ggplot output via vdiffr.

26.3 Orientation

Analysis code rarely gets tested, and the reason is a respectable one. Software is written to run on inputs its author has never seen; an analysis script is written to run on one dataset, once, and then to be thrown away. Testing something that runs once against known data can look like ceremony.

The reasoning has two holes. The first is that the script does not run once. It runs every time a data extract is refreshed, an exclusion is revised, or a reviewer asks for a sensitivity analysis, and by then its author has forgotten which of its branches were ever exercised. The second is more serious: the claim ‘it runs on one known dataset’ is doing no work, because the failures that matter in analysis code do not stop it running. A filter that removes the wrong rows, a join that duplicates a patient, a unit conversion applied twice, all produce output. The question a test answers is not ‘did it run’ but ‘would I have noticed if it were wrong’, and for most analysis scripts the honest answer is no.

We bring in this chapter the habits of software engineering to bear on that question.

The framework used in this chapter is testthat (3rd edition), because it is the framework the reader is most likely to meet in the wider R ecosystem and in the companion textbook (chapter 20 of Statistical Computing in the Age of AI, which covers it in the package-development context). This chapter focuses on applying it to non-package analysis code.

NoteA note on tinytest

The house standard for projects scaffolded with zzcollab is tinytest, not testthat. The two are close cousins: tinytest has no dependencies, keeps each test file runnable on its own, and uses the same expect_*() vocabulary. Everything in this chapter transfers with small substitutions, principally expect_equal() and expect_true() in place of their testthat namesakes, and tinytest::test_package() or run_test_dir() in place of testthat::test_local(). Where a project already depends on testthat, stay with it; where you are starting fresh under zzcollab, prefer tinytest. The discipline the chapter teaches, testing statistical content rather than plumbing, is identical under either framework.

26.4 The statistician’s contribution

Writing a test is five minutes of work. Writing a test that would actually fail if the code were wrong is the part that requires thought:

Test the statistical content, not just the plumbing. A test that ‘returns a tibble of the right shape’ is a structural check, not a test of correctness. A test that the regression coefficient on a known synthetic dataset matches the closed-form answer is a test of correctness.

Edge cases are where bugs live. Empty input, single observation, all-NA, perfect collinearity, boundary value (zero variance, all-zero counts, class with one observation). Test these deliberately.

Test the pipeline, not just the helpers. Unit tests catch helper-function bugs. Integration tests on synthetic data catch pipeline bugs, including ones where helpers work individually but the composition fails.

Don’t test what you don’t trust. Mocking out a function call inside a test of that function makes the test tautological. The test passes because you made it pass, not because the code is correct.

These judgments are what make tests useful rather than performative.

26.5 Why test analyses?

The Orientation argued the general case. Three specific occasions make it concrete, and each is a moment that arrives in ordinary projects rather than unusual ones.

The first is refactoring. You restructure a cleaning script to use dplyr where it previously used base::merge, and the output should be identical. Note what is peculiar about this situation: you have a precise expectation, you have no way to check it by looking, and the consequence of being wrong is a silently different cohort. A test that compares the output against a saved baseline settles the question in a second, and without one the honest position after a refactor is that you do not know whether it changed anything.

The second is upstream drift, and it is the case tests catch that nothing else does. Your script reads an extract from an institutional warehouse, and one quarter the warehouse begins reporting a duration in seconds where it previously reported days. Nothing errors. The column is still numeric, the join still succeeds, the model still fits, and every downstream number is wrong by a factor of 86,400. A test asserting that the column falls in a plausible range catches this on the first run after the change; no amount of reading your own code ever will, because your code is not what changed.

The third is documentation. A test of the form ‘given input X, this function returns Y’ is a statement about intent that cannot rot, because it is executed. Comments drift away from the code they describe; assertions cannot, since a false one fails the build.

Against these the cost is a few minutes per test, weighed against a silent miscalculation in something you have published under your name.

Tests come in layers, and the layers catch different failures. Figure 26.1 is the shape, with the important qualification that the pyramid familiar from software engineering is somewhat flatter for analysis code: the integration layer, a small synthetic dataset run through the whole pipeline, does more work here than it does in most applications.

flowchart TD
  U["<b>Unit</b><br/><i>one function, known input</i><br/>catches: helper bugs<br/>cost: seconds"]
  I["<b>Integration</b><br/><i>whole pipeline, 20-line<br/>synthetic dataset</i><br/>catches: composition bugs,<br/>silent filter and join errors<br/>cost: seconds"]
  E["<b>End to end</b><br/><i>raw data to final figures</i><br/>catches: everything left<br/>cost: minutes to hours"]
  U --> I --> E
Figure 26.1: The three test layers and what each one catches. Unit tests are fast and narrow; end-to-end tests are slow and broad; the integration layer, where a synthetic dataset is run through the whole cleaning pipeline, has the best cost-benefit ratio for analysis code and is the layer most often missing.

26.6 Unit tests with testthat

For a function in R/clean.R:

# R/clean.R
#' Compute age groups from numeric age
#' @param age numeric vector
age_group <- function(age) {
  cut(age,
      breaks = c(0, 18, 40, 65, Inf),
      right = FALSE,
      labels = c("under-18", "18-39", "40-64", "65+"))
}

Tests in tests/testthat/test-clean.R:

test_that("age_group bins typical adult ages correctly", {
  expect_equal(as.character(age_group(c(20, 50, 75))),
               c("18-39", "40-64", "65+"))
})

test_that("age_group respects boundary values", {
  expect_equal(as.character(age_group(c(18, 40, 65))),
               c("18-39", "40-64", "65+"))
  expect_equal(as.character(age_group(c(17, 39, 64))),
               c("under-18", "18-39", "40-64"))
})

test_that("age_group handles NA", {
  expect_true(is.na(age_group(NA)))
})

test_that("age_group errors on non-numeric input", {
  expect_error(age_group("twenty"))
})

Each test_that() is one named test with one or more expect_*() assertions. Run all tests:

testthat::test_local()
# or, if the project is a package:
devtools::test()

For non-package analysis projects, place tests in tests/testthat/ and source the relevant scripts at the top of the test file.

26.7 Integration tests on synthetic data

Unit tests catch helper-function bugs. They do not catch bugs where the helpers work individually but the pipeline composition fails.

Integration tests run the full pipeline on a small synthetic dataset:

# tests/testthat/test-pipeline.R

test_that("cleaning pipeline produces expected analytic dataset", {
  # synthetic data that exercises every branch
  raw <- tibble::tibble(
    patient_id = c(1, 2, 3, 4),
    age        = c(25, 50, 75, 17),     # one excluded
    bp_v1      = c(120, 140, 160, NA),
    bp_v2      = c(118, 138, NA,  NA),
    sex        = c("M", "F", "F", "M"),
    treatment  = c("placebo", "active", "placebo", "active")
  )

  clean <- run_pipeline(raw)

  # structural assertions
  expect_s3_class(clean, "tbl_df")
  expect_named(clean, c("patient_id", "age_group", "sex",
                        "treatment", "visit", "bp"))

  # row count: 3 adults x 2 visits, minus NAs
  expect_equal(nrow(clean), 5)    # patient 4 excluded; patient 3 has 1 NA bp

  # specific values
  expect_equal(clean[clean$patient_id == 1 & clean$visit == 1, "bp"][[1]],
               120)
})

The synthetic dataset is small (4 patients, 2 visits) but exercises each cleaning rule: age exclusion (patient 4), missing values (patient 3 visit 2). A test that runs in a fraction of a second covers the pipeline’s logic.

For each cleaning rule (filter, derive, pivot, recode), the synthetic dataset should have at least one row that exercises it.

Question. Your test asserts expect_equal(nrow(result), 100). The test passes. Does this mean the result is correct?

Answer.

No. The test asserts only that the result has 100 rows. It says nothing about which rows, what values, or whether the right rows were filtered. A bug that produces 100 rows of the wrong patients passes the test silently. To test correctness:

expect_equal(nrow(result), 100)            # structural
expect_setequal(result$patient_id,         # correctness
                expected_ids)
expect_equal(result$age, expected_ages)

The combination tests structure and content. The structural test alone is a false positive generator. This pattern, ‘looks right at the surface, wrong underneath’, is the most common testing failure mode in data analysis.

26.8 Snapshot tests

For complex outputs that are hard to assert literally, snapshots:

test_that("summary table renders correctly", {
  fit <- lm(mpg ~ wt + hp, data = mtcars)
  expect_snapshot(summary(fit))
})

test_that("Table 1 has expected structure", {
  expect_snapshot(gtsummary::tbl_summary(d, by = treatment))
})

On first run, the output is captured to tests/testthat/_snaps/. On subsequent runs, output is compared to the snapshot; if it differs, the test fails and you decide whether to accept (testthat::snapshot_accept()) or investigate.

Snapshots are useful when:

  • Output is complex and tedious to assert literally.
  • Output is human-readable formatting (tables, printed summaries, plot text).
  • You care more about ‘is this what was reviewed’ than ‘is this exactly the right value’.

Snapshots are not useful when:

  • Output changes for trivial reasons (whitespace, random seed, system locale). These produce false positives.
  • The snapshot is too large to read manually.

26.9 Visual regression testing with vdiffr

For ggplot output, vdiffr saves an SVG of the rendered plot and compares on subsequent runs:

library(vdiffr)

test_that("regression plot is unchanged", {
  fit <- lm(mpg ~ wt, data = mtcars)
  p <- ggplot(broom::augment(fit), aes(.fitted, .resid)) +
        geom_point() +
        geom_hline(yintercept = 0)
  expect_doppelganger("regression-residual-plot", p)
})

vdiffr handles cross-platform rendering issues (font differences, anti-aliasing) better than raw image comparison. Useful for catching unintended changes to plot code (an axis label change, a geom swap).

26.10 Continuous integration

usethis::use_github_action("check-standard") sets up GitHub Actions to run R CMD check on every push and pull request, on multiple OS / R-version combinations. (Older material calls this use_github_action_check_standard(), a function usethis has since removed in favor of naming the workflow as an argument.) For a non-package analysis project, the equivalent is a workflow that runs tests:

# .github/workflows/test.yaml
name: tests

on: [push, pull_request]

jobs:
  test:
    runs-on: ubuntu-latest
    steps:
      - uses: actions/checkout@v4
      - uses: r-lib/actions/setup-r@v2
      - uses: r-lib/actions/setup-r-dependencies@v2
        with:
          packages: |
            local::.
            any::testthat
      - name: Run tests
        run: Rscript -e 'testthat::test_local()'

Push the workflow; on every commit, GitHub runs the tests and reports pass/fail. Regressions caught at push time are far cheaper to fix than regressions caught in production.

26.11 Test smells

Common patterns that look like tests but produce false confidence:

Tautological tests.

test_that("foo works", {
  expect_equal(foo(1), foo(1))    # always passes
})

The test cannot fail because it tests foo() against itself. Useless.

Tests that mock the tested function.

test_that("fit_model returns a list", {
  with_mocked_bindings(
    code = expect_type(fit_model(d), "list"),
    fit_model = function(d) list(...)
  )
})

You replaced the function under test with a stub; the test passes because the stub returns a list, not because fit_model works.

Brittle snapshots.

test_that("output looks right", {
  expect_snapshot(do_thing(today()))     # date in output
})

The snapshot includes the current date; the test fails every day. Either remove the date from the output or use expect_snapshot_value with a date-tolerant comparison.

Overly permissive tolerances.

expect_equal(result, 0.5, tolerance = 0.5)   # passes for any value

A 100% tolerance is no test.

Tests that test the test framework, not the code.

expect_true(TRUE)

These are placeholders that should be filled in. Watch for them in code reviews.

26.12 Worked example: testing a cleaning pipeline

# R/cleaning.R
clean_visits <- function(raw) {
  raw |>
    janitor::clean_names() |>
    dplyr::filter(age >= 18) |>
    tidyr::pivot_longer(
      cols = dplyr::starts_with("bp_"),
      names_to = "visit",
      names_prefix = "bp_v",
      values_to = "bp"
    ) |>
    dplyr::filter(!is.na(bp))
}

# tests/testthat/test-cleaning.R
source("../../R/cleaning.R")

test_that("cleaning excludes minors", {
  raw <- tibble::tibble(
    patient_id = 1:3,
    age        = c(25, 17, 30),
    bp_v1      = c(120, 110, 130),
    bp_v2      = c(122, 112, 128)
  )
  result <- clean_visits(raw)
  expect_setequal(result$patient_id, c(1, 3))
})

test_that("cleaning drops NA blood pressures", {
  raw <- tibble::tibble(
    patient_id = c(1, 1, 2),
    age        = c(25, 25, 30),
    bp_v1      = c(120, NA, 130),
    bp_v2      = c(NA,  118, 128)
  )
  result <- clean_visits(raw)
  expect_equal(nrow(result), 4)        # 6 cells, 2 NA → 4
  expect_true(all(!is.na(result$bp)))
})

test_that("cleaning produces long format", {
  raw <- tibble::tibble(
    patient_id = 1,
    age        = 25,
    bp_v1      = 120,
    bp_v2      = 118
  )
  result <- clean_visits(raw)
  expect_named(result,
               c("patient_id", "age", "visit", "bp"),
               ignore.order = TRUE)
  expect_equal(nrow(result), 2)
})

The tests cover: the age filter, the NA filter, the pivot. Each test is a single named scenario; each asserts specific values, not just structure. Running the test suite catches any regression to the cleaning logic.

26.13 Collaborating with an LLM on tests

LLMs draft tests well; the judgment about what to test needs human input.

Prompt 1: drafting tests. Paste the function and ask: ‘write testthat tests covering happy path, edge cases, and error handling. Make the assertions specific to values, not just structure.’

What to watch for. Default LLM tests tend to be shape-checks. Push for value-checks. Edge cases: empty, NA, single observation, type mismatch.

Verification. Introduce a bug into the function and re-run the tests. If the bug is caught, the tests are useful. If not, add a more specific test.

Prompt 2: integration test on synthetic data. Describe the pipeline; ask the LLM to generate a small synthetic dataset that exercises every branch.

What to watch for. The LLM may generate a dataset that exercises the happy path but not the edges. Verify each branch (filter rules, missing-data paths, factor levels) is covered.

Verification. Add a deliberate bug to one branch; the test should fail. If it does not, the synthetic data does not exercise that branch.

Prompt 3: diagnosing a flaky test. Paste the test and the failure; ask: ‘what’s the source of non-determinism?’

What to watch for. Common causes: random number generation without seed, parallel processing without seed-aware RNG, system-time-dependent output, locale issues. The LLM should know these.

Verification. Apply the fix and run the test multiple times. If it stops being flaky, fixed.

26.14 Principle in use

Three habits separate a suite that catches bugs from one that reports green:

  1. Test the math, not just the shape. Structural checks alone produce false confidence.
  2. Integration tests on synthetic data. A 20-line synthetic dataset that exercises every branch catches pipeline bugs that unit tests miss.
  3. CI as a tripwire. GitHub Actions running tests on every push catches regressions at the commit, not at the publication.

26.15 Exercises

  1. Add unit tests for two functions from an existing analysis of yours. Aim for at least: a happy-path test, an edge-case test, and an error-path test.
  2. Write an integration test: a 20-line synthetic dataset that exercises every branch of your pipeline. Run the full pipeline in under a second.
  3. Set up GitHub Actions with a workflow that runs testthat::test_local() on every push. Push and verify the workflow runs green.
  4. Introduce a deliberate bug into your pipeline. Verify the test suite catches it. If it does not, add a test that does.
  5. Apply vdiffr::expect_doppelganger to one of your project’s plots. Modify the plot’s theme slightly; verify the test fails and you can accept or reject the change.
  6. Review a colleague’s test suite rather than their code. For each test, decide whether it checks statistical content or only structure, and find at least one test that would pass even if the function under test were subtly wrong. Propose the value-level assertion that would catch the bug, and deliver the review in writing as a pull-request comment (Chapter 8).

26.16 Further reading

  • (Wickham & Bryan, 2023) testing chapters, the canonical testthat reference.
  • vdiffr documentation on CRAN, visual regression testing for ggplot2.
  • The testthat package vignettes.

26.17 Prerequisites answers

  1. Even a one-off analysis benefits from tests: they catch regressions when you refactor mid-analysis, they serve as executable documentation of expected behavior, and they detect silent failures when upstream data changes format. The cost of writing a handful of tests is small; the cost of silently mis-analyzing data is much larger.
  2. A unit test exercises a single function with known inputs and checks the output. An integration test runs a pipeline segment (multiple functions) on synthetic data and checks the pipeline’s output. An end-to-end test runs the full analysis pipeline from raw data to final figures/tables. Unit tests are fast and narrow; end-to-end tests are slow and broad; integration tests sit in between, with good cost-benefit for analytic pipelines.
  3. expect_snapshot() captures the output of an expression (printed text, numerical result, plot) to a file on first run. On subsequent runs it compares the new output against the saved snapshot. Use it when the output is complex (a printed object, a plot, a formatted table) and not easily expressed as a literal value. Review changes manually when they occur; do not auto-accept.