8  Causal Panel Designs

8.1 Fixed Effects

Longitudinal datasets have become increasingly common in the 20th and 21st centuries and come in a variety of forms. Some of these datasets are purpose-built for tracking individuals or other types of units over time, like the National Longitudinal Survey of Youth (NSLY) or the Panel Income of Survey Dynamics, and are often called panel data. Others are harmonized from existing administrative datasets corresponding to geographic aggregates like the FBI’s Uniform Crime Reports. Still others are representative surveys repeated periodically, like the Censuses, which are called repeated cross-sectional datasets. These do not follow the same people over time, but rather provide an ongoing snapshot of a representative population over multiple periods. These datasets, with their rich structure, have opened the door to new analytical approaches and have been used in a variety of insightful ways.

In this chapter, I’ll focus on two models for analyzing panel data: the fixed effects estimator and the pooled OLS (POLS) estimator, which disregards the panel structure of the data. My treatment of these models will take a traditional approach by assuming that treatment effects are homogeneous—constant across individuals and time. I acknowledge this assumption may feel at odds with the potential outcomes framework we discussed earlier, but panel data methods historically emerged from outside the potential outcomes tradition. But that being said, these traditions have tended to converge, with recent work being focused intently on the role that heterogeneous treatment effects in panel data have on identification and estimation. However, for now and in this chapter, I’ll stick to the classic setup and relax the common treatment effect assumption in a later chapter.

Matrix Notation

When we use the term “panel data,” what do we mean? We mean a dataset where we observe the same units (individuals, firms, countries, schools, etc.) over more than one time period. Often our outcome variable depends on several factors, some of which are observed and some of which are unobserved in our data, and insofar as the unobserved variables are correlated with the treatment variable, then the treatment variable is endogenous and correlations are not estimates of a causal effect. This chapter focuses on the conditions under which a correlation between \(D\) and \(Y\) reflects a causal effect even with unobserved variables that are correlated with the treatment variable. Specifically, if these omitted variables are constant over time, then even if they are heterogeneous across units, we can use panel data estimators to consistently estimate the effect of our treatment variable on outcomes.

There are several different kinds of estimators for panel data, but we will in this chapter only cover two: pooled ordinary least squares (POLS) and fixed effects (FE).1

First we need to set up our notation. With some exceptions, panel methods are usually based on the traditional notation and not the potential outcomes notation.

