14  Lab 12: Elasticity — What the Coefficient Means, and What It Costs to Get Wrong

Run this one live: the executable lab.

14.1 What we are doing

A Bayesian regression coefficient is a distribution, and a price elasticity is the one number a category manager will actually act on. Today you will estimate one, and then find three separate ways it can be wrong while every diagnostic still says the model is fine.

Two of the three are not statistical failures at all. One is a routine preprocessing step you have used in every lab since Lecture 6.

Suggested time: ~45 minutes · You will need: pandas, numpy, matplotlib, pymc, arviz, and data/processed/m5_weekly.csv. Four model fits, most under 30 seconds.

14.2 Setup

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import pymc as pm
import arviz as az
import warnings

warnings.filterwarnings("ignore")
RANDOM_STATE = 42

weekly = pd.read_csv("../data/processed/m5_weekly.csv", parse_dates=["week_start"])
weekly["series"] = weekly.store_id + "_" + weekly.cat_id

df = weekly.dropna(subset=["units", "avg_price"]).copy()
df = df[df.avg_price > 0]
df["ly"]    = np.log(df.units.clip(lower=1))
df["lp"]    = np.log(df.avg_price)
df["snap"]  = (df.snap_days >= 5).astype(float)
df["event"] = (df.event_days > 0).astype(float)
df["yr"]    = (df.week_start - df.week_start.min()).dt.days / 365.25

codes, uniq = pd.factorize(df.series)
df["idx"] = codes
J = len(uniq)
y = df.ly.values

print(f"{len(df):,} rows | {J} series | sd(log price) = {df.lp.std():.4f}")

With units and price both in logs, the price coefficient is an elasticity: the percent change in demand per percent change in price. Hold on to that sentence — Step 3 is about the moment it stops being true.

One model function, with switches for the three things we want to vary.

def fit(price_col, controlled=True, standardize_price=False, seed=RANDOM_STATE):
    x = df[price_col].values.astype(float)
    if standardize_price:
        x = (x - x.mean()) / x.std()
    with pm.Model():
        s  = pm.Exponential("sigma", 1.0)
        bp = pm.Normal("b_price", 0, 1)
        bs = pm.Normal("b_snap",  0, 1)
        be = pm.Normal("b_event", 0, 1)
        if controlled:
            a  = pm.Normal("a", y.mean(), 1, shape=J)
            bt = pm.Normal("b_trend", 0, 1)
            mu = (a[df.idx.values] + bp * x + bs * df.snap.values
                  + be * df.event.values + bt * df.yr.values)
        else:
            a  = pm.Normal("a", y.mean(), 1)
            mu = a + bp * x + bs * df.snap.values + be * df.event.values
        pm.Normal("obs", mu=mu, sigma=s, observed=y)
        return pm.sample(1000, tune=1500, chains=4, progressbar=False,
                         random_seed=seed, target_accept=0.9)

def show(idata, controlled=True, label=""):
    vs = ["b_price", "b_snap", "b_event"] + (["b_trend"] if controlled else [])
    su = az.summary(idata, var_names=vs, ci_kind="hdi", ci_prob=0.94, round_to=4)
    print(label)
    print(su[["mean", "hdi94_lb", "hdi94_ub", "r_hat", "ess_bulk"]].to_string())
    return su
8,310 rows | 30 series | sd(log price) = 0.2444

Both units and price are in logs, and that is what makes the price coefficient an elasticity — percent per percent. Say that sentence carefully now, because Step 3 is entirely about the moment it stops being true.

And note 0.2444, the standard deviation of log price across the whole panel. Keep it in view. In Step 3 a different, much smaller number takes its place, and the gap between them is where the damage comes from.

14.3 Step 1 — The naive elasticity

Regress log units on log price, SNAP and events. Nothing else.

idata_naive = fit("lp", controlled=False)
_ = show(idata_naive, controlled=False, label="Naive")

An elasticity near \(-2.8\): a 1% price rise costs almost 3% of volume. Diagnostics are clean, \(\hat{R}\) is 1.00, ESS is in the thousands, and the interval is tight.

