Each box summarizes a conditional distribution. The group means are the cell-mean estimates of \mathbf{E}[\log(\textit{spend}) \mid \textit{broadband}].
Reading a boxplot
Element
Meaning
Thick middle line
Median
Box
Interquartile range (IQR): first to third quartile
Whiskers
Most extreme points within 1.5 IQR of the box
Dots beyond the whiskers
Outliers
Households with broadband spend more on average, and the whole distribution shifts up, not only its centre.
Checkpoint: sample or population?
The histogram, the boxplots and the group means were computed from n = 10{,}000 households.
Are these sample objects or population objects?
What population object does each one estimate?
What would change if we drew a fresh sample of 10,000 households?
The linear model and OLS
A line as an approximation to the conditional mean
Linear regression assumes a linear conditional mean
Here \mathbf{E}[Y_i \mid X_i] = 1 + 3X_i by construction.
Increasing X by one unit raises Y by 3 on average. The red line is the conditional mean.
Cell means recover the slope, with noise
table(X)
## X
## 0 1 2 3 4
## 181 193 78 36 12
mean(Y[X ==2]) -mean(Y[X ==1])
## [1] 2.874353
The difference is close to 3 but not equal to it: each cell mean carries sampling error, and the cells with large X are small.
Why cell means do not scale:
With a continuous covariate, no two units share exactly the same x. Every cell has one observation.
With many covariates, the number of cells explodes and most of them are empty.
A linear model pools information across cells: two parameters instead of one mean per cell.
OLS chooses the best-fitting line
The ordinary least squares (OLS) estimator solves
\min_{b} \sum_{i=1}^n (Y_i - X_i'b)^2.
Each unit contributes its squared approximation error; the minimizer is the best fit in total. The first-order conditions give a closed form involving a matrix inverse:
The first level, Dominick’s, will be the baseline category in regressions. The three brands have equal numbers of rows, but how often each is featured differs.
Each brand occupies its own price range
p_price <-ggplot(oj,aes(x = brand, y =log(price), fill = brand)) +geom_boxplot(show.legend =FALSE) +scale_fill_manual(values = smu_colors) +labs(x =NULL, y ="log(price)")
Tropicana is the premium brand, Minute Maid sits in the middle, and Dominick’s is the store brand.
Prices also vary a lot within brand across stores and weeks. That variation identifies the price response.
Logs straighten the demand curve
Raw scale: occasional bursts of sales at low prices, consistent with stocking up during promotions. Log scale: a roughly linear law of demand, sales fall as price rises.
model.matrix() shows the regressors R actually used: an intercept column, log price, and one indicator per non-baseline brand.
glm() with its default Gaussian family gives the same estimates as lm(). We use it here so that the interface is familiar when the family changes next lecture.
The fitted model
summary(OLS_model)
##
## Call:
## glm(formula = log(sales) ~ log(price) + brand, data = oj)
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 10.82882 0.01453 745.04 <2e-16 ***
## log(price) -3.13869 0.02293 -136.89 <2e-16 ***
## brandminute.maid 0.87017 0.01293 67.32 <2e-16 ***
## brandtropicana 1.52994 0.01631 93.81 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for gaussian family taken to be 0.6296804)
##
## Null deviance: 30079 on 28946 degrees of freedom
## Residual deviance: 18225 on 28943 degrees of freedom
## AIC: 68765
##
## Number of Fisher Scoring iterations: 2
Read the coefficients
Term
Estimate
Interpretation
log(price)
-3.14
A 1% higher price goes with about 3.1% lower sales, holding brand fixed
brandminute.maid
0.87
Log sales are 0.87 higher than Dominick’s at the same price
brandtropicana
1.53
Log sales are 1.53 higher than Dominick’s at the same price
The price coefficient in a log-log model is an elasticity.
Every coefficient is a comparison that holds the other regressors fixed. It is not the raw difference between brands in the data.
Interpreting a logged outcome
Model
Change in D
Change in Y
\log Y on \log D
1% higher D
about \beta%, an elasticity
\log Y on continuous D
D up by one unit
about 100\beta%
\log Y on a 0/1 dummy D
D from 0 to 1
(e^{\beta} - 1)\times 100%
The first two rows use the calculus approximation \Delta\log Y \approx \Delta Y / Y, accurate for small changes.
For a dummy, no approximation is needed: \log(Y_1) - \log(Y_0) = \beta gives \dfrac{Y_1 - Y_0}{Y_0} = e^{\beta} - 1.
(exp(coef(OLS_model)["brandtropicana"]) -1) *100
## brandtropicana
## 361.7913
At the same price, Tropicana sells about 362% more units than Dominick’s, not 153% more.
log(price) * brand * feat includes all main effects, all pairwise products, and the three-way product. The 12 coefficients give an intercept and a slope for each of the six brand-by-promotion cells.
Elasticities with and without promotion
Each cell’s slope is the baseline slope plus every interaction that is switched on:
Promotion makes every brand’s demand more elastic, and the ranking across brands changes.
Why does promotion raise elasticity?
More elastic with ads. A feature ad reaches shoppers beyond the brand loyalists, and these marginal customers respond more to price. Promotions also tend to coincide with price cuts, and stocking up amplifies the response.
The ranking changed. Without feat, Minute Maid looked as price sensitive as Dominick’s (-3.32 against -3.38) while Tropicana stood apart (-2.71). With feat, Minute Maid and Tropicana look alike within each promotion state.
Something correlated with brand and with sales was missing from the earlier model.
Two-thirds of Minute Maid’s units are sold in featured weeks, against 40% of Tropicana’s.
Omitting promotion biased the brand elasticities
Minute Maid’s sales came disproportionately from featured weeks, which combine low prices with high volumes.
Without feat, the price coefficient for Minute Maid absorbs the promotion effect, so Minute Maid looks as price sensitive as the store brand.
Omitted-variable bias. A left-out variable that is correlated with an included regressor and affects the outcome shifts the regressor’s coefficient. The coefficient then describes a different comparison from the intended one.
We return to this logic with the housing data, where it appears in a two-regressor model.
Predict for new observations
Given X_{n+1}, the prediction is X_{n+1}'\hat\beta. predict() builds the row of the model matrix for us:
new_dat <-data.frame(price =c(3, 3.5, 2),brand =c("tropicana", "minute.maid", "dominicks"),feat =c(0, 0, 1))exp(predict(OLS_model, newdata = new_dat)) # predicted units, not log units
## 1 2 3
## 7409.805 2361.282 5728.715
new_dat must contain every variable in the model’s formula, with brand labels that match the factor levels.
OLS_model does not use feat, so that column is ignored here; the full model would use it.
The model predicts log sales, so we exponentiate to get units sold.
Checkpoint: read the orange-juice models
In OLS_model, the Tropicana coefficient is 1.53. Is Tropicana’s sales advantage 153%?
In the interaction model, which coefficient answers “do Tropicana buyers respond to price differently from Dominick’s buyers”?
Why did adding feat change the Minute Maid elasticity so much more than the Tropicana elasticity?
Which model should you use to predict next week’s sales during a promotion?
One-bedroom condos in Los Angeles
Confounding, multicollinearity, and interactions with few observations
## LOCATION n
## 1 C42 - Downtown L.A. 35
## 2 Downtown L.A. 80
## 3 Westwood - Century City 18
One-bedroom condos for sale, downloaded from Redfin and filtered to the three locations with the most listings.
Variable
Meaning
PRICE
Listing price ($)
SQUARE_FEET
Interior size
BATHS
Bathrooms: 1, 1.5, 2, or 2.5
LOCATION
Neighbourhood label from the listing service
Price against size and bathrooms
Price rises with size in a roughly linear way, with substantial spread. Bathrooms take four values, so the right panel is really four conditional distributions.
Start with bathrooms alone
fit_baths <-lm(PRICE ~ BATHS, data = redfin)printCoefmat(coef(summary(fit_baths)))
The bathrooms coefficient falls from about $500,000 to about $113,000 and is no longer statistically distinguishable from zero. Its standard error rises.
cor(redfin$BATHS, redfin$SQUARE_FEET)
## [1] 0.6159194
Bathrooms and size are strongly correlated. The first model credited bathrooms with the effect of size.
Two consequences of correlated regressors
Omitted-variable bias, a bias problem. Size affects price and is correlated with bathrooms. Leaving size out makes the bathrooms coefficient absorb part of the size effect.
Multicollinearity, a precision problem. With both included, the data contain little independent variation in bathrooms at a given size, so the bathrooms coefficient is estimated less precisely.
Model
BATHS estimate
BATHS std. error
PRICE ~ BATHS
500,430
83,542
PRICE ~ BATHS + SQUARE_FEET
112,732
91,016
A coefficient is defined relative to the other regressors in the model. Changing the model changes the question it answers.
Signs can flip: Simpson’s paradox
Regress first-year college GPA on SAT score for one high-school graduating class.
Across the whole class, the relationship is negative.
Within each college the students attend, it is positive.
Students with higher scores enter more selective colleges, where grades are harder to earn. College is a confounder that reverses the sign.
Conditioning on a third variable can change a coefficient’s size, its precision, or its sign.
Location as a categorical regressor
fit_location <-lm(PRICE ~ LOCATION, data = redfin)printCoefmat(coef(summary(fit_location)))
LOCATION is a character column; lm() converts it to a factor and creates two indicators. The baseline is the first level alphabetically, “C42 - Downtown L.A.”.
Westwood listings are about $296,000 more expensive on average than the baseline; the two Downtown labels are indistinguishable.
Size and location: parallel slopes
fit_location_size <-lm(PRICE ~ SQUARE_FEET + LOCATION, data = redfin)printCoefmat(coef(summary(fit_location_size)))
grid_la has 300 rows: 100 sizes for each location. predict() applies the location-specific intercept and slope to each row.
The plot uses the same recipe as the orange-juice slides: points from redfin, lines from grid_la.
The Westwood line rests on 18 listings. Interactions double the parameters spent on location.
Checkpoint: which model prices a condo?
A client asks what a 900 square-foot, one-bathroom Westwood condo should list for.
Which of the fitted models would you use, and what does it predict?
In the interaction model, the Westwood indicator is about -748{,}000 dollars. Is Westwood cheaper?
What would you check before trusting the Westwood slope?
Regression estimates conditional means
Define the target.\mathbf{E}[Y \mid X] is the prediction. OLS approximates it with a line and remains the best linear predictor when the truth is curved.
Build the design. Indicators encode categories, logs turn slopes into elasticities, and interactions let slopes differ across groups.
Read coefficients conditionally. A coefficient is a comparison holding the other regressors fixed. Adding or removing a regressor changes the question.
Keep correlation and causation apart. Omitted promotions and omitted size both shifted coefficients. A regression describes an association in the data at hand.
Reading and data sources
Taddy (2019), Business Data Science, McGraw-Hill, chapters 1 and 2: the online-spending and orange-juice examples.
Common in economic data: price dispersion grows with home value, earnings dispersion grows with education, and revenue dispersion grows with firm size.
A simulated example
set.seed(123)n <-500X <-sample(0:3, n, replace =TRUE)eps <-rnorm(n, mean =0, sd =sqrt(0.5+ X))Y <-0.5+3* X + epshet <-data.frame(X, Y)
The error standard deviation is \sqrt{0.5 + X}, so the spread of Y around the line fans out as X grows.
Consequences of heteroskedasticity
Quantity
Under heteroskedasticity
OLS point estimates
Still unbiased and consistent
summary() standard errors
Wrong
t-tests and confidence intervals
Unreliable
Heteroskedasticity does not affect the validity of \hat\beta. It affects inference about \beta.
Use heteroskedasticity-robust standard errors whenever heteroskedasticity is plausible. In economic data, that is most of the time.
Robust standard errors in R
The sandwich package supplies the robust covariance matrix; coeftest() from lmtest reports the table.
Without an intercept, each brand coefficient is that brand’s own intercept. With Tropicana as the baseline, the other two brands’ coefficients become negative shifts.
Where the OLS formula comes from
Differentiate the sum of squared residuals with respect to b and set the result to zero:
The inverse exists when no regressor is a linear combination of the others. That is why a factor with an intercept needs one fewer indicator than levels.
Check the formula by hand
model.matrix() returns the n \times d matrix whose rows are the X_i'; crossprod() computes the sums.
