Comparison with Stata’s ‘margins’ command

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.

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

library("margins")
x <- lm(mpg ~ cyl + hp + wt, data = mtcars)
summary(margins(x))

##  factor     AME     SE       z      p   lower   upper
##     cyl -0.9416 0.5509 -1.7092 0.0874 -2.0214  0.1382
##      hp -0.0180 0.0119 -1.5188 0.1288 -0.0413  0.0052
##      wt -3.1670 0.7406 -4.2764 0.0000 -4.6185 -1.7155

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

x <- lm(mpg ~ cyl + hp * wt, data = mtcars)
summary(margins(x))

##  factor     AME     SE       z      p   lower   upper
##     cyl -0.3652 0.5086 -0.7181 0.4727 -1.3621  0.6316
##      hp -0.0253 0.0105 -2.4046 0.0162 -0.0459 -0.0047
##      wt -3.8376 0.6731 -5.7014 0.0000 -5.1568 -2.5183

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

x <- lm(mpg ~ factor(cyl) + hp + wt, data = mtcars)
summary(margins(x))

##  factor     AME     SE       z      p   lower   upper
##    cyl6 -3.3590 1.4017 -2.3964 0.0166 -6.1062 -0.6118
##    cyl8 -3.1859 2.1705 -1.4678 0.1422 -7.4399  1.0682
##      hp -0.0231 0.0120 -1.9344 0.0531 -0.0465  0.0003
##      wt -3.1814 0.7196 -4.4211 0.0000 -4.5918 -1.7710

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

x <- lm(mpg ~ cyl + hp + I(hp^2) + wt, data = mtcars)
summary(margins(x))

##  factor     AME     SE       z      p   lower   upper
##     cyl -0.3696 0.6164 -0.5997 0.5487 -1.5776  0.8384
##      hp -0.0429 0.0178 -2.4054 0.0162 -0.0779 -0.0079
##      wt -2.8736 0.7301 -3.9357 0.0001 -4.3046 -1.4425

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

x <- lm(mpg ~ cyl + I(hp^2) + wt, data = mtcars)
summary(margins(x))

##  factor     AME     SE       z      p   lower   upper
##     cyl -1.2192 0.5031 -2.4235 0.0154 -2.2052 -0.2332
##      hp -0.0082 0.0081 -1.0124 0.3114 -0.0241  0.0077
##      wt -3.2186 0.7571 -4.2514 0.0000 -4.7025 -1.7348

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

x <- glm(am ~ cyl + hp + wt, data = mtcars, family = binomial)
# AME
summary(margins(x, type = "response"))

##  factor     AME     SE       z      p   lower   upper
##     cyl  0.0215 0.0470  0.4567 0.6479 -0.0706  0.1135
##      hp  0.0014 0.0006  2.3197 0.0204  0.0002  0.0026
##      wt -0.4025 0.1154 -3.4880 0.0005 -0.6287 -0.1764

# AME and MEM equivalent on "link" scale
summary(margins(x, type = "link"))

##  factor     AME     SE       z      p    lower   upper
##     cyl  0.4876 1.0716  0.4550 0.6491  -1.6127  2.5879
##      hp  0.0326 0.0189  1.7280 0.0840  -0.0044  0.0696
##      wt -9.1495 4.1533 -2.2029 0.0276 -17.2898 -1.0091

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

x <- glm(am ~ factor(cyl) + hp + wt, data = mtcars, family = binomial)
# Log-odds
summary(margins(x, type = "link"))

##  factor      AME      SE       z      p    lower   upper
##    cyl6   2.7658  3.1568  0.8761 0.3810  -3.4214  8.9529
##    cyl8  -8.3890 13.1671 -0.6371 0.5240 -34.1960 17.4181
##      hp   0.1032  0.0961  1.0744 0.2826  -0.0851  0.2915
##      wt -10.6760  5.4418 -1.9618 0.0498 -21.3417 -0.0102

# Probability with continuous factors
summary(margins(x, type = "response"))

##  factor     AME     SE       z      p   lower   upper
##    cyl6  0.1198 0.1063  1.1271 0.2597 -0.0885  0.3281
##    cyl8 -0.3479 0.2068 -1.6825 0.0925 -0.7531  0.0574
##      hp  0.0033 0.0030  1.1144 0.2651 -0.0025  0.0092
##      wt -0.3441 0.1189 -2.8953 0.0038 -0.5771 -0.1112

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

x <- glm(am ~ cyl + hp * wt, data = mtcars, family = binomial)

## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred

# AME
summary(margins(x, type = "response"))

##  factor     AME     SE       z      p   lower  upper
##     cyl  0.0216 0.0493  0.4377 0.6616 -0.0750 0.1181
##      hp  0.0027 0.0023  1.1596 0.2462 -0.0018 0.0072
##      wt -0.5158 0.2685 -1.9209 0.0547 -1.0421 0.0105

# AME and MEM equivalent on "link" scale
summary(margins(x, type = "link"))

##  factor      AME     SE       z      p    lower  upper
##     cyl   0.5156 1.1695  0.4409 0.6593  -1.7765 2.8077
##      hp   0.0515 0.0357  1.4429 0.1490  -0.0185 0.1215
##      wt -12.2426 7.6784 -1.5944 0.1108 -27.2920 2.8067

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

x <- glm(am ~ cyl + hp * wt, data = mtcars, family = binomial(link="probit"))

## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred

# AME (log-odds)
summary(margins(x, type = "link"))

##  factor     AME     SE       z      p    lower  upper
##     cyl  0.2975 0.6472  0.4596 0.6458  -0.9711 1.5660
##      hp  0.0278 0.0184  1.5075 0.1317  -0.0083 0.0639
##      wt -6.6269 3.9095 -1.6951 0.0901 -14.2894 1.0355

# AME (probability)
summary(margins(x, type = "response"))

##  factor     AME     SE       z      p   lower   upper
##     cyl  0.0226 0.0492  0.4598 0.6456 -0.0738  0.1190
##      hp  0.0026 0.0021  1.2244 0.2208 -0.0015  0.0067
##      wt -0.5088 0.2479 -2.0523 0.0401 -0.9948 -0.0229

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

x <- glm(carb ~ cyl + hp * wt, data = mtcars, family = poisson)
# AME (linear/link)
summary(margins(x, type = "link"))

##  factor     AME     SE       z      p   lower  upper
##     cyl -0.0994 0.1479 -0.6720 0.5016 -0.3893 0.1905
##      hp  0.0067 0.0024  2.7468 0.0060  0.0019 0.0114
##      wt  0.1225 0.2035  0.6019 0.5472 -0.2764 0.5214

# AME (probability)
summary(margins(x, type = "response"))

##  factor     AME     SE       z      p   lower  upper
##     cyl -0.2795 0.4170 -0.6703 0.5026 -1.0968 0.5378
##      hp  0.0176 0.0067  2.6189 0.0088  0.0044 0.0308
##      wt  0.2075 0.4860  0.4271 0.6693 -0.7449 1.1600