It is also badly wrong, and Lecture 12’s DAG says why before any model is fitted. Our thirty series are different categories in different stores. Categories differ systematically in both typical price and typical volume — HOBBIES items cost more per unit and sell fewer of them. So series identity is an arrow into price and into units: a textbook confounder.

Comparing a hobby item to a food item is not measuring a price response. It is measuring the difference between two product categories.

Log units on log price, SNAP, and events. Nothing else.

Naive
           mean  hdi94_lb  hdi94_ub   r_hat   ess_bulk
b_price -2.8046   -2.8531   -2.7522  1.0010  2296.4356

An elasticity of minus 2.8. A one percent price rise costs almost three percent of volume. That is a dramatic, highly elastic result — and it would change how a category is priced.

Look at the diagnostics. R-hat 1.001. ESS 2,296. The interval runs from −2.85 to −2.75 — a width of one tenth. This model is certain. It is also badly wrong, and the useful part is that we knew that before fitting it.

Lecture 12’s DAG says why. Our thirty series are different categories in different stores, and categories differ systematically in both typical price and typical volume — hobby items cost more per unit and sell fewer of them. So series identity is an arrow into price and an arrow into units. That is the textbook shape of a confounder.

Comparing a hobby item to a food item is not measuring a price response. It is measuring the difference between two product categories, and calling it elasticity.

14.4 Step 2 — Control for the confounder

Add per-series intercepts and a trend. Center log price within each series, so the model uses only price variation inside a series over time — which is the variation a manager could actually create.

df["lpc"] = df.lp - df.groupby("series").lp.transform("mean")

idata_ctrl = fit("lpc", controlled=True)
su_ctrl = show(idata_ctrl, label="Controlled")
b_naive = idata_naive.posterior["b_price"].values.ravel()
b_ctrl  = idata_ctrl.posterior["b_price"].values.ravel()

print(f"  naive       {b_naive.mean():+.3f}   width of 94% HDI "
      f"{np.diff(az.hdi(b_naive, prob=0.94))[0]:.3f}")
print(f"  controlled  {b_ctrl.mean():+.3f}   width of 94% HDI "
      f"{np.diff(az.hdi(b_ctrl, prob=0.94))[0]:.3f}")
print(f"  the naive estimate is {b_naive.mean()/b_ctrl.mean():.1f}x too large")
Precision is not accuracy

The naive estimate is nearly six times too large and its interval is half the width. It is not noisy. It is confidently, precisely wrong.

A tight credible interval says the model is sure. It says nothing about whether the model is answering your question — and no diagnostic in ArviZ will ever tell you that it isn’t. Only the DAG does.

# ---- your turn: uncomment and answer in the string ----
# direction = "..."   # <- the naive estimate is too NEGATIVE. Which arrow in the DAG
#                     #    makes the bias run that way rather than the other?
# print(direction)

Add per-series intercepts and a trend, and center log price within each series — so the model uses only price variation inside a series over time, which is the only variation a manager could actually create.

Controlled
           mean  hdi94_lb  hdi94_ub
b_price -0.4921   -0.5949   -0.3872
  naive       -2.805   width of 94% HDI 0.100
  controlled  -0.492   width of 94% HDI 0.207
  the naive estimate is 5.7x too large

Minus 0.49 rather than minus 2.8. The naive estimate was 5.7 times too large.

Now the line that matters, and it is the second column. The naive interval is 0.100 wide. The correct one is 0.207 — twice as wide.

So the wrong model was not noisy. It was more precise than the right one, by a factor of two, while being nearly six times off. It was confidently, tightly, precisely wrong.

That is worth dwelling on, because it inverts an instinct almost everybody has. A narrow credible interval feels like a good sign. All it says is that the model is sure about the answer to the question it was actually asked — and the naive model was asked “how do hobby items differ from food items?”, which the data answers very precisely indeed.

No diagnostic in ArviZ will ever tell you the question was wrong. Only the DAG does, and the DAG is drawn before any code is written.

Try it

Answer the question in the cell: the naive estimate is too negative. Which arrow in the DAG makes the bias run that way rather than the other?

14.5 Step 3 — The preprocessing step that silently destroys the answer

