---
title: "Choosing a balancing method"
output: rmarkdown::html_vignette
vignette: >
  %\VignetteIndexEntry{Choosing a balancing method}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---

```{r, include = FALSE}
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)
options(balancing.quiet = TRUE)
set.seed(2024)
```

balancing provides six methods. They fall into two families that differ in what they optimize and in the guarantees they offer. This article describes the families, the exposures and estimands each method supports, the difference between exact and tolerance-based balance, the constraint options in `balance_terms()`, and the handling of continuous exposures.

```{r setup}
library(balancing)
```

## Two families of methods

Every method solves a convex optimization problem, but the problems come in two shapes.

The **estimating-equation family** fits weights that solve smooth estimating equations. `bw_entropy()` reweights to the maximum-entropy solution that matches the covariate targets; `bw_ipt()` fits a propensity model through a tilted moment condition; `bw_cbps()` chooses propensity parameters to satisfy covariate-balancing moment conditions. Because these weights solve estimating equations, their uncertainty can be propagated analytically into effect standard errors through `ipw()`. With mean balance and a logit link, all three produce the same ATT weights for a binary exposure.

The **quadratic-program family** minimizes a dispersion or distance objective subject to balance constraints. `bw_energy()` minimizes the energy distance between the reweighted groups; `bw_cfd()` minimizes a kernel distance; `bw_sbw()` minimizes the weight dispersion subject to a balance tolerance. These methods do not carry estimating equations, so inference for them uses resampling rather than the M-estimation route `ipw()` takes.

The family a method belongs to is reflected in its class hierarchy: the estimating-equation methods subclass `estimating_equation_method` and the quadratic-program methods subclass `quadratic_program_method`.

## Which exposures and estimands each method supports

The methods differ in the exposure types and estimands they support. The table below is built from the methods themselves, so it always reflects what the installed version accepts.

```{r support-matrix, echo = FALSE}
methods <- list(
  bw_entropy = bw_entropy(),
  bw_ipt = bw_ipt(),
  bw_cbps = bw_cbps(),
  bw_energy = bw_energy(),
  bw_cfd = bw_cfd(),
  bw_sbw = bw_sbw()
)

exposure_types <- c("binary", "categorical", "continuous")

support_row <- function(name, method) {
  supported <- supported_exposure_types(method)
  cells <- vapply(exposure_types, function(type) {
    if (!type %in% supported) {
      return("-")
    }
    paste(supported_estimands(method, type), collapse = ", ")
  }, character(1))
  c(method = name, cells)
}

support <- as.data.frame(
  do.call(rbind, Map(support_row, names(methods), methods)),
  stringsAsFactors = FALSE
)
rownames(support) <- NULL
names(support) <- c("method", exposure_types)

knitr::kable(support)
```

A dash marks an unsupported exposure type. The estimand vocabulary follows propensity: `"ate"` (average treatment effect), `"att"` (average treatment effect on the treated), `"atu"` (on the untreated, accepted as `"atc"`), and `"ato"` (the overlap estimand). The overlap estimand is available only for `bw_cbps()` with a binary exposure. Continuous exposures permit only `"ate"`.

## Exact versus tolerance-based balance

The families also differ in the balance they guarantee.

The estimating-equation methods achieve **exact** balance on the requested moments. Their weighted covariate means equal the target means to solver precision, which is why the printed summary reports a largest imbalance at the level of numerical noise.

```{r exact}
n <- 600
age <- rnorm(n)
score <- rnorm(n)
exposure <- rbinom(n, 1, plogis(0.5 * age - 0.5 * score))
study <- data.frame(exposure, age, score)

exact <- balance(study, exposure, c(age, score), method = bw_entropy())
exact
```

Exact balance is not always achievable or desirable. When the groups overlap poorly, forcing exact balance can drive the weights to extremes and collapse the effective sample size. A **tolerance** relaxes the constraints: each weighted statistic is held within a band rather than pinned to a point. For entropy balancing, a positive tolerance selects an inexact problem that trades a little residual imbalance for less variable weights. The printed summary records the tolerance the fit used and the largest imbalance it leaves.

```{r tolerance}
relaxed <- balance(
  study,
  exposure,
  c(age, score),
  method = bw_entropy(),
  estimand = "att",
  constraints = balance_terms(tolerance = 0.05)
)
relaxed
```

For `bw_sbw()` the tolerance is not optional. Stable balancing weights minimize weight dispersion subject to a balance band, so the tolerance is the method's central tuning parameter and a fit without a positive tolerance is refused.

```{r sbw}
sbw_fit <- balance(
  study,
  exposure,
  c(age, score),
  method = bw_sbw(),
  constraints = balance_terms(tolerance = 0.02)
)
sbw_fit
```

For `bw_energy()` and `bw_cfd()`, the objective itself drives balance, so the tolerance relaxes any moment constraints you add rather than selecting an inexact solver. A tolerance supplied without moment constraints has nothing to relax and is ignored with a warning.

