--- title: "Comparison with Stata's 'margins' command" date: "`r Sys.Date()`" output: html_document: fig_caption: false toc: true toc_float: collapsed: false smooth_scroll: false toc_depth: 2 vignette: > %\VignetteIndexEntry{Comparison with Stata's 'margins' command} %\VignetteEngine{knitr::rmarkdown} %\VignetteEncoding{UTF-8} --- **margins** is intended as a port of (some of) the features of Stata's `margins` command. This vignette compares output from Stata's `margins` command for linear models against the output of **margins**. ```{r} library("margins") options(width = 100) ``` ## OLS marginal effects ### Stata ``` . quietly reg mpg cyl hp wt . margins, dydx(*) Average marginal effects Number of obs = 32 Model VCE : OLS Expression : Linear prediction, predict() dy/dx w.r.t. : cyl hp wt ------------------------------------------------------------------------------ | Delta-method | dy/dx Std. Err. t P>|t| [95% Conf. Interval] -------------+---------------------------------------------------------------- cyl | -.9416166 .5509165 -1.71 0.098 -2.070118 .1868846 hp | -.0180381 .0118763 -1.52 0.140 -.0423655 .0062893 wt | -3.166973 .740576 -4.28 0.000 -4.683975 -1.649972 ------------------------------------------------------------------------------ ``` ### R ```{r} library("margins") x <- lm(mpg ~ cyl + hp + wt, data = mtcars) summary(margins(x)) ``` --- ## OLS marginal effects with interaction ### Stata ``` . quietly reg mpg cyl c.hp##c.wt . margins, dydx(*) Average marginal effects Number of obs = 32 Model VCE : OLS Expression : Linear prediction, predict() dy/dx w.r.t. : cyl hp wt ------------------------------------------------------------------------------ | Delta-method | dy/dx Std. Err. t P>|t| [95% Conf. Interval] -------------+---------------------------------------------------------------- cyl | -.3652391 .5086204 -0.72 0.479 -1.408842 .6783638 hp | -.0252715 .0105097 -2.40 0.023 -.0468357 -.0037073 wt | -3.837584 .6730996 -5.70 0.000 -5.21867 -2.456498 ------------------------------------------------------------------------------ ``` ### R ```{r} x <- lm(mpg ~ cyl + hp * wt, data = mtcars) summary(margins(x)) ``` --- ## OLS marginal effects with factor term ### Stata ``` . quietly reg mpg i.cyl hp wt . margins, dydx(*) Average marginal effects Number of obs = 32 Model VCE : OLS Expression : Linear prediction, predict() dy/dx w.r.t. : 6.cyl 8.cyl hp wt ------------------------------------------------------------------------------ | Delta-method | dy/dx Std. Err. t P>|t| [95% Conf. Interval] -------------+---------------------------------------------------------------- cyl | 6 | -3.359024 1.40167 -2.40 0.024 -6.235014 -.4830353 8 | -3.185884 2.170476 -1.47 0.154 -7.639332 1.267564 | hp | -.0231198 .0119522 -1.93 0.064 -.0476437 .0014041 wt | -3.181404 .7196011 -4.42 0.000 -4.657904 -1.704905 ------------------------------------------------------------------------------ Note: dy/dx for factor levels is the discrete change from the base level. ``` ### R ```{r} x <- lm(mpg ~ factor(cyl) + hp + wt, data = mtcars) summary(margins(x)) ``` --- ## OLS marginal effects with squared term ### Stata ``` . quietly reg mpg cyl c.hp##c.hp wt . margins, dydx(*) Average marginal effects Number of obs = 32 Model VCE : OLS Expression : Linear prediction, predict() dy/dx w.r.t. : cyl hp wt ------------------------------------------------------------------------------ | Delta-method | dy/dx Std. Err. t P>|t| [95% Conf. Interval] -------------+---------------------------------------------------------------- cyl | -.3696041 .6163571 -0.60 0.554 -1.634264 .8950561 hp | -.0429018 .0178353 -2.41 0.023 -.0794969 -.0063066 wt | -2.873553 .7301251 -3.94 0.001 -4.371646 -1.37546 ------------------------------------------------------------------------------ ``` ### R ```{r} x <- lm(mpg ~ cyl + hp + I(hp^2) + wt, data = mtcars) summary(margins(x)) ``` --- ## OLS marginal effects with squared term (but no first-order term) ### Stata ``` . gen hp2 = hp^2 . quietly reg mpg cyl hp2 wt . margins, dydx(*) Average marginal effects Number of obs = 32 Model VCE : OLS Expression : Linear prediction, predict() dy/dx w.r.t. : cyl hp2 wt ------------------------------------------------------------------------------ | Delta-method | dy/dx Std. Err. t P>|t| [95% Conf. Interval] -------------+---------------------------------------------------------------- cyl | -1.21919 .5030753 -2.42 0.022 -2.249693 -.1886869 hp2 | -.000028 .0000276 -1.01 0.320 -.0000846 .0000286 wt | -3.218637 .7570747 -4.25 0.000 -4.769435 -1.66784 ------------------------------------------------------------------------------ ``` ### R ```{r} x <- lm(mpg ~ cyl + I(hp^2) + wt, data = mtcars) summary(margins(x)) ``` --- ## Logit effects on log-odds and probability scales ### Stata ``` . quietly logit am cyl hp wt . margins, dydx(*) Average marginal effects Number of obs = 32 Model VCE : OIM Expression : Pr(am), predict() dy/dx w.r.t. : cyl hp wt ------------------------------------------------------------------------------ | Delta-method | dy/dx Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- cyl | .0214527 .0469746 0.46 0.648 -.0706157 .1135212 hp | .0014339 .0006182 2.32 0.020 .0002224 .0026455 wt | -.4025475 .1154098 -3.49 0.000 -.6287466 -.1763484 ------------------------------------------------------------------------------ . quietly logit am cyl hp wt . margins, dydx(*) predict(xb) Average marginal effects Number of obs = 32 Model VCE : OIM Expression : Linear prediction (log odds), predict(xb) dy/dx w.r.t. : cyl hp wt ------------------------------------------------------------------------------ | Delta-method | dy/dx Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- cyl | .4875978 1.071621 0.46 0.649 -1.612741 2.587936 hp | .0325917 .0188611 1.73 0.084 -.0043753 .0695587 wt | -9.14947 4.153326 -2.20 0.028 -17.28984 -1.009101 ------------------------------------------------------------------------------ ``` ### R ```{r} x <- glm(am ~ cyl + hp + wt, data = mtcars, family = binomial) # AME summary(margins(x, type = "response")) # AME and MEM equivalent on "link" scale summary(margins(x, type = "link")) ``` --- ## Logit effects with factor variable on log-odds and probability scales ### Stata ``` . quietly logit am i.cyl hp wt . margins, dydx(*) predict(xb) Average marginal effects Number of obs = 32 Model VCE : OIM Expression : Linear prediction (log odds), predict(xb) dy/dx w.r.t. : 6.cyl 8.cyl hp wt ------------------------------------------------------------------------------ | Delta-method | dy/dx Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- cyl | 6 | 2.765754 3.156829 0.88 0.381 -3.421517 8.953025 8 | -8.388958 13.16745 -0.64 0.524 -34.1967 17.41878 | hp | .103209 .0960655 1.07 0.283 -.0850759 .2914939 wt | -10.67598 5.441998 -1.96 0.050 -21.3421 -.0098575 ------------------------------------------------------------------------------ Note: dy/dx for factor levels is the discrete change from the base level. . margins, dydx(*) Average marginal effects Number of obs = 32 Model VCE : OIM Expression : Pr(am), predict() dy/dx w.r.t. : 6.cyl 8.cyl hp wt ------------------------------------------------------------------------------ | Delta-method | dy/dx Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- cyl | 6 | .1197978 .1062873 1.13 0.260 -.0885214 .3281171 8 | -.3478575 .2067542 -1.68 0.092 -.7530883 .0573732 | hp | .0033268 .0029852 1.11 0.265 -.0025241 .0091777 wt | -.3441297 .1188604 -2.90 0.004 -.5770919 -.1111675 ------------------------------------------------------------------------------ Note: dy/dx for factor levels is the discrete change from the base level. ``` ### R ```{r} x <- glm(am ~ factor(cyl) + hp + wt, data = mtcars, family = binomial) # Log-odds summary(margins(x, type = "link")) # Probability with continuous factors summary(margins(x, type = "response")) ``` --- ## Logit with interaction on probability and Log-Odds scales ### Stata ``` . quietly logit am cyl c.hp##c.wt . margins, dydx(*) Average marginal effects Number of obs = 32 Model VCE : OIM Expression : Pr(am), predict() dy/dx w.r.t. : cyl hp wt ------------------------------------------------------------------------------ | Delta-method | dy/dx Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- cyl | .0215633 .0492676 0.44 0.662 -.0749994 .1181261 hp | .0026673 .0023004 1.16 0.246 -.0018414 .007176 wt | -.5157922 .2685806 -1.92 0.055 -1.042201 .0106162 ------------------------------------------------------------------------------ . margins, dydx(*) predict(xb) Average marginal effects Number of obs = 32 Model VCE : OIM Expression : Linear prediction (log odds), predict(xb) dy/dx w.r.t. : cyl hp wt ------------------------------------------------------------------------------ | Delta-method | dy/dx Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- cyl | .5156396 1.169458 0.44 0.659 -1.776456 2.807735 hp | .0515116 .035699 1.44 0.149 -.0184571 .1214804 wt | -12.24264 7.678428 -1.59 0.111 -27.29208 2.806807 ------------------------------------------------------------------------------ ``` ### R ```{r} x <- glm(am ~ cyl + hp * wt, data = mtcars, family = binomial) # AME summary(margins(x, type = "response")) # AME and MEM equivalent on "link" scale summary(margins(x, type = "link")) ``` ## Probit effects on latent and probability scales ### Stata ``` . quietly probit am cyl c.hp##c.wt . margins, dydx(*) predict(xb) Average marginal effects Number of obs = 32 Model VCE : OIM Expression : Linear prediction, predict(xb) dy/dx w.r.t. : cyl hp wt ------------------------------------------------------------------------------ | Delta-method | dy/dx Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- cyl | .2974758 .6629205 0.45 0.654 -1.001825 1.596776 hp | .0277713 .0193121 1.44 0.150 -.0100797 .0656223 wt | -6.626949 4.096208 -1.62 0.106 -14.65537 1.401471 ------------------------------------------------------------------------------ . margins, dydx(*) Average marginal effects Number of obs = 32 Model VCE : OIM Expression : Pr(am), predict() dy/dx w.r.t. : cyl hp wt ------------------------------------------------------------------------------ | Delta-method | dy/dx Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- cyl | .022611 .0498253 0.45 0.650 -.0750447 .1202667 hp | .0025769 .0022607 1.14 0.254 -.001854 .0070077 wt | -.508829 .2625404 -1.94 0.053 -1.023399 .0057408 ------------------------------------------------------------------------------ ``` ### R ```{r} x <- glm(am ~ cyl + hp * wt, data = mtcars, family = binomial(link="probit")) # AME (log-odds) summary(margins(x, type = "link")) # AME (probability) summary(margins(x, type = "response")) ``` ## Poisson effects on latent and probability scales ### Stata ``` . quietly poisson carb cyl c.hp##c.wt . margins, dydx(*) predict(xb) Average marginal effects Number of obs = 32 Model VCE : OIM Expression : Linear prediction, predict(xb) dy/dx w.r.t. : cyl hp wt ------------------------------------------------------------------------------ | Delta-method | dy/dx Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- cyl | -.0993854 .1478936 -0.67 0.502 -.3892516 .1904808 hp | .0066519 .0024217 2.75 0.006 .0019054 .0113984 wt | .1225051 .2035185 0.60 0.547 -.2763837 .521394 ------------------------------------------------------------------------------ . margins, dydx(*) Average marginal effects Number of obs = 32 Model VCE : OIM Expression : Predicted number of events, predict() dy/dx w.r.t. : cyl hp wt ------------------------------------------------------------------------------ | Delta-method | dy/dx Std. Err. z P>|z| [95% Conf. Interval] -------------+---------------------------------------------------------------- cyl | -.2795214 .4169931 -0.67 0.503 -1.096813 .53777 hp | .0175935 .0067179 2.62 0.009 .0044267 .0307604 wt | .2075447 .4859868 0.43 0.669 -.7449719 1.160061 ------------------------------------------------------------------------------ ``` ### R ```{r} x <- glm(carb ~ cyl + hp * wt, data = mtcars, family = poisson) # AME (linear/link) summary(margins(x, type = "link")) # AME (probability) summary(margins(x, type = "response")) ```