Let \(Y\) and \(D\equiv(D_1, D_2, \dots, D_k)\) be observable random variables and \(u\) be an unobservable random variable. We are interested in the partial effects of variable \(D_j\) in the population regression function: \[ E\big[Y\mid D_1,D_2, \dots , D_k, u\big] \tag{8.1}\] We observe a sample of \(i=1,2,\dots,N\) cross-sectional units for \(t=1,2,\dots,T\) time periods (a balanced panel). For each unit \(i\), we denote the observable variables for all time periods as \(\{(Y_{it},D_{it}):t=1,2,\dots,T\}\).2 Let \(D_{it}\equiv(D_{it1},D_{it2}, \dots, D_{itk})\) is a \(1\times K\) vector. We typically assume that the actual cross-sectional units (e.g., individuals in a panel) are identical and independent draws from the population in which case \(\{Y_i, D_i, u_i\}^N_{i=1}\sim i.i.d.\), or cross-sectional independence. We describe the main observables, then, as \(Y_i \equiv (Y_{i1}, Y_{i2}, \dots, Y_{iT})'\) and \(D_i \equiv (D_{i1},D_{i2}, \dots, D_{iT})\).

It’s helpful now to illustrate the actual stacking of individual units across their time periods. A single unit \(i\) will have multiple time periods \(t\) \[ Y_i=\left( \begin{array}{c} Y_{i1} \\ \vdots \\ Y_{it} \\ \vdots \\ Y_{iT} \end{array} \right)_{T\times 1} \qquad D_i = \left( \begin{array}{ccccc} D_{i,1,1} & D_{i,1,2} & D_{i,1,j} & \dots & D_{i,1,K} \\ \vdots & \vdots & \vdots & & \vdots \\ D_{i,t,1} & D_{i,t,2} & D_{i,t,j} & \dots & D_{i,t,K} \\ \vdots & \vdots & \vdots & & \vdots \\ D_{i,T,1} & D_{i,T,2} & D_{i,T,j} & \dots & D_{i,T,K} \end{array} \right)_{T\times K} \tag{8.2}\]

And the entire panel itself with all units included will look like this: \[ Y = \left( \begin{array}{c} Y_{1} \\ \vdots \\ Y_{i} \\ \vdots \\ Y_{N} \end{array} \right)_{NT\times 1} \qquad D = \left( \begin{array}{c} D_1 \\ \vdots \\ D_i \\ \vdots \\ D_N \end{array} \right)_{NT \times K} \tag{8.3}\] For a randomly drawn cross-sectional unit \(i\), the model is given by \[ Y_{it}=\delta D_{it}+u_i+\varepsilon_{it}, \quad t=1,2,\dots,T \tag{8.4}\]

As always, we use our schooling-earnings example for motivation. Let \(Y_{it}\) be log earnings for a person \(i\) in year \(t\). Let \(D_{it}\) be schooling for person \(i\) in year \(t\). Let \(\delta\) be the returns to schooling. Let \(u_i\) be the sum of all time-invariant person-specific characteristics, such as unobserved ability. This is, as I said earlier, often called the unobserved heterogeneity. And let \(\varepsilon_{it}\) be the time-varying unobserved factors that determine a person’s wage in a given period. This is often called the idiosyncratic error. We want to know what happens when we regress \(Y_{it}\) on \(D_{it}\).

Pooled OLS Regression

The first estimator we will discuss is the pooled ordinary least squares or POLS estimator. When we ignore the panel structure and regress \(Y_{it}\) on \(D_{it}\) we get \[ Y_{it}=\delta D_{it} + \eta_{it}; \quad t=1,2,\dots,T \tag{8.5}\] with composite error \(\eta_{it} \equiv c_i + \varepsilon_{it}\). The main assumption necessary to obtain consistent estimates for \(\delta\) is: \[ E\big[\eta_{it}\mid D_{i1},D_{i2}, \dots, D_{iT}\big] = E\big[\eta_{it}\mid D_{it}\big] = 0 \quad \text{for $t=1,2,\dots, T$} \tag{8.6}\] While our DAG did not include \(\varepsilon_{it}\), this would be equivalent to assuming that the unobserved heterogeneity, \(c_i\), was uncorrelated with \(D_{it}\) for all time periods.

But this is not an appropriate assumption in our case because our DAG explicitly links the unobserved heterogeneity to both the outcome and the treatment in each period. Or using our schooling-earnings example, schooling is likely based on unobserved background factors, \(u_i\), and therefore without controlling for it, we have omitted variable bias and \(\widehat{\delta}\) is biased. No correlation between \(D_{it}\) and \(\eta_{it}\) necessarily means no correlation between the unobserved \(u_i\) and \(D_{it}\) for all \(t\), and that is just probably not a credible assumption. An additional problem is that \(\eta_{it}\) is serially correlated for unit \(i\) since \(u_i\) is present in each \(t\) period. And thus heteroskedastic robust standard errors are also likely too small.

Fixed Effects Regression

Let’s rewrite our unobserved effects model so that this is still firmly in our minds: \[ Y_{it} = \delta D_{it} + u_i + \varepsilon_{it}; \quad t=1,2,\dots,T \tag{8.7}\] If we have data on multiple time periods, we can think of \(u_i\) as fixed effects to be estimated. OLS estimation with fixed effects yields \[ \big(\widehat{\delta}, \widehat{u}_1, \dots, \widehat{u}_N\big) = \underset{b,m_1,\dots,m_N}{\arg\,\min} \sum_{i=1}^N\sum_{t=1}^T (Y_{it}-D_{it}b- m_i)^2, \tag{8.8}\] which amounts to including \(N\) individual dummies in regression of \(Y_{it}\) on \(D_{it}\).

The first-order conditions (FOC) for this minimization problem are: \[ \sum_{i=1}^N \sum_{t=1}^T D_{it}' \big(Y_{it} - D_{it}\widehat{\delta} - \widehat{u}_i\big)=0 \tag{8.9}\] and \[ \sum_{t=1}^T\big(Y_{it} - D_{it}\widehat{\delta} - \widehat{u}_i\big)=0 \tag{8.10}\] for \(i=1,\dots,N\).

Therefore, for \(i=1, \dots, N\), \[ \widehat{u}_i = \dfrac{1}{T} \sum_{t=1}^T\big(Y_{it}-D_{it}\widehat{\delta}\big)= \overline{Y}_i-\overline{D}_i\widehat{\delta} \tag{8.11}\] where \[ \overline{D}_i \equiv \dfrac{1}{T}\sum_{t=1}^TD_{it}; \bar{Y}_i \equiv \dfrac{1}{T} \sum_{t=1}^T Y_{it} \tag{8.12}\] Plug this result into the first FOC to obtain: \[ \begin{gather*} \widehat{\delta} = \bigg( \sum_{i=1}^N \sum_{t=1}^T (D_{it} - \overline{D}_i)'(D_{it} - \overline{D}_i) \bigg)^{-1} \bigg( \sum_{i=1}^N \sum_{t=1}^T (D_{it} - \overline{D}_i)'(Y_{it} - \overline{Y})\bigg) \\ \widehat{\delta} = \bigg(\sum_{i=1}^N \sum_{t=1}^T \ddot{D}_{it}'\ddot{D}_{it} \bigg)^{-1} \bigg( \sum_{i=1}^N \sum_{t=1}^T \ddot{D}_{it}' \ddot{D}_{it} \bigg) \end{gather*} \tag{8.13}\] with time-demeaned variables \(\ddot{D}_{it} \equiv D_{it}-\overline{D},\ddot{Y}_{it} \equiv Y_{it} - \overline{Y}_i\).

In case it isn’t clear, though, running a regression with the time-demeaned variables \(\ddot{Y}_{it}\equiv Y_{it} - \overline{Y}_i\) and \(\ddot{D}_{it} \equiv D_{it} - \overline{D}\) is numerically equivalent to a regression of \(Y_{it}\) on \(D_{it}\) and unit specific dummy variables. Hence, this is sometimes called the “within” estimator, and sometimes called the “fixed effects” estimator. And when year fixed effects are included, the “two-way fixed effects” estimator. They are the same thing.3

Even better, the regression with the time demeaned variables is consistent for \(\delta\) even when \(C[D_{it},u_i]\ne 0\) because time-demeaning eliminates the unobserved effects. Let’s see this now: \[ \begin{align*} Y_{it} &=\delta D_{it} + u_i + \varepsilon_{it} \\ \overline{Y}_i &=\delta \overline{D}_i + u_i + \overline{\varepsilon}_{i} \\ (Y_{it} - \overline{Y}_i) &=(\delta D_{it} - \delta\overline{D}) + (u_i - u_i) + (\varepsilon_{it} - \overline{\varepsilon}_i) \\ \ddot{Y}_{it}&=\delta\ddot{D}_{it}+\ddot{\varepsilon}_{it} \end{align*} \tag{8.14}\] Where’d the unobserved heterogeneity go?! It was deleted when we time demeaned the data. And as we said, including individual fixed effects does this time demeaning automatically so that you don’t have to go to the actual trouble of doing it yourself manually.4

So, how do we precisely do this form of estimation? There are three ways to implement the fixed effects (within) estimator. They are:

  1. Demean and regress \(\ddot{Y}_{it}\) on \(\ddot{D}_{it}\) (need to correct degrees of freedom)

  2. Regress \(Y_{it}\) on \(D_{it}\) and unit dummies (dummy variable regression)

  3. Regress \(Y_{it}\) on \(D_{it}\) with canned fixed effects routine in Stata or R

Later in this chapter, I will review an example from my research. We will estimate a model involving sex workers (!) using pooled OLS, a FE and a demeaned OLS model. I’m hoping this exercise will help you see the similarities and differences across all three approaches.

Identifying Assumptions

We kind of reviewed the assumptions necessary to identify \(\delta\) with our fixed effects (within) estimator when we walked through that original DAG, but let’s supplement that DAG intuition with some formality. The main identification assumptions are:

  1. \(E[\varepsilon_{it}\mid D_{i1}, D_{i2}, \dots, D_{iT}, u_i]=0; t=1,2,\dots,T\)

    • This means that the regressors are strictly exogenous conditional on the unobserved effect. This allows \(D_{it}\) to be arbitrarily related to \(u_i\), though. It only concerns the relationship between \(D_{it}\) and \(\varepsilon_{it}\), not \(D_{it}\)’s relationship to \(u_i\).
  2. \(rank\bigg( \sum_{t=1}^T E[\ddot{D}_{it}'\ddot{D}_{it}]\bigg) = K\)

    • It shouldn’t be a surprise to you by this point that we have a rank condition, because even when we were working with the simpler linear models, the estimated coefficient was always a scaled covariance, where the scaling was by a variance term. Thus regressors must vary over time for at least some \(i\) and not be collinear in order that \(\widehat{\delta} \approx \delta\).

The properties of the estimator under assumptions 1–2 are that \(\widehat{\delta}_{FE}\) is consistent (\(\underset{N\rightarrow \infty}{p\lim}\ \widehat{\delta}_{FE,N}=\delta\)) and \(\widehat{\delta}_{FE}\) is unbiased conditional on D

I only briefly mention inference. But the standard errors in this framework must be “clustered” by panel unit (e.g., individual) to allow for correlation in the \(\varepsilon_{it}\)’s for the same person \(i\) over time. This yields valid inference so long as the number of clusters is “large.”5

Do Not Use Fixed Effects For These Things

Reverse Causality

But, there are still things that fixed effects (within) estimators cannot solve. For instance, let’s say we regressed crime rates onto police spending per capita. Becker (1968) argues that increases in the probability of arrest, usually proxied by police per capita or police spending per capita, will reduce crime. But at the same time, police spending per capita is itself a function of crime rates. This kind of reverse causality problem shows up in most panel models when regressing crime rates onto police. For instance, see Cornwell and Trumbull (1994). I’ve reproduced a portion of this in Table 8.1, column 2. The dependent variable is crime rates by county in North Carolina for a panel, and they find a positive correlation between police and crime rates, which is the opposite of what Becker predicts. Does this mean that more police in an area cause higher crime rates? Or does it likely reflect the reverse causality problem?

Table 8.1: Panel Estimates of Police on Crime
Dependent variable Between Within 2SLS (FE) 2SLS (no FE)
Police 0.364 0.413 0.504 0.419
(0.060) (0.027) (0.617) (0.218)
Controls Yes Yes Yes Yes

North Carolina county level data. Standard errors in parenthesis.

So, one situation in which you wouldn’t want to use panel fixed effects is if you have reverse causality or simultaneity bias. And specifically when that reverse causality is very strong in observational data. This would technically violate the DAG, though, that we presented at the start of the chapter. Notice that if we had reverse causality, then \(Y \rightarrow D\), which is explicitly ruled out by this theoretical model contained in the DAG. But obviously, in the police—crime example, that DAG would be inappropriate, and any amount of reflection on the problem should tell you that that DAG is inappropriate. Thus it requires, as I’ve said repeatedly, some careful reflection, and writing out exactly what the relationship is between the treatment variables and the outcome variables in a DAG can help you develop a credible identification strategy.

Time-Variant Unobserved Heterogeneity

The second situation in which panel fixed effects don’t buy you anything is if the unobserved heterogeneity is time varying. In this situation, the demeaning has simply demeaned an unobserved time-variant variable, which is then moved into the composite error term, and which since time demeaned \(\ddot{u}_{it}\) correlated with \(\ddot{D}_{it}\), \(\ddot{D}_{it}\) remains endogenous. Again, look carefully at the DAG—panel fixed effects is only appropriate if \(u_i\) is unchanging. Otherwise it’s just another form of omitted variable bias. So, that said, don’t just blindly use fixed effects and think that it solves your omitted variable bias problem—in the same way that you shouldn’t use matching just because it’s convenient to do. You need a DAG, based on an actual economic model, which will allow you to build the appropriate research design. Nothing substitutes for careful reasoning and economic theory, as they are the necessary conditions for good research design.

When might this be true? Let’s use an example from Cornwell and Rupert (1997) in which the authors attempt to estimate the causal effect of marriage on earnings. It’s a well-known stylized fact that married men earn more than unmarried men, even controlling for observables. But the question is whether that correlation is causal, or whether it reflects unobserved heterogeneity (i.e., selection bias).

So, let’s say that we had panel data on individuals. These individuals \(i\) are observed for four periods \(t\). We are interested in the following equation:6 \[ Y_{it}=\alpha+\delta M_{it}+\beta X_{it}+A_{i} +\gamma_i+\varepsilon_{it} \tag{8.15}\] Let the outcome be their wage \(Y_{it}\) observed in each period, and which changes each period. Let wages be a function of marriage, which changes over time \(M_{it}\), other covariates that change over time \(X_{it}\), race and gender, which do not change over the panel period \(A_i\), and an unobserved variable we call unobserved ability \(\gamma_i\). This could be intelligence, non-cognitive ability, motivation, or some other unobserved confounder. The key here is that it is unit-specific, unobserved, and time-invariant. The \(\varepsilon_{it}\) is the unobserved determinants of wages that are assumed to be uncorrelated with marriage and other covariates.

Cornwell and Rupert (1997) estimate both a feasible generalized least squares model (FGLS) and three fixed effects models (each of which includes different time-varying controls). The authors call the fixed effects regression a “within” estimator, because it uses the within unit variation for eliminating the confounding. Their estimates are presented in Table 8.2.

Table 8.2: Estimated Wage Regressions
Dependent variable FGLS Within Within Within
Married 0.083 0.056 0.051 0.033
(0.022) (0.026) (0.026) (0.028)
Education controls Yes No No No
Tenure No No Yes Yes
Quadratics in years married No No No Yes

Standard errors in parentheses.

Notice that the FGLS (column 1) finds a strong marriage premium of around 8.3%. But, once we begin estimating fixed effects models, the effect gets smaller and less precise. The inclusion of marriage characteristics, such as years married and job tenure, causes the coefficient on marriage to fall by around 60% from the FGLS estimate, and is no longer statistically significant at the 5% level.

Sex Work and Risky Sex Data Exercise

Next I’d like to introduce a Stata exercise based on data collection for my own research: a survey of sex workers. You may or may not know this—the internet has had a profound effect on sex markets. It has moved women indoors from the streets while simultaneously breaking the link with pimps. It has increased safety and anonymity, too, which has had the effect of causing new entrants. The marginal sex worker has more education and better outside options than traditional US sex workers (Cunningham and Kendall 2011, 2014, 2017). The internet, in sum, caused the marginal sex worker demographic to shift towards women more sensitive to detection, harm, andarrest.

In 2008 and 2009, I surveyed (with Todd Kendall) approximately 700 US internet-mediated sex workers. The survey was a basic labor market survey; I asked them about their illicit and legal labor market experiences, and demographics. The survey had two parts: a “static” provider-specific section and a “panel” section. The panel section asked respondents to share information about each of the last four sessions with clients.7

I have created a shortened version of the dataset and uploaded it to Github. It includes a few time-invariant provider characteristics, such as race, age, marital status, years of schooling, and body mass index, as well as several time-variant, session-specific characteristics, including the log of the hourly price, the log of the session length (in hours), characteristics of the client himself, whether a condom was used in any capacity during the session, whether the client was a”regular,” etc.

In this exercise, you will estimate three types of models: a pooled OLS model, a fixed effects (FE), and a demeaned OLS model. The model will be of the following form: \[ \begin{align*} Y_{is}&=\beta_i X_i + \gamma_{is} Z_{is} + u_i + \varepsilon_{is} \\ \ddot{Y}_{is} &= \gamma_{is} \ddot{Z}_{is} + \ddot \eta_{is} \end{align*} \tag{8.16}\] where \(u_i\) is both unobserved and correlated with \(Z_{is}\).

The first regression model will be estimated with pooled OLS and the second model will be estimated using both fixed effects and OLS. In other words, I’m going to have you estimate the model using canned routines in Stata and R with individual fixed effects, as well as demean the data manually and estimate the demeaned regressionusing OLS.

Notice that the second regression has a different notation on the dependent and independent variable; it represents the fact that the variables are columns of demeaned variables. Thus \(\ddot Y_{is} = Y_{is} - \overline{Y}_i\). Secondly, notice that the time-invariant \(X_i\) variables are missing from the second equation. Do you understand why that is the case? These variables have also been demeaned, but since the demeaning is across time, and since these time-invariant variables do not change over time, the demeaning deletes them from the expression. Notice, also, that the unobserved individual specific heterogeneity, \(u_i\), has disappeared. It has disappeared for the same reason that the \(X_i\) terms are gone—because the mean of \(u_i\) over time is itself, and thus the demeaning deletes it.

Let’s examine these models using the following R and Stata programs.

sasp.do

Code
use https://github.com/scunning1975/mixtape/raw/master/sasp_panel.dta, clear
tsset id session
foreach x of varlist lnw age asq bmi hispanic black other asian schooling cohab married divorced separated age_cl unsafe llength reg asq_cl appearance_cl provider_second asian_cl black_cl hispanic_cl othrace_cl hot massage_cl {
    drop if `x'==.
}
    bysort id: gen s=_N
    keep if s==4


    foreach x of varlist lnw age asq bmi hispanic black other asian schooling cohab married divorced separated age_cl unsafe llength reg asq_cl appearance_cl provider_second asian_cl black_cl hispanic_cl  othrace_cl hot massage_cl {
    egen mean_`x'=mean(`x'), by(id)
    gen demean_`x'=`x' - mean_`x'
    drop mean*
}

xi: reg lnw  age asq bmi hispanic black other asian schooling cohab married divorced separated age_cl unsafe llength reg asq_cl appearance_cl provider_second asian_cl black_cl hispanic_cl othrace_cl hot massage_cl, robust

xi: xtreg lnw  age asq bmi hispanic black other asian schooling cohab married divorced separated age_cl unsafe llength reg asq_cl appearance_cl provider_second asian_cl black_cl hispanic_cl othrace_cl hot massage_cl, fe i(id) robust

reg demean_lnw demean_age demean_asq demean_bmi demean_hispanic demean_black demean_other demean_asian demean_schooling demean_cohab demean_married demean_divorced demean_separated demean_age_cl demean_unsafe demean_llength demean_reg demean_asq_cl demean_appearance_cl demean_provider_second demean_asian_cl demean_black_cl demean_hispanic_cl demean_othrace_cl demean_hot demean_massage_cl, robust cluster(id)

sasp.R

Code
library(tidyverse)
library(haven)
library(estimatr)
library(plm)

read_data <- function(df)
{
  full_path <- paste("https://raw.github.com/scunning1975/mixtape/master/", 
                     df, sep = "")
  df <- read_dta(full_path)
  return(df)
}

sasp <- read_data("sasp_panel.dta")

#-- Delete all NA
sasp <- na.omit(sasp)

#-- order by id and session 
sasp <- sasp %>% 
  arrange(id, session)

#Balance Data
balanced_sasp <- make.pbalanced(sasp, 
                                balance.type = "shared.individuals")

#Demean Data
balanced_sasp <- balanced_sasp %>% mutate( 
  demean_lnw = lnw - ave(lnw, id),
  demean_age = age - ave(age, id),
  demean_asq = asq - ave(asq, id),
  demean_bmi = bmi - ave(bmi, id),
  demean_hispanic = hispanic - ave(hispanic, id),
  demean_black = black - ave(black, id),
  demean_other = other - ave(other, id),
  demean_asian = asian - ave(asian, id),
  demean_schooling = schooling - ave(schooling, id),
  demean_cohab = cohab - ave(cohab, id),
  demean_married = married - ave(married, id),
  demean_divorced = divorced - ave(divorced, id),
  demean_separated = separated - ave(separated, id),
  demean_age_cl = age_cl - ave(age_cl, id),
  demean_unsafe = unsafe - ave(unsafe, id),
  demean_llength = llength - ave(llength, id),
  demean_reg = reg - ave(reg, id),
  demean_asq_cl = asq_cl - ave(asq_cl, id),
  demean_appearance_cl = appearance_cl - ave(appearance_cl, id),
  demean_provider_second = provider_second - ave(provider_second, id),
  demean_asian_cl = asian_cl - ave(asian_cl, id),
  demean_black_cl = black_cl - ave(black_cl, id),
  demean_hispanic_cl = hispanic_cl - ave(hispanic_cl, id),
  demean_othrace_cl = othrace_cl - ave(lnw, id),
  demean_hot = hot - ave(hot, id),
  demean_massage_cl = massage_cl - ave(massage_cl, id)
  )

#-- POLS
ols <- lm_robust(lnw ~ age + asq + bmi + hispanic + black + other + asian + schooling + cohab + married + divorced + separated + 
           age_cl + unsafe + llength + reg + asq_cl + appearance_cl + provider_second + asian_cl + black_cl + hispanic_cl + 
           othrace_cl + hot + massage_cl, data = balanced_sasp)
summary(ols)


#-- FE
formula <- as.formula("lnw ~ age + asq + bmi + hispanic + black + other + asian + schooling + 
                      cohab + married + divorced + separated + 
                      age_cl + unsafe + llength + reg + asq_cl + appearance_cl + 
                      provider_second + asian_cl + black_cl + hispanic_cl + 
                      othrace_cl + hot + massage_cl")

model_fe <- lm_robust(formula = formula,
                  data = balanced_sasp, 
                  fixed_effect = ~id, 
                  se_type = "stata")

summary(model_fe)

#-- Demean OLS
dm_formula <- as.formula("demean_lnw ~ demean_age + demean_asq + demean_bmi + 
                demean_hispanic + demean_black + demean_other +
                demean_asian + demean_schooling + demean_cohab + 
                demean_married + demean_divorced + demean_separated +
                demean_age_cl + demean_unsafe + demean_llength + demean_reg + 
                demean_asq_cl + demean_appearance_cl + 
                demean_provider_second + demean_asian_cl + demean_black_cl + 
                demean_hispanic_cl + demean_othrace_cl +
                demean_hot + demean_massage_cl")

ols_demean <- lm_robust(formula = dm_formula, 
                data = balanced_sasp, clusters = id,
                se_type = "stata")

summary(ols_demean)

A few comments about this analysis. Some of the respondents would leave certain questions blank, most likely due to concerns about anonymity and privacy. So, I dropped anyone who had missing values for the sake of this exercise. This leaves us with a balanced panel. I have organized the output into Table 8.3. There’s a lot of interesting information in these three columns, some of which may surprise you if only for the novelty of the regressions. So, let’s talk about the statistically significant ones. The pooled OLS regressions, recall, do not control for unobserved heterogeneity, because by definition those are unobservable. So, these are potentially biased by the unobserved heterogeneity, which is a kind of selection bias, but we will discuss them anyhow.

Table 8.3: POLS, FE, and Demeaned OLS Estimates of the Determinants of Log Hourly Price for a Panel of Sex Workers
Depvar: POLS FE Demeaned OLS
Unprotected sex with client of any kind 0.013 0.051* 0.051*
(0.028) (0.028) (0.026)
Ln(Length) −0.308*** −0.435*** −0.435***
(0.028) (0.024) (0.019)
Client was a regular −0.047* −0.037** −0.037**
(0.028) (0.019) (0.017)
Age of client −0.001 0.002 0.002
(0.009) (0.007) (0.006)
Age of client squared 0.000 −0.000 −0.000
(0.000) (0.000) (0.000)
client attractiveness (scale of 1 to 10) 0.020*** 0.006 0.006
(0.007) (0.006) (0.005)
Second provider involved 0.055 0.113* 0.113*
(0.067) (0.060) (0.048)
Asian client −0.014 −0.010 −0.010
(0.049) (0.034) (0.030)
Black client 0.092 0.027 0.027
(0.073) (0.042) (0.037)
Hispanic client 0.052 −0.062 −0.062
(0.080) (0.052) (0.045)
Other ethnicity client 0.156** 0.142*** 0.142***
(0.068) (0.049) (0.045)
Met client in hotel 0.133*** 0.052* 0.052*
(0.029) (0.027) (0.024)
Gave client a massage −0.134*** −0.001 −0.001
(0.029) (0.028) (0.024)
Age of provider 0.003 0.000 0.000
(0.012) (.) (.)
Age of provider squared −0.000 0.000 0.000
(0.000) (.) (.)
Body mass index −0.022*** 0.000 0.000
(0.002) (.) (.)
Hispanic −0.226*** 0.000 0.000
(0.082) (.) (.)
Black 0.028 0.000 0.000
(0.064) (.) (.)
Other −0.112 0.000 0.000
(0.077) (.) (.)
Asian 0.086 0.000 0.000
(0.158) (.) (.)
Imputed years of schooling 0.020** 0.000 0.000
(0.010) (.) (.)
Cohabiting (living with a partner) but unmarried −0.054 0.000 0.000
(0.036) (.) (.)
Currently married and living with your spouse 0.005 0.000 0.000
(0.043) (.) (.)
Divorced and not remarried −0.021 0.000 0.000
(0.038) (.) (.)
Married but not currently living with your spouse −0.056 0.000 0.000
(0.059) (.) (.)
N 1,028 1,028 1,028
Mean of dependent variable 5.57 5.57 0.00

Heteroskedastic robust standard errors in parenthesis clustered at the provider level. *\(p<0.10\), **\(p<0.05\), ***\(p<0.01\)

First, a simple scan of the second and third column will show that the fixed effects regression, which included (not shown) dummies for the individual herself is equivalent to a regression on the demeaned data. This should help persuade you that the fixed effects and the demeaned (within) estimators are yielding the same coefficients.

But second, let’s dig into the results. One of the first things we observe is that in the pooled OLS model, there is not a compensating wage differential detectable on having unprotected sex with a client.8 But, notice that in the fixed effects model, unprotected sex has a premium. This is consistent with Rosen (1986) who posited the existence of risk premia, as well as Gertler, Shah, and Bertozzi (2005) who found risk premia for sex workers using panel data. Gertler, Shah, and Bertozzi (2005), though, find a much larger premia of over 20% for unprotected sex, whereas I am finding only a mere 5%. This could be because a large number of the unprotected instances are fellatio, which carry a much lower risk of infection than unprotected receptive intercourse. Nevertheless, it is interesting that unprotected sex, under the assumption of strict exogeneity, appears to cause wages to rise by approximately 5%, which is statistically significant at the 10% level. Given an hourly wage of $262, this amounts to a mere $13 additional dollars per hour. The lack of a finding in the pooled OLS model seems to suggest that the unobserved heterogeneity was masking the effect.

Next, we look at the session length. Note that I have already adjusted the price the client paid for the length of the session so that the outcome is a log wage, as opposed to a log price. As this is a log-log regression, we can interpret the coefficient on log length as an elasticity. When we use fixed effects, the elasticity increases from \(-0.308\) to \(-0.435\). The significance of this result, in economic terms, though, is that there appears to be “volume discounts” in sex work. That is, longer sessions are more expensive, but at a decreasing rate. Another interesting result is whether the client was a “regular,” which meant that she had seen him before in another session. In our pooled OLS model, regulars paid 4.7% less, but this shrinks slightly in our fixed effects model to 3.7% reductions. Economically, this could be lower because new clients pose risks that repeat customers do not pose. Thus, if we expect prices to move closer to marginal cost, the disappearance of some of the risk from the repeated session should lower price, which it appears to do.

Another factor related to price is the attractiveness of the client. Interestingly, this does not go in the direction we may have expected. One might expect that the more attractive the client, the less he pays. But in fact it is the opposite. Given other research that finds beautiful people earn more money (Hamermesh and Biddle 1994), it’s possible that sex workers are price discriminating. That is, when they see a handsome client, they deduce he earns more, and therefore charge him more. This result does not hold up when including fixed effects, though, suggesting that it is due to unobserved heterogeneity, at least in part.

Similar to unprotected sex, a second provider present has a positive effect on price, which is only detectable in the fixed effects model. Controlling for unobserved heterogeneity, the presence of a second provider increases prices by 11.3%. We also see that she discriminates against clients of “other” ethnicity who pay 14.2% more than White clients. There’s a premium associated with meeting in a hotel that is considerably smaller when controlling for provider fixed effects by almost a third. This positive effect, even in the fixed effects model, may simply represent the higher costs associated with meeting in a hotel room. The other coefficients are not statistically significant.

Concluding Remarks

In conclusion, we have been exploring the usefulness of panel data for estimating causal effects. We noted that the fixed effects (within) estimator is a very useful method for addressing a very specific form of endogeneity, with some caveats. First, it will eliminate any and all unobserved and observed time-invariant covariates correlated with the treatment variable. So long as the treatment and the outcome varies over time, and strict exogeneity holds, then the fixed effects (within) estimator will identify the causal effect of the treatment on some outcome.

But this came with certain qualifications. For one, the method couldn’t handle time variant unobserved heterogeneity. It’s thus the burden of the researcher to determine which type of unobserved heterogeneity problem they face, but if they face the latter, then the panel methods reviewed here are not unbiased and consistent. Second, when there are strong reverse causality pathways, then panel methods are biased. Thus, we cannot solve the problem of simultaneity, such as Wright faced when estimating the price elasticity of demand, using the fixed effects (within) estimator. Most likely, we are going to have to move into a different framework when facing that kind of problem.

Still, many problems in the social sciences may credibly be caused by a time-invariant unobserved heterogeneity problem, in which case the fixed effects (within) panel estimator is useful and appropriate.


  1. A common third type of panel estimator is the random effects estimator, but in my experience, I have used it less often than fixed effects, so I decided to omit it. Again, this is not because it is unimportant. It is important. I just have chosen to do fewer things in more detail based on whether I think they qualify as the most common methods used in the present period by applied empiricists. See Wooldridge (2010) for a more comprehensive treatment, though, of all panel methods including random effects.↩︎

  2. For simplicity, I’m ignoring the time-invariant observations \(X_i\) from our DAG for reasons that will hopefully soon be made clear.↩︎

  3. One of the things you’ll find over time is that things have different names, depending on the author and tradition, and those names are often completely uninformative.↩︎

  4. Though feel free to do it if you want to convince yourself that they are numerically equivalent, probably just starting with a bivariate regression for simplicity.↩︎

  5. In my experience, when an econometrician is asked how large is large, they say “the size of your data.” But that said, there is a body of small clusters literature, and usually it’s thought that fewer than 30 clusters is too small (as a rule of thumb). So, it may be that having around 30–40 clusters is sufficient for the approaching of infinity. This will usually hold in most panel applications such as US states or individuals in the NSLY, etc.↩︎

  6. We use the same notation as used in their paper, as opposed to the \(\ddot{Y}\) notation presented earlier.↩︎

  7. Technically, I asked them to share about the last five sessions, but for this exercise, I have dropped the fifth due to low response rates on the fifth session.↩︎

  8. There were three kinds of sexual encounter—vaginal receptive sex, anal receptive sex, and fellatio. Unprotected sex is coded as any sex act without a condom.↩︎