## Constraints in depth

`balance_terms()` records which covariate functions a method should balance. It has four ingredients.

**Moments** set the highest power of each numeric covariate to balance. A scalar applies to every covariate; a named integer vector sets powers per covariate.

```{r moments}
# Balance means and variances of both covariates
balance_terms(moments = 2)

# Balance the mean of age but the mean and variance of score
balance_terms(moments = c(age = 1, score = 2))
```

Powers above one are ignored for binary indicator columns, since squaring a zero/one column reproduces it.

**Interactions** add all pairwise products of the base columns when set to `TRUE`, excluding products of two indicators of the same factor. This balances the joint distribution of pairs of covariates, not just their marginals, and adds the product terms to the constraint count the printed summary reports. The word names a constraint here rather than a causal quantity: `interactions = TRUE` expands the set of covariate functions the weights must balance, equating them across the exposure groups of a discrete exposure and driving them to zero correlation with a continuous one, and says nothing about two treatments acting together, which is an effect `ipw()` reports for a joint exposure rather than a constraint requested here.

```{r interactions}
with_interactions <- balance(
  study,
  exposure,
  c(age, score),
  method = bw_entropy(),
  constraints = balance_terms(moments = 1, interactions = TRUE)
)
with_interactions
```

**Quantiles** add an indicator column at each requested probability, so that mean balance on the indicator is quantile balance on the covariate. This targets the shape of the distribution rather than its moments, and each added indicator becomes another constraint term. Quantile constraints apply to discrete exposures only.

```{r quantiles}
with_quantiles <- balance(
  study,
  exposure,
  c(age, score),
  method = bw_entropy(),
  constraints = balance_terms(
    moments = 1,
    quantiles = c(0.25, 0.5, 0.75)
  )
)
with_quantiles
```

**Tolerance** sets the per-constraint balance band, as described above. It can be a scalar or a named vector giving a different tolerance per source covariate, with derived columns inheriting their source's tolerance.

```{r per-cov-tolerance}
balance_terms(tolerance = c(age = 0.01, score = 0.1))
```

When `constraints` is left at its default, the estimating-equation methods and stable balancing weights balance first moments, while energy and characteristic function distance balancing add no moment constraints and let their objective drive balance.

## Continuous exposures

For a continuous exposure there are no groups to balance against one another. Instead, the target is independence between the exposure and the covariates: the weighted exposure-covariate correlations should be zero. `bw_entropy()`, `bw_cbps()`, `bw_energy()`, and `bw_sbw()` support continuous exposures; the estimand is always `"ate"`.

```{r continuous}
dose <- 0.7 * age - 0.5 * score + rnorm(n)
dose_study <- data.frame(dose, age, score)

dose_fit <- balance(
  dose_study,
  dose,
  c(age, score),
  method = bw_entropy()
)
dose_fit
```

For a continuous exposure the printed summary reports the largest weighted exposure-covariate correlation rather than a standardized mean difference.

`bw_energy()` is the exception to that target. Its objective minimizes the weighted distance covariance, which asks for full distributional independence rather than for zero correlations, and it leaves a residual weighted correlation of roughly 0.1 to 0.3 at ordinary sample sizes. Pass `constraints = balance_terms(moments = 1)` to add the correlation constraints outright, which costs effective sample size; `?bw_energy` describes the trade and the part `weight_penalty` plays in it.

For `bw_entropy()` and `bw_energy()`, the `distribution_moments` argument on the method constructor sets how many moments of the exposure and covariate marginals are held equal to the sample, extending the identifying conditions beyond the correlation constraints. Both measure that sample under the base measure: any sampling weights, times the `base_weights` of `bw_entropy()`, which `bw_energy()` does not carry. With neither in play the marginals are held equal to the unweighted sample. For `bw_entropy()` it is raised automatically to at least the moments the constraints require; for `bw_energy()` it stands alone, since there `moments` asks for correlation constraints rather than marginal ones.

```{r distribution-moments}
dose_moments <- balance(
  dose_study,
  dose,
  c(age, score),
  method = bw_entropy(distribution_moments = 2, convergence_tolerance = 1e-8)
)
dose_moments
```

The marginal-moment conditions are enforced inside the solver rather than added to the balance table, so they leave the constraint count in the printed summary unchanged. They also make the problem stiffer, so this fit relaxes the `convergence_tolerance` to a looser value to reach the requested balance.

## A note on solvers

Each method has a default numerical backend: the estimating-equation family uses a damped Newton solver, and the quadratic-program family uses an operator-splitting solver, with stable balancing weights falling back automatically to an interior-point solver when the default certifies a problem infeasible. The backend that produced a fit is recorded on the result. These defaults are set to be reliable across a wide range of problems, so a typical analysis does not need to change them.