jkcshea / ivmte

An R package for implementing the method in Mogstad, Santos, and Torgovitsky (2018, Econometrica).
GNU General Public License v3.0
18 stars 2 forks source link

Custom weight functions must have arguments (else produces uninformative error) #189

Closed johnnybonney closed 4 years ago

johnnybonney commented 4 years ago

The package requires custom weighting functions to have arguments. I personally think this is an unnecessary requirement. However, if the requirement is maintained, then it would be helpful to have an informative error message. Here is an example that illustrates the behavior and current error message.

library(data.table)
library(ivmte)

set.seed(1001)
N <- 10000

# setup
dt <- data.table(
  id = 1:N,
  u = runif(N),
  z = sample(1:4, N, replace = T),
  x = sample(0:1, N, replace = T),
  epsilon = rnorm(N, sd = .1)
)

dt[, p := (z == 1)*0.12 + (z == 2)*0.29 + (z == 3)*0.48 + (z == 4)*0.78]
dt[, d := as.integer(u <= p)]
dt[, y0 := 0.9 - 1.1*u + 0.3*u^2]
dt[, y1 := 0.35 - 0.3*u - 0.05*u^2]
dt[, y := y1*d + y0*(1 - d) + epsilon]

# custom weights to reproduce ATE
weight1 <- function() 1
weight0 <- function() -1
knot1 <- function() 0
knot2 <- function() 1

args_manual <- list(
  data = dt,
  propensity = d ~ z,
  ivlike = c(y ~ d | factor(z),
             y ~ d),
  components = l(),
  m1 = ~uSplines(knots = c(0.2, 0.6, 0.8), degree = 0),
  m0 = ~uSplines(knots = c(0.2, 0.6, 0.8), degree = 0),
  target.knots0 = c(knot1, knot2),
  target.knots1 = c(knot1, knot2),
  target.weight0 = c(0, weight0, 0),
  target.weight1 = c(0, weight1, 0),
  lpsolver = "gurobi"
)

(test <- do.call(ivmte, args_manual))
Obtaining propensity scores...

Generating target moments...
Integrating non-spline terms for control group...
Error in if (wKnotVarsU[1] == "...") { : argument is of length zero

It took me a while to figure out what was causing the issue, which can be resolved by adding arguments to the weight/knot functions:

# redefine functions to have (unnecessary) argument
weight1 <- function(x) 1
weight0 <- function(x) -1
knot1 <- function(x) 0
knot2 <- function(x) 1

args_manual$target.knots0 = c(knot1, knot2)
args_manual$target.knots1 = c(knot1, knot2)
args_manual$target.weight0 = c(0, weight0, 0)
args_manual$target.weight1 = c(0, weight1, 0)

(test2 <- do.call(ivmte, args_manual))
Bounds on the target parameter: [-0.3788724, 0.1789346]
Audit terminated successfully after 1 round
jkcshea commented 4 years ago

Ah, but you don't need to declare functions! But this is mentioned only in the reference manual and not the vignette example. I'll update the vignette, as I realize now users would've benefited from this information, as well as some more details on how these weight functions work. (Edit: Sorry, the vignette does mention that scalars can be passed, but I'll add in the details regarding what the arguments of the functions should be.)

Here is what that information is. Each argument of the function is supposed to correspond to a variable name (but not the unobserved term u). In the vignette, the function argument is x because the variable that determines the weight is x. If you wanted the weights to depend on variables x1, x2, x3, then you need to include those arguments in the weight functions. Note that these variables need not be part of the model, as demonstrated by your example. But they do need to be in the data object passed to ivmte.

If you want the knots and weights to be constant, you can simply declare them to be constants.

> ## Declare the knots and weights as constants
> args_manual$target.knots0 = c(0, 1)
> args_manual$target.knots1 = c(0, 1)
> args_manual$target.weight0 = c(0, -1, 0)
> args_manual$target.weight1 = c(0, 1, 0)
> (test3 <- do.call(ivmte, args_manual))

Bounds on the target parameter: [-0.3788724, 0.1789346]
Audit terminated successfully after 1 round
johnnybonney commented 4 years ago

I see. That makes sense.

Personally, I still don't like that ivmte can't handle constant functions, and I really don't like that I get the confusing error message Error in if (wKnotVarsU[1] == "...") { : argument is of length zero. My recommendation would be at the very least to add a check that the knots/weighting functions have arguments (using length(formals(FUNC_HERE))) and return an informative error if the number of arguments is 0.

Moreover, this line in the vignette seems to suggest that scalars are converted to constant functions: (A scalar is interpreted as a constant function.) But a constant function itself is not OK? Seems odd to me.

jkcshea commented 4 years ago

You're absolutely right, sorry about that, I should've made those changes. Let me do that.

jkcshea commented 4 years ago

Done! Thanks for the suggestion. The first example in your post should run just fine now, i.e. weight and knot functions without arguments are interpreted as constant functions.

johnnybonney commented 4 years ago

Great! Thanks!