# Is delta method better than bootstrap to generate standard error for marginal effects?

I read here, here, here, here, and elsewhere that "

• Parametric bootstrap closely related to objective Bayes. (That’s why it’s a good importance sampling choice.)
• When it applies, parboot approach has both computational and interpretational advantages over MCMC/Gibbs. (...)
• (...) the use of parametric bootstrap sampling to carry out Bayesian inference calculations. (...) possible in a subset of those problems amenable to MCMC analysis, but when feasible the bootstrap approach offers both computational and theoretical advantages."

On the other hand, the delta method is widely used to estimate standard error (SE) for marginal effects (ME) as in the example below taken from here:

set.seed(12345)
n <- 20000
dat <- data.frame(x = rnorm(n))
dat$$y <- ((2 + 2 * dat$$x + + rlogis(n)) > 0) + 0
library(margins)
fit.logit <- glm(y ~ x, binomial, dat)
summary(margins(fit.logit, variables = "x"))

My questions are:

1) Are the SE from the delta method more "robust" than those that would be obtained via the parboot approach?

2) If not, how would the parboot approach be implemented to estimate SE for ME in R?

Any help would be much appreciated.

• Delta method is a first order Taylor series approximation that relies on normality assumptions. Bootstrap is less restrictive in that regard so should be preferable. With your actual example and not this contrived one, fit your model and use something like the syntax below for bootstrapping: summary(margins(fit.logit, vce = "bootstrap", iterations = 5000)). margins offers bootstrap. – Heteroskedastic Jim Nov 14 '18 at 5:41
• Thank you very much @HeteroskedasticJim. I ran the code using the sample above but I get an error. I wanted to know if this is because of this example being contrived? Here is the error: Error in [[<-.data.frame(*tmp*, "z", value = c(Var_dydx_x = NA_real_, : replacement has 5000 rows, data has 1 – Krantz Nov 14 '18 at 6:18
• I noted not to use the contrived example. One, it's 20k cases so takes too long. Two, it fails for some reason with one predictor. – Heteroskedastic Jim Nov 14 '18 at 6:21
• Now used this example x <- lm(mpg ~ cyl + hp * wt, data = mtcars) summary(margins(x, variables = "hp", vce = "bootstrap", iterations = 20000)) from cran.r-project.org/web/packages/margins/vignettes/…. But the error is the same: Error in [[<-.data.frame(*tmp*, "z", value = c(Var_dydx_hp = NA_real_, : replacement has 20000 rows, data has 1. Any thoughts? – Krantz Nov 14 '18 at 6:39
• Do not specify any variable to pull out. Just get everything and interpret the choice effects. The package manual is a good reference for working examples. Though I must admit that the package has some strange behaviours. – Heteroskedastic Jim Nov 14 '18 at 7:09