---
title: "Descriptive Statistics Layers"
output:
  rmarkdown::html_vignette:
    toc: true
vignette: >
  %\VignetteIndexEntry{Descriptive Statistics Layers}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---

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

# Introduction

Descriptive statistics tables are among the most common outputs in clinical trial reporting. Whether you are summarizing demographics in a Table 14.1 or lab parameters across visits, the pattern is the same: compute summary statistics for a continuous variable, then present them in a formatted, publication-ready layout grouped by treatment arm.

In tplyr2, descriptive statistics layers are created with `group_desc()`. The core of your control over the output comes from the `format_strings` parameter within `layer_settings()`. Format strings let you specify exactly which statistics appear, what row label each statistic gets, and how numbers are formatted -- all in one place.

Let's start with a typical example. Using the built-in `tplyr_adsl` dataset, we will summarize age by treatment group.

```{r intro-example}
spec <- tplyr_spec(
  cols = "TRT01P",
  layers = tplyr_layers(
    group_desc("AGE",
      by = "Age (years)",
      settings = layer_settings(
        format_strings = list(
          "n"          = f_str("xx", "n"),
          "Mean (SD)"  = f_str("xx.x (xx.xx)", "mean", "sd"),
          "Median"     = f_str("xx.x", "median"),
          "Q1, Q3"     = f_str("xx.x, xx.x", "q1", "q3"),
          "Min, Max"   = f_str("xx, xx", "min", "max"),
          "Missing"    = f_str("xx", "missing")
        )
      )
    )
  )
)

result <- tplyr_build(spec, tplyr_adsl)
kable(result[, !grepl("^ord", names(result))])
```

A few things to note about this example:

- The `format_strings` parameter is a named list. Each **name** becomes the row label in the output (e.g., "Mean (SD)"), and each **value** is an `f_str()` object that controls the numeric format.
- Inside `f_str()`, the first argument is the format template. The characters `x` define the display width: `xx.x` means two integer digits and one decimal place. The remaining arguments are strings naming the statistics to plug into each format group.
- The `by` argument `"Age (years)"` does not match any column in the data, so tplyr2 treats it as a text label. It appears as an additional `rowlabel` column, which is useful for distinguishing blocks of statistics when multiple layers are combined.

# Built-in Summaries

tplyr2 provides a set of built-in summary statistics that cover the most common needs in clinical reporting. These are computed automatically for every `group_desc()` layer -- you simply reference them by name in your format strings.

| Statistic | Description | Details |
|-----------|-------------|---------|
| `n` | Non-missing count | `length(x[!is.na(x)])` |
| `n_records` | Records assessed (incl. missing analysis values) | `length(x)` |
| `mean` | Arithmetic mean | `mean(x, na.rm = TRUE)` |
| `sd` | Standard deviation | `sd(x, na.rm = TRUE)` |
| `median` | Median | `median(x, na.rm = TRUE)` |
| `var` | Variance | `var(x, na.rm = TRUE)` |
| `min` | Minimum | `min(x)` of finite values |
| `max` | Maximum | `max(x)` of finite values |
| `iqr` | Interquartile range | `IQR(x, type = ...)` |
| `q1` | First quartile (25th percentile) | `quantile(x, 0.25, type = ...)` |
| `q3` | Third quartile (75th percentile) | `quantile(x, 0.75, type = ...)` |
| `missing` | Missing count | `sum(is.na(x))` |
| `total` | Denominator for `pct` | The column population, per `denoms_by` |
| `pct` | Percent of the population with data | `n / total * 100` |

A few important notes about these built-in summaries:

- All statistics use `na.rm = TRUE` by default, so missing values are excluded from calculations (except for `missing` itself, which counts them).
- `total` and `pct` give a desc layer an `n (%)` row -- the number of subjects who contributed data and what share of the arm that is. They respect `denoms_by` and `denom_where` just as a count layer would. See `vignette("format_strings")`.
- `min` and `max` operate on finite values only. If all values in a group are `NA`, the result is `NA_real_` rather than `Inf` or `-Inf`. This avoids formatting issues where infinity symbols would appear in your output.
- The `n` statistic counts non-missing observations, while `missing` counts the `NA` values. Together they sum to the total number of rows in that group.

