Optimal treatment regime by exhaustive tree search

Laber & Zhao's argument is that a regime is only useful if a clinician can read it, and that the usual two-step approach – regress the outcome, then take the argmax – optimizes the wrong thing: a small error in the outcome model where the two arms are close flips the recommendation. Their tree instead maximizes the ESTIMATED VALUE directly, \(V(d) = sum Y_i 1{A_i = d(W_i)} / pi_i / sum 1{A_i = d(W_i)} / pi_i\), over axis-aligned splits. This implementation enumerates every split of every covariate at every observed value and recurses.

Usage

Opttre(y, A, W, pi = NULL, max_depth = 2L, min_leaf = 1L)

Arguments

y

Outcome, larger is better.

A

Observed binary treatment, 0/1.

W

Covariates.

pi

Propensity P(A = A_i | W_i), or NULL for the marginal randomization probability.

max_depth

Maximum tree depth; 0 gives one constant rule.

min_leaf

Minimum observations in a leaf.

Value

List with estimate, value, value_all_treated, value_all_control, rule, split_var, split_point, n_leaves, depth, n.

Details

The candidate scan runs COVARIATE-MAJOR then observation, and an improvement must be STRICT, so ties keep the first candidate in that fixed order – R scans a matrix column-major and Python row-major, and an unpinned scan would return a different-but-equally-optimal split per language.

References

Laber, E. B. & Zhao, Y.-Q. (2015). Tree-based methods for individualized treatment regimes. Biometrika, 102(3), 501-514. doi:10.1093/biomet/asv028

Examples

set.seed(6)
n <- 400
W <- matrix(rnorm(n * 2), n, 2)
A <- rbinom(n, 1, 0.5)
y <- W[, 1] + A * ifelse(W[, 1] > 0, 1, -1) + rnorm(n)
t <- Opttre(y, A, W, max_depth = 1L)
c(t$split_var, round(t$split_point, 2))   # splits on W1 near 0
#> [1]  0.00 -0.24