Synth: predictors.op is applied to the treated unit but not to the controls

A bug report for the Synth R package (Hainmueller & Diamond), affecting 1.1-10 (current CRAN) and 1.2-0 (current development). R/dataprep.R is byte-identical between the two.

This has been submitted upstream: j-hai/Synth#2

This repository just holds a standalone reproduction and patch. It is not a fork and not a package.

Correction. An earlier version of this README said Synth "has no public issue tracker". That was wrong — j-hai/Synth exists and is active. I had checked the URL field in Synth's DESCRIPTION, which lists only the author's Stanford homepage, instead of trying the maintainer's GitHub account directly. Apologies for the error.

The problem

dataprep() aggregates predictors over time.predictors.prior separately for the treated unit and the controls. The treated unit is aggregated with the operator the user supplied:

# R/dataprep.R:216
X1 <- apply(X1, 2, paste(predictors.op), na.rm = TRUE)

The controls are aggregated with a hard-coded mean:

# R/dataprep.R:257
X0 <- sapply(X0, apply, 2, mean, na.rm = TRUE, simplify = TRUE)

So for any predictors.op other than "mean", X1 and X0 are computed under different aggregation rules, and the predictor-balance matrix that synth() optimises over is internally inconsistent.

Two things suggest this is an oversight rather than a deliberate asymmetry:

1. The message printed in the control branch a few lines above refers to the argument by name — "We ignore (na.rm = TRUE) all missing values for predictors.op" (R/dataprep.R:252) — in the same block that hard-codes mean.
2. special.predictors operators are applied to both the treated unit and the controls (R/spec.pred.func.R), so the package is inconsistent with itself.

The documentation (?dataprep) describes predictors.op as "the method to be used when aggregating the predictor variables", with no treated/control distinction.

Reproduction

repro.R in this repository. Output on Synth 1.1-10, using the bundled basque data with predictors = "invest", predictors.op = "median", and time.predictors.prior = 1964:1969:

TREATED unit (regionno 17)
  dataprep X1        : 24.41825
  hand median        : 24.41825   <- matches, operator honoured
  hand mean          : 24.64738

CONTROL units
  max |X0 - hand median| : 2.061043   <- should be 0 if predictors.op were honoured
  max |X0 - hand mean|   : 0.000000   <- controls received the MEAN

Every control's predictor is off by up to 2.06 investment-share points relative to the operator that was requested. Asking for predictors.op = "median" changes only the treated unit's row of the balance matrix; the controls' rows are byte-identical to the "mean" call.

Patch

predictors-op.patch — one line, mirroring what line 216 already does for the treated unit:

     X0 <- split(X0, X0[,dim(X0)[2]])
-    X0 <- sapply(X0, apply, 2, mean, na.rm = TRUE, simplify = TRUE)
+    X0 <- sapply(X0, apply, 2, paste(predictors.op), na.rm = TRUE, simplify = TRUE)
     X0 <- as.matrix(X0[-dim(X0)[1],])

With the patch applied, max |X0 - hand median| becomes 0.000000.

Synth's own test suite is unaffected: FAIL 0 | PASS 46 before and after, run via the package's own harness (cd tests && Rscript testthat.R).

What I have not established

Whether the asymmetry is intended. It is possible the authors deliberately want control units summarised by their mean regardless of predictors.op; nothing in the documentation says so, and the two points above argue against it, but only the authors can settle that. If it is intended, the fix is a documentation change rather than a code change.

Licence

Synth is GPL (>= 2). The patch here is offered under the same terms.