## Quantile Algorithm

By default, tplyr2 uses R's default quantile algorithm (Type 7) for computing `q1`, `q3`, and `iqr`. This is fine for many applications, but clinical trial reporting often needs to match SAS output, which uses a different algorithm (closest to R's Type 3).

You can change the quantile algorithm globally using `tplyr2_options()`:

```{r quantile-type7}
# Default Type 7 (R default)
spec <- tplyr_spec(
  cols = "TRT01P",
  layers = tplyr_layers(
    group_desc("AGE",
      settings = layer_settings(
        format_strings = list(
          "Q1, Q3" = f_str("xx.x, xx.x", "q1", "q3")
        )
      )
    )
  )
)

result_type7 <- tplyr_build(spec, tplyr_adsl)
kable(result_type7[, !grepl("^ord", names(result_type7))],
      caption = "Type 7 (R default)")
```

```{r quantile-type3}
# Type 3 (matches SAS PROC UNIVARIATE default)
tplyr2_options(quantile_type = 3)

result_type3 <- tplyr_build(spec, tplyr_adsl)
kable(result_type3[, !grepl("^ord", names(result_type3))],
      caption = "Type 3 (SAS-like)")

# Reset to default
tplyr2_options(quantile_type = 7)
```

Notice how the quartile values differ between the two algorithms. The difference is typically small but can matter when you need to produce outputs that match SAS exactly. Type 3 uses the nearest even order statistic and is the closest match to SAS's default behavior.

# Custom Summaries

The built-in summaries cover most standard needs, but clinical reporting sometimes calls for statistics that are not part of the default set -- geometric means, coefficients of variation, trimmed means, and so on. tplyr2 handles this through custom summaries.

## Layer-Level Custom Summaries

You can define custom summaries directly in `layer_settings()` using the `custom_summaries` parameter. Each custom summary is a named expression that uses `.var` as a placeholder for the target variable's values.

Here is an example computing a geometric mean alongside the standard mean:

```{r custom-layer}
spec <- tplyr_spec(
  cols = "TRT01P",
  layers = tplyr_layers(
    group_desc("AGE",
      settings = layer_settings(
        format_strings = list(
          "n"              = f_str("xx", "n"),
          "Mean (SD)"      = f_str("xx.x (xx.xx)", "mean", "sd"),
          "Geometric Mean" = f_str("xx.xx", "geo_mean")
        ),
        custom_summaries = list(
          geo_mean = quote(exp(mean(log(.var[.var > 0]), na.rm = TRUE)))
        )
      )
    )
  )
)

result <- tplyr_build(spec, tplyr_adsl)
kable(result[, !grepl("^ord", names(result))])
```

The key points about custom summaries:

- Use `quote()` to wrap the expression. This delays evaluation until build time, when `.var` is replaced with the actual data vector.
- `.var` refers to the values of the target variable for the current group. It behaves like a numeric vector, so you can apply any R function to it.
- If a custom summary expression throws an error (e.g., trying to take the log of negative values), tplyr2 catches the error and returns `NA_real_` for that group, so your table build will not fail.

## Session-Level Custom Summaries

If you find yourself using the same custom summary across many tables in a study, you can register it at the session level using `tplyr2_options()`. Once registered, the custom statistic is available by name in any `format_strings` specification, just like the built-in summaries.