You have standardized predictors in nearly every lab since Lecture 6. It makes Normal(0, 1) priors mean the same thing across variables, and it helps samplers. It is a good habit.

Do it to log price here and watch what happens.

idata_std = fit("lpc", controlled=True, standardize_price=True)
b_std = idata_std.posterior["b_price"].values.ravel()

print(f"  unstandardized  b_price = {b_ctrl.mean():+.4f}")
print(f"  standardized    b_price = {b_std.mean():+.4f}")
print(f"  sd(centered log price)   = {df.lpc.std():.4f}")
print(f"  ratio                   = {b_ctrl.mean()/b_std.mean():.4f}")

The two fits agree perfectly — they are the same model. But the coefficient is no longer the same quantity. Standardizing divided log price by its standard deviation, so the coefficient is now the change in log units per standard deviation of log price, not per percent.

Both diagnostics are clean. Both intervals exclude zero. Nothing anywhere reports a problem. The only thing that changed is what the number means, and the model has no way to tell you.

Now watch it become a business error.

def scenario(b, pct=5.0):
    """Percent change in units, and in revenue, from a `pct` price rise."""
    units = 100 * (np.exp(b * np.log(1 + pct/100)) - 1)
    rev   = 100 * ((1 + pct/100) * np.exp(b * np.log(1 + pct/100)) - 1)
    return units, rev

for label, b in [("correct (elasticity)", b_ctrl), ("standardized", b_std)]:
    u, r = scenario(b)
    print(f"  {label:<22s} 5% price rise -> units {u.mean():+.2f}%  "
          f"revenue {r.mean():+.2f}%")
This is the lab’s point

Read the two lines. Using the standardized coefficient, a manager is told a 5% price rise costs about a tenth of a percent of volume. The truth is more than twenty times that — because the coefficient was divided by the standard deviation of centered log price, which is only 0.0455. The revenue projection is inflated too, from \(+2.5\%\) to \(+4.9\%\).

The smaller the standard deviation you divide by, the larger the distortion. Here the within-series price variation is genuinely small, which makes the error enormous — and entirely invisible.

The model was right. The sampler was right. The diagnostics were right. The error entered through a preprocessing convenience three cells earlier, and it survived all the way to a stocking decision because nothing downstream carries units.

The habit that prevents it: write down what one unit of each variable is before you rescale it, and check that the coefficient still answers the question you asked.

# ---- your turn: uncomment and fill in ----
# print("To recover the elasticity from the standardized fit, multiply it by ___")

You have standardized predictors in nearly every lab since Lecture 6. It makes Normal(0,1) priors mean the same thing across variables, it helps samplers, and it is a genuinely good habit that I am not about to tell you to abandon. Do it to log price here.

  unstandardized  b_price = -0.4921
  standardized    b_price = -0.0224
  sd(centered log price)   = 0.0455
  ratio                   = 21.9384

The two fits agree perfectly — they are the same model, and both have clean diagnostics and intervals excluding zero. But the coefficient went from −0.49 to −0.022, a factor of nearly twenty-two.

Nothing is broken. Standardizing divided log price by its standard deviation, so the coefficient is now the change in log units per standard deviation of log price rather than per percent. Same model, same fit, different unit — and the number no longer means “elasticity.”

Now look at what it was divided by. 0.0455.

Remember the setup printed 0.2444 for log price across the whole panel? This is not that number. This is the standard deviation of log price within a series, after centering — and it is five times smaller, because prices move very little inside a single store-category over time.

And here is the irony worth naming. That within-series variation is small because of the good decision we made in Step 2. Centering within series is what removed the confounding. It is also what makes the standard deviation tiny, and the smaller the number you divide by, the larger the distortion. The correct modeling choice amplified the preprocessing error.

  correct (elasticity)   5% price rise -> units -2.37%  revenue +2.51%
  standardized           5% price rise -> units -0.11%  revenue +4.89%

Read those two lines as advice given to a person.

Using the standardized coefficient, a manager is told that a 5% price rise costs about a tenth of a percent of volume — essentially free. The truth is −2.37%, more than twenty times larger. And the revenue projection nearly doubles, from +2.5% to +4.9%. That is a real decision, made from a real model, that passed every check.

