16  Elasticity

Warning

Slopes and elasticities can only be calculated for continuous numeric variables. The slopes() functions will automatically revert to comparisons() for binary or categorical variables.

In some contexts, it is useful to interpret the results of a regression model in terms of elasticity or semi-elasticity. One strategy to achieve that is to estimate a log-log or a semilog model, where the left and/or right-hand side variables are logged. Another approach is to note that \(\frac{\partial ln(x)}{\partial x}=\frac{1}{x}\), and to post-process the marginal effects to transform them into elasticities or semi-elasticities.

For example, say we estimate a linear model of this form:

\[y = \beta_0 + \beta_1 x_1 + \beta_2 x_2 + \varepsilon\]

Let \(\hat{y}\) be the adjusted prediction made by the model for some combination of covariates \(x_1\) and \(x_2\). The slope with respect to \(x_1\) (or “marginal effect”) is:

\[\frac{\partial \hat{y}}{\partial x_1}\]

We can estimate the slope (“dydx”), semi-elasticities (“dyex” and “eydx”), or elasticity (“eyex”) with respect to \(x_1\) as follows:

\[\eta_1=\frac{\partial \hat{y}}{\partial x_1}\] \[\eta_2=\frac{\partial \hat{y}}{\partial x_1}\cdot x_1\] \[\eta_3=\frac{\partial \hat{y}}{\partial x_1}\cdot \frac{1}{\hat{y}}\] \[\eta_4=\frac{\partial \hat{y}}{\partial x_1}\cdot \frac{x_1}{\hat{y}}\]

with interpretations as follows:

  1. “dydx”: Increasing the \(x_1\) variable by 1 unit is associated with a change of \(\eta_1\) units in \(y\).
  2. “dyex”: Increasing the \(x_1\) variable by 1% of its baseline value is associated with a change of \(\eta_2/100\) units in \(y\).
  3. “eydx”: Increasing the \(x_1\) variable by 1 unit is associated with a proportional change of \(\eta_3\) in \(y\), that is, a change of \(100\cdot\eta_3\)% of \(y\)’s baseline value.
  4. “eyex”: Increasing \(x_1\) by 1% of its baseline value is associated with a change of \(\eta_4\)% in the baseline value of \(y\).

These finite-change interpretations are exact when the prediction is linear in the focal predictor, as in the examples below. For nonlinear prediction functions, they are first-order approximations valid for sufficiently small changes in the predictor.

With the marginaleffects package, these quantities are easy to compute with the slopes() function and the slope argument.

The avg_slopes() function averages the observation-specific quantities. In particular, average proportional changes generally differ from proportional changes in the average prediction.

library(marginaleffects)
dat <- transform(mtcars, y = mpg, x = wt)
mod <- lm(y ~ x, data = dat)

pct_change_x <- dat$x / 100
p0 <- predict(mod)
p1 <- predict(mod, newdata = transform(dat, x = x + pct_change_x))
mean(p1 - p0)
#> [1] -0.171945

s <- avg_slopes(mod, slope = "dyex")
s$estimate
#> [1] -17.1945
import numpy as np
import statsmodels.formula.api as smf
from marginaleffects import avg_slopes
from statsmodels.datasets import get_rdataset

mtcars = get_rdataset("mtcars").data

dat = mtcars.assign(
    y=mtcars["mpg"],
    x=mtcars["wt"]
)

mod = smf.ols("y ~ x", data=dat).fit()

pct_change_x = dat["x"] / 100

p0 = mod.predict(dat)
p1 = mod.predict(dat.assign(x=dat["x"] + pct_change_x))

print(np.mean(p1 - p0))
#> -0.17194501167342024

s_dyex = avg_slopes(mod, slope="dyex")
print(float(s_dyex["estimate"][0]))
#> -17.194501167338508

The avg_slopes() estimate is the average semi-elasticity, -17.195. Increasing each observation’s x by 1% of its baseline value is associated with an average change of -0.172 units in predicted y. The avg_slopes() estimate is 100 times the average change computed manually.

yhat <- predict(mod)
p0 <- predict(mod, newdata = dat)
p1 <- predict(mod, newdata = transform(dat, x = x + 1))
mean((p1 - p0) / yhat)
#> [1] -0.29217

s <- avg_slopes(mod, slope = "eydx")
s$estimate
#> [1] -0.29217
yhat = mod.predict(dat)

p0 = mod.predict(dat)
p1 = mod.predict(dat.assign(x=dat["x"] + 1))

print(np.mean((p1 - p0) / yhat))
#> -0.29216997514688453

s_eydx = avg_slopes(mod, slope="eydx")
print(float(s_eydx["estimate"][0]))
#> -0.29216997514682524

Increasing each observation’s x by 1 unit is associated with an average proportional change of -0.292 in predicted y, or -29.2%. Each observation’s change is expressed relative to its own baseline prediction before averaging.

pct_change_x <- 1 / dat$x
pct_change_y <- (p1 - p0) / yhat
mean(pct_change_y / pct_change_x)
#> [1] -1.038292

s <- avg_slopes(mod, slope = "eyex")
s$estimate
#> [1] -1.038292
pct_change_x = 1 / dat["x"]
pct_change_y = (p1 - p0) / yhat

print(np.mean(pct_change_y / pct_change_x))
#> -1.038292137479023

s_eyex = avg_slopes(mod, slope="eyex")
print(float(s_eyex["estimate"][0]))
#> -1.0382921374787804

Increasing each observation’s x by 1% of its baseline value is associated with an average percentage change of -1.038% in predicted y, with each change expressed relative to that observation’s baseline prediction.