```{r custom-session}
# Register a coefficient of variation summary for the session
tplyr2_options(
  custom_summaries = list(
    cv = quote(sd(.var, na.rm = TRUE) / mean(.var, na.rm = TRUE) * 100)
  )
)

spec <- tplyr_spec(
  cols = "TRT01P",
  layers = tplyr_layers(
    group_desc("AGE",
      settings = layer_settings(
        format_strings = list(
          "n"          = f_str("xx", "n"),
          "Mean (SD)"  = f_str("xx.x (xx.xx)", "mean", "sd"),
          "CV (%)"     = f_str("xx.x", "cv")
        )
      )
    )
  )
)

result <- tplyr_build(spec, tplyr_adsl)
kable(result[, !grepl("^ord", names(result))])

# Clean up
tplyr2_options(custom_summaries = NULL)
```

## Overriding Built-in Summaries

Custom summaries can even overwrite built-in statistics. If you name a custom summary `"mean"`, it replaces the built-in mean calculation. This is useful when your study requires a non-standard definition of a standard statistic, such as using a trimmed mean instead of the arithmetic mean.

```{r custom-override}
spec <- tplyr_spec(
  cols = "TRT01P",
  layers = tplyr_layers(
    group_desc("AGE",
      settings = layer_settings(
        format_strings = list(
          "Trimmed Mean" = f_str("xx.x", "mean")
        ),
        custom_summaries = list(
          mean = quote(mean(.var, trim = 0.1, na.rm = TRUE))
        )
      )
    )
  )
)

result <- tplyr_build(spec, tplyr_adsl)
kable(result[, !grepl("^ord", names(result))])
```

Layer-level custom summaries always take priority over session-level custom summaries, and both take priority over built-in statistics. This layered precedence gives you fine-grained control: set sensible defaults at the session level, then override on a per-layer basis when needed.

# Multi-Variable Descriptive Statistics

It is common to summarize several continuous variables in a single table -- for example, a demographics table that includes age, height, and weight. Rather than creating separate layers for each variable, you can pass a character vector of variable names to `group_desc()`.

```{r multi-target}
spec <- tplyr_spec(
  cols = "TRT01P",
  layers = tplyr_layers(
    group_desc(c("AGE", "HEIGHTBL", "WEIGHTBL"),
      settings = layer_settings(
        format_strings = list(
          "n"          = f_str("xx", "n"),
          "Mean (SD)"  = f_str("xx.x (xx.xx)", "mean", "sd"),
          "Median"     = f_str("xx.x", "median"),
          "Q1, Q3"     = f_str("xx.x, xx.x", "q1", "q3"),
          "Min, Max"   = f_str("xx, xx", "min", "max"),
          "Missing"    = f_str("xx", "missing")
        )
      )
    )
  )
)

result <- tplyr_build(spec, tplyr_adsl)
kable(result[, !grepl("^ord", names(result))])
```

When multiple target variables are specified:

- Each variable gets its own block of summary rows. The variable name appears as an additional `rowlabel` column (here, `rowlabel1`), with the statistic labels in the next column (`rowlabel2`).
- The same `format_strings` are applied to every variable. Each variable's statistics are computed independently, so differences in scale (e.g., age in years vs. height in centimeters) are handled naturally.
- Ordering is preserved: the first variable's rows appear first, followed by the second, and so on.

You can also combine multi-variable descriptive layers with the `by` parameter. If you add a text label through `by`, it becomes the outermost row label, followed by the variable name, then the statistic label.

```{r multi-target-by}
spec <- tplyr_spec(
  cols = "TRT01P",
  layers = tplyr_layers(
    group_desc(c("AGE", "WEIGHTBL"),
      by = "Demographics",
      settings = layer_settings(
        format_strings = list(
          "n"         = f_str("xx", "n"),
          "Mean (SD)" = f_str("xx.x (xx.xx)", "mean", "sd"),
          "Median"    = f_str("xx.x", "median")
        )
      )
    )
  )
)

result <- tplyr_build(spec, tplyr_adsl)
kable(result[, !grepl("^ord", names(result))])
```

# Records Assessed: `n` Versus `n_records`

Two counts are available and they answer different questions. `n` is the number of *non-missing* values -- the analysis count that pairs naturally with mean/SD. `n_records` is the number of records *assessed* (non-missing plus missing), which some tables report as an "N assessed" line. When there are no missing values the two are identical; they diverge only when the target variable has `NA` values.