The model was right. The sampler was right. The diagnostics were right. The error entered through a preprocessing convenience three cells earlier, and it survived all the way to a stocking decision because nothing downstream carries units. A float is a float. -0.0224 looks exactly as authoritative as -0.4921.

Try it

Fill in the blank: to recover the elasticity from the standardized fit, what do you multiply by?

And the habit that prevents this, which is the one line to take from the lab: write down what one unit of each variable is before you rescale it, and check that the coefficient still answers the question you asked. Keep standardizing. Just know what you are holding afterwards.

14.6 Step 4 — Ask the posterior a business question

The posterior is samples, so any question is answered by counting. Use the correct fit.

print(f"  P(elasticity < 0)     = {(b_ctrl < 0).mean():.3f}")
print(f"  P(elasticity < -0.5)  = {(b_ctrl < -0.5).mean():.3f}")
print(f"  P(-0.5 < e < 0)       = {((b_ctrl > -0.5) & (b_ctrl < 0)).mean():.3f}")

The first number settles the direction: every posterior sample is negative, so demand falls when price rises. There is no doubt at all.

The second is the interesting one. Elasticity \(-0.5\) is a meaningful business threshold, and the posterior puts us within a few points of a coin flip. The data cannot settle whether demand is elastic or inelastic, and that is the honest report. A significance test would have announced the direction and said nothing about the question the manager actually asked.

fig, ax = plt.subplots(figsize=(7, 3.6))
ax.hist(b_ctrl, bins=60, color="steelblue", alpha=0.8)
ax.axvline(-0.5, color="crimson", ls="--", lw=1.5, label="elastic / inelastic boundary")
ax.axvline(0, color="gray", ls=":", lw=1)
ax.set_xlabel("price elasticity"); ax.set_yticks([])
ax.legend(); plt.tight_layout(); plt.show()

Now the decision the manager actually faces.

u5, r5 = scenario(b_ctrl, 5.0)
print(f"  5% price rise:  units {u5.mean():+.2f}%  "
      f"[{np.percentile(u5,3):+.2f}, {np.percentile(u5,97):+.2f}]")
print(f"                  revenue {r5.mean():+.2f}%  "
      f"[{np.percentile(r5,3):+.2f}, {np.percentile(r5,97):+.2f}]")
print(f"  P(revenue rises) = {(r5 > 0).mean():.3f}")

Revenue rises with probability 1.00 across every posterior sample — even though the elasticity might well be past \(-0.5\). Both facts are true at once, and only the posterior lets you say so.

# ---- your turn: uncomment and answer in the strings ----
# recommend = "..."   # <- two sentences to the manager, with a probability attached
# caveat    = "..."   # <- name one DAG arrow that, if wrong, invalidates this
# print(recommend); print(caveat)

Back to the correct fit. The posterior is samples, so any question is answered by counting.

  P(elasticity < 0)     = 1.000
  P(elasticity < -0.5)  = 0.446
  P(-0.5 < e < 0)       = 0.554

The first line settles the direction completely. Every posterior sample is negative — demand falls when price rises, with no doubt whatsoever.

The second line is the interesting one, and it is the reason to have fitted a Bayesian model at all. Minus one half is a meaningful business boundary: past it, demand is elastic and a price rise loses revenue. The posterior puts us at 0.446 — within a few points of a coin flip.

So the honest report is: the data cannot settle whether demand is elastic or inelastic. Not “it is inelastic.” Not “the effect is significant.” A significance test here would have announced the direction — which was never in question — and said nothing at all about the boundary the manager actually cares about.

  5% price rise:  units -2.37%  [-2.86, -1.87]
                  revenue +2.51%  [+2.00, +3.04]
  P(revenue rises) = 1.000

And now the resolution, which is the most useful thing in the lab.

Revenue rises with probability 1.000 — across every posterior sample — even though the elasticity may well be past −0.5. Both of those are true at once. The threshold question is genuinely unresolved; the decision question is not.

