When I first time learn multiple linear regression, I remember the interpretation of the regression coefficient is that: the marginal contribution of a specific predictor.
Now I am rethinking this interpretation and create an example to ask myself a question: is this interpretation always correct?
Suppose I have 2 dummy variables in my data set $X_i$ and $Z_i$, to make it concrete, they are defined as the following \begin{equation*} X_i=\begin{cases} &1,~\text{if people $i$ go to school before}\\ &0,~\text{otherwise} \end{cases},~~~Z_i=\begin{cases} &1,~\text{if people $i$ drink milk}\\ &0,~\text{otherwise} \end{cases} \end{equation*} Let $Y_i$ be the response variable, say the wage of people $i$. There are 4 combinations of the value of $(X_i,Z_i)$. If we can FORCE all people in the population have a specific value of $(X_i,Z_i)$, then we will be able to get a wage distribution, denoted by random variable $Y^Z_{Xi}$. So we have $Y^1_{1i}$, $Y^1_{0i}$, $Y^0_{0i}$ and $Y^0_{1i}$. Note that, without loss of generality, we can always assume that \begin{equation*} Y^1_{1i}=\mu^1_1+\epsilon^1_{1i},~Y^1_{0i}=\mu^1_0+\epsilon^1_{0i},~Y^0_{1i}=\mu^0_1+\epsilon^0_{1i},~Y^0_{0i}=\mu^0_0+\epsilon^0_{0i}\\ \end{equation*} where $\mu^X_Z$ is the mean of the wage if the people in the population ALL have value $X$ and $Z$. After rearrangement, it is easily to write $Y_i$, which is the observed wage in the real world, as \begin{align*} Y_i&=\mu^0_0+X_i(\mu^0_1-\mu^0_0)+Z_i(\mu^1_0-\mu^0_0)+X_iZ_i(\mu^1_1-\mu^0_1-\mu^1_0+\mu^0_0)\\ &+\underset{\epsilon_i}{\underbrace{\epsilon^0_{0i}+X_i(\epsilon^0_{1i}-\epsilon^0_{0i})+Z_i(\epsilon^1_{0i}-\epsilon^0_{0i})+X_iZ_i(\epsilon^1_{1i}-\epsilon^0_{1i}-\epsilon^1_{0i}+\epsilon^0_{0i})}} \end{align*} To write it clearly, we have \begin{equation} Y_i=a+bX_i+cZ_i+dX_iZ_i+\epsilon_i,\label{eq:2D} \end{equation} where $a=\mu^0_0, b=\mu^0_1-\mu^0_0, c=\mu^1_0-\mu^0_0, d=\mu^1_1-\mu^0_1-\mu^1_0+\mu^0_0$ Note that this equation decribes the TRUTH of the world (without making ANY assumption) and all of these coefficients have real meaning, e.g., $b$ is the increment of mean wage if the population have no milk but change from no school to all go to school, i.e., the effect of schooling when $Z_i=0$ (no milk). Note that whether we can consistently estimate these coefficient is another story.
Now, what if when you see the data set $\{(Y_i,X_i,Z_i)\}_{i=1}^n$, you build up the following model \begin{equation} Y_i=a'+b'X_i+c'Z_i+\epsilon'_i, \end{equation} i.e., you didn't incorporate the term $X_iZ_i$ in your model, which means that you dump this term into the error term, then here comes my questions:
(1) if you don't include the term $X_iZ_i$, is it always true that the error term will be correlated with the predictors? Intuitively, it looks yes, since how could the random variable $X$ and $XZ$ being uncorrelated? Then this could implies that OLS cannot consistently estimate $a',b',c'$. If for some reason that I don't know, $\epsilon'_i$ is uncorrelated with $X$ and $Z$, then we know that the omitted variable $XZ$ does not affect the OLS, so in this case, $a'=a,b'=b,c'=c$, but then what is the correct interpretation of these coefficients? $b'$ is still: "the increment of mean wage if the population have no milk but change from no school to all go to school"?
(2) somewhat related to (1), if initially I only have data $\{(Y_i,X_i)\}_{i=1}^n$, then I build the model $Y_i=a'' + b'' X_i+\epsilon''_i$, where $a''=\mathbb{E}[Y_{0i}]$, i.e., the mean wage of population if nobody go to school regardless of whether they drink milk and $b''=\mathbb{E}[Y_{1i}-Y_{0i}]$. Now somehow I observe the information about whether each individual drink milk or not. Then I project the error term $\epsilon''_i$ onto $Z_i$, i.e., $\epsilon''_i=c'' Z_i+\nu_i$. Then I run $Y_i=a'' + b'' X_i + c''_i Z_i + \nu_i$, what is the meaning of $a''_i$ and $b''_i$ in the new model? Are they the same as in the old model (both in the theoretical value and the interpretation)?
(3) For any dataset we have (with a response variable $Y$ and a predictors vector $W$), we can always write $Y_i=\mathbb{E}[Y_i|W_i]+e_i$, where $e_i$ is the residual, which is defined as $e_i\triangleq Y_i-\mathbb{E}[Y_i|W_i]$. By definition, $e_i$ is uncorrelated with $W_i$. In this post, $W_i$ are all dummies, $W_i=(X_i,Z_i)$, then $\mathbb{E}[Y_i|W_i]$ is linear, i.e., you could write $\mathbb{E}[Y_i|W_i]$ as $\mathbb{E}[Y_i|X_i,Z_i]=\alpha+\beta X_i+\gamma Z_i+\delta X_iZ_i+e_i$. Note that those $\alpha,\beta,\gamma,\delta$ have nothing to do with the aforementioned $a,b,c,d$. Here those Greek letters have completely different interpretation, e.g., here $\beta$ means the impact of unit change in $X$ on $\mathbb{E}[Y_i|X_i,Z_i]$ when $Z=0$, which is a pure description about the relationship between the predictors and the conditional expectation, that has nothing to do with the counterfactual effect like the interpretation of $b$. So in this case, if we omit the variable $X_iZ_i$, and model the linear functional $\mathbb{E}[Y_i|X_i,Z_i]$ as $\mathbb{E}[Y_i|X_i,Z_i]=\alpha'+\beta'X_i+\gamma'Z_i+e'_i$, then $\beta'$ just means the impact of unit change in $X$ on $Y_i$ (no "when $Z=0$"). In empirical study, which interpretation of the parameters should we adopt? the $\beta$ one or the $b$ one? Note that as $e_i$ is the residual (NOT the $\epsilon_i$ aforementioned), which is defined as $e_i\triangleq Y_i-\mathbb{E}[Y_i|W_i]$, then by definition, $e_i$ is uncorrelated with $W_i$. Hence, in this case, the OLS is always consistent, in the sense that, OLS estimator will converge to those $\beta$'s (the Greek letter), which does not necessarily equal to $b$ (the English letter). I recall that when I learned econometrics, the consistency of OLS estimator is a big chunk of lectures, so I assume the correct interpretation of the coefficients in the linear regression model should be the $b$ one, rather than the $\beta$ one, since otherwise, the OLS estimator is always consistent. Did I miss something important?