```{r n-records}
spec <- tplyr_spec(
  cols = "TRT01P",
  layers = tplyr_layers(
    group_desc("AGE",
      settings = layer_settings(
        format_strings = list(
          "N assessed"   = f_str("xxx", "n_records"),
          "n (analyzed)" = f_str("xxx", "n"),
          "Mean (SD)"    = f_str("xx.x (xx.xx)", "mean", "sd")
        )
      )
    )
  )
)

result <- tplyr_build(spec, tplyr_adsl)
kable(result[, !grepl("^ord", names(result))])
```

# Statistics as Columns

The default descriptive layout puts statistics in rows. Some tables instead want each statistic as its own column. Set `stats_as_columns = TRUE` to transpose the block: without a `by` variable, the treatment groups become the rows and each statistic becomes a column.

```{r stats-as-columns}
spec <- tplyr_spec(
  cols = "TRT01P",
  layers = tplyr_layers(
    group_desc("AGE",
      settings = layer_settings(
        stats_as_columns = TRUE,
        format_strings = list(
          "Mean (SD)" = f_str("xx.x (xx.xx)", "mean", "sd"),
          "Median"    = f_str("xx.x", "median"),
          "Min, Max"  = f_str("xx, xx", "min", "max")
        )
      )
    )
  )
)

result <- tplyr_build(spec, tplyr_adsl)
kable(result[, !grepl("^ord", names(result))])
```

When a `by` variable is present, the by-groups stay as rows and each column becomes a treatment-by-statistic combination (labelled `"<arm> | <stat>"`), so the by dimension is preserved rather than collapsed.

# Adding a Comparison p-value

A demographics table often carries a p-value comparing a continuous characteristic across arms -- an ANOVA, a Kruskal-Wallis test, or a t-test. Attach one with `assoc_test()` in omnibus mode: it runs a function you supply once over the raw source rows and lands the result in a trailing `pval` column. This is the same mechanism count layers use, so the continuous and categorical rows of a demographics table can share one p-value column.

```{r desc-assoc}
spec <- tplyr_spec(
  cols = "TRT01P",
  layers = tplyr_layers(
    group_desc("AGE",
      settings = layer_settings(
        format_strings = list(
          "Mean (SD)" = f_str("xx.x (xx.xx)", "mean", "sd"),
          "Median"    = f_str("xx.x", "median")
        ),
        assoc_test = assoc_test(
          fn = function(.data) anova(lm(AGE ~ TRT01P, .data))[["Pr(>F)"]][1],
          format = f_str("x.xxx", "p"),
          label = "ANOVA p"
        )
      )
    )
  )
)

result <- tplyr_build(spec, tplyr_adsl)
kable(result[, !grepl("^ord", names(result))])
```

The p-value sits on the first statistic row of each `by` group (here there is no `by`, so it appears once, on the first row). The `fn` receives the by-group's raw data frame and returns a scalar; a character return is passed through verbatim. `vignette("binding-statistics")` covers this in full, including how to bind externally computed model results.

# Where to Go From Here

This vignette covered the fundamentals of descriptive statistics layers in tplyr2: built-in summaries, custom summaries, quantile algorithms, and multi-variable analysis. But there is more to explore when it comes to controlling how your numbers look on the page.

- `vignette("format_strings")` -- the format string grammar, the complete statistic keyword reference, rounding, and what happens when a statistic is `NA` (including the `empty` argument of `f_str()`).
- `vignette("precision_alignment")` -- **auto-precision**, which lets the data set the decimal width via `precision_by`, `precision_on`, `precision_data`, and `precision_cap`, and **parenthesis hugging**, which closes the gap between a number and its delimiter.
- `vignette("display_conventions")` -- **`stats_as_columns`**, which transposes the statistics into columns, plus the other display rules a shell may impose.