Only a posterior lets you say both. A point estimate forces you to pick one number and then argue about which side of the line it falls on. The distribution lets you answer the question that was actually asked — should we raise the price? — with a probability attached, while being honest that a different question remains open.

Try it

Write the two strings: your two-sentence recommendation to the manager with a probability attached, and one DAG arrow that, if wrong, invalidates it.

14.7 Before you leave

Discuss (the three failures): This lab produced a wrong elasticity three ways — an omitted confounder, a rescaled coefficient, and a threshold question the data cannot answer. Which of the three would any ArviZ diagnostic have caught? What does your answer imply about what diagnostics are for?

Discuss (the habit): Standardizing is genuinely good practice and you should keep doing it. State the rule that lets you keep the habit without repeating Step 3’s error.

Discuss (the decision): You must advise on the 5% rise using only Step 4’s output. Write the recommendation. Then say what you would want measured before a 20% rise, and why extrapolating this elasticity that far is a different kind of claim.

Discuss (the course): You have now estimated the same kind of quantity with OLS, XGBoost, and a Bayesian regression. Only one of them naturally answers “how likely is it that revenue rises?” Explain what about the Bayesian output makes that question cheap, and name one forecasting task from this course where you would still not reach for it.

Solutions posted to Canvas after class.

Try it

Four questions, and the last one closes the course.

The three failures. You produced a wrong elasticity three ways — an omitted confounder, a rescaled coefficient, and a threshold the data cannot resolve. Which would any ArviZ diagnostic have caught?

that the computation did what the code said. Everything about whether the code said the right thing is yours.*

The habit. Standardizing is good practice and you should keep doing it. State the rule that lets you keep the habit without repeating Step 3.

The decision. Advise on the 5% rise using only Step 4’s output. Then say what you would want measured before a 20% rise — and why extrapolating this elasticity that far is a different kind of claim.

that actually occurred within series, and the standard deviation of that is 0.0455 — a few percent. A 20% move is far outside the range the data has ever seen. That is extrapolation, not estimation, and it is the same limit you met in Lab 4 when a tree could not predict beyond its training range.*

The course. You have now estimated this kind of quantity with OLS, with XGBoost, and with a Bayesian regression. Only one of them naturally answers “how likely is it that revenue rises?” Explain what about the Bayesian output makes that question cheap — and name one forecasting task from this course where you would still not reach for it.

And that is the last lab. Fifteen weeks, one dataset, a dozen model families, and a scoreboard that ends with a 46-coefficient linear model at the top.

The through-line was never which method wins. It was that every one of these models will return a confident number whether or not it deserves one — a tree that cannot extrapolate, a network whose result moves 230 RMSE with the seed, a forecast interval that is honest and useless, a hierarchy pooling series that never belonged together, and today a coefficient that stopped being an elasticity three cells before anybody read it.

None of those announced themselves. Every one of them was caught by holding something out and checking. That habit is the course. Solutions go up on Canvas after the deadline. See you in Week 16 for the synthesis.

14.8 What your run should produce

Your numbers will differ where a model is stochastic; these are the values this lab produces on the pinned environment.

random_seed=42; four fits, whole lab ≈ 2 minutes.

Quantity Value
Panel 8,310 rows · 30 series · sd(log price) = 0.2444
Naive b_price −2.8046 · HDI [−2.8531, −2.7522] · R-hat 1.001 · ESS 2,296
Controlled b_price −0.4921 · HDI [−0.5949, −0.3872] · R-hat 1.002 · ESS 2,736
Bias naive is 5.7× too large
Interval widths naive 0.100 vs controlled 0.207 — the wrong model is twice as precise
Standardized b_price −0.0224
sd(centered log price) 0.0455 — five times smaller than the panel-wide 0.2444
Rescaling ratio 21.94 ( = 1 / 0.0455 )
5% rise, correct units −2.37% [−2.86, −1.87] · revenue +2.51% [+2.00, +3.04]
5% rise, standardized units −0.11% · revenue +4.89% — 21× understated
P(elasticity < 0) 1.000
P(elasticity < −0.5) 0.446 — near a coin flip
P(revenue rises at +5%) 1.000