Calculating odds ratio when there is a confounder: contingency table vs logistic regression

Question

For a 2x2 contingency table,

$\begin{array}{ccc} \hline & Y=0 & Y=1 \\ \hline T=0 & n_{00} & n_{01} \\ \hline T=1 & n_{10} & n_{11} \\ \hline \end{array}$

calculating the odds ratio is straightforward: $\displaystyle \text{OR}=\frac{n_{11} / n_{10}}{n_{01}/n_{00}}$.

This odds ratio is also the same as that one calculated using logistic regression: one can prove that $\displaystyle e^{\beta_0} = \frac{n_{01}}{n_{00}}$ and $\displaystyle e^{\beta_T}=\frac{n_{11} / n_{10}}{n_{01}/n_{00}}$. For details, please see Appendix 2 below.

However, I am not sure if such an agreement holds when there is a confounder:

$\begin{array}{ccc} \hline & & Y=0 & Y=1 \\ \hline Z=0 & T=0 & n_{000} & n_{001} \\ \hline & T=1 & n_{010} & n_{011} \\ \hline Z=1 & T=0 & n_{100} & n_{101} \\ \hline & T=1 & n_{110} & n_{111} \\ \hline \end{array}$

In this case, my thought of calculating the odds ratio using the table is in the following way:

Conditional on $Z=0$, $\displaystyle \text{OR}_0 = \frac{n_{011}/n_{010}}{n_{001}/n_{000}}$. This can be thought as the "individual treatment effect" for people of $Z=0$.

Conditional on $Z=1$, $\displaystyle \text{OR}_1 = \frac{n_{111}/n_{110}}{n_{101}/n_{100}}$. This can be thought as the "individual treatment effect" for people of $Z=1$.

Then the "average treatment effect" is $\Pr(Z=0) \, \text{OR}_0 + \Pr(Z=1) \, \text{OR}_1$, where $\Pr(Z=0)$ can be estimated as the sum of the first two rows over the sum of all rows, and similarly for $\Pr(Z=1)$.

However, such an estimate of odds ratio does not agree with the logistic regression estimate (which by the way has no analytical solution to $e^{\beta}$ in this case). I am demonstrating my calculations on the "kidney stone data" at Appendix 1.

My questions:

1. Is my way of calculating the odds ratio correct for the 4x2 contingency table with a confounder?

2. If indeed there is no agreement between the two ways of calculating the odds ratio when there is a confounder, what justifications do we have for using logistic regression to "adjust" or "correct" for confounders and estimate the "true" causal effect of the treatment variable?

3. Is there an intuitive explanation why logistic regression agrees with contingency table at all for the 2x2 table? To me, logistic regression is a parametric model with rigid assumptions while contingency table seems nonparametric with minimal assumptions.

---

Appendix 1: Kidney stone data

$\begin{array}{ccc} \hline & & Y=0 & Y=1 \\ \hline Z=0 & T=0 & 6 & 81 \\ \hline & T=1 & 36 & 234 \\ \hline Z=1 & T=0 & 71 & 192 \\ \hline & T=1 & 25 & 55 \\ \hline \end{array}$

Ignoring the confounder $Z$, we have the typical 2x2 table:

$\begin{array}{cc} \hline & Y = 0 & Y = 1 \\ \hline T=0 & 77 & 273 \\ \hline T=1 & 61 & 289 \\ \hline \end{array}$

And both ways of calculating the odds ratio give 1.336.

But for the original 4x2 table, my way of calculating the odds ratio gives 0.328, but logistic regression gives 0.700 (which agrees with numerically solving the system of nonlinear equations at the end of the post). (In both cases, the odds ratio decreases from above 1 to below 1, hence the "Simpson's paradox".)

import numpy as np
from sklearn.linear_model import LogisticRegression


X = np.array([(0,0)]*87 + [(0,1)]*270 + [(1,0)]*263 + [(1,1)]*80)
y = np.array([0]*6+[1]*81 + [0]*36+[1]*234 + [0]*71+[1]*192 + [0]*25+[1]*55)

f = LogisticRegression(C=1e9).fit(X, y)

np.exp(f.intercept_), np.exp(f.coef_)
# (array([9.91206217]), array([[0.28349989, 0.69962134]]))

---

Appendix 2: Derivation

For the 2x2 contingency table: $\displaystyle \hat{y} = \frac{1}{1+e^{-(\beta_0 + \beta_T T)}}$

$\displaystyle \hat{\beta} = \arg \min_{\beta} \hat{R} = \arg \min_{\beta} \sum_i L(y_i, \hat{y}_i) = \arg \min_{\beta} \sum_i \Big(y_i \log \hat{y}_i + (1-y_i) \log (1-\hat{y}_i) \Big) \\= \arg \min_{\beta} \Big( n_{00} \big(1 - \frac{1}{1+e^{-\beta_0}}\big) + n_{01} \frac{1}{1+e^{-\beta_0}} + n_{10} \big(1 - \frac{1}{1+e^{-(\beta_0+\beta_T)}}\big) + n_{11} \frac{1}{1+e^{-(\beta_0+\beta_T)}} \Big)$

Setting the partial derivative of empirical risk $\hat{R}$ with respect to parameters to zero, one has:

$\displaystyle \frac{\partial \hat{R}}{\partial \beta_0} = 0 \Rightarrow \frac{-n_{00}}{1+e^{-\beta_0}} + \frac{n_{01}}{1+e^{\beta_0}} + \frac{ -n_{10}}{1+e^{-(\beta_0+\beta_T)}} + \frac{n_{11}}{1+e^{\beta_0+\beta_T}} = 0$

$\displaystyle \frac{\partial \hat{R}}{\partial \beta_T} = 0 \Rightarrow \frac{- n_{10}}{1+e^{-(\beta_0+\beta_T)}} + \frac{n_{11}}{1+e^{\beta_0+\beta_T}} = 0$

Solving them, one has: $\displaystyle e^{\beta_0} = \frac{n_{01}}{n_{00}}$ and $\displaystyle e^{\beta_T}=\frac{n_{11} / n_{10}}{n_{01}/n_{00}}$.

For the 4x2 contingency table with a confounder: $\displaystyle \hat{y} = \frac{1}{1+e^{-(\beta_0+\beta_Z Z + \beta_T T)}}$. Repeating the same procedure one has the following system of nonlinear equations:

Let $\displaystyle A_1 = \frac{-n_{000}}{1+e^{-\beta_0}} + \frac{n_{001}}{1+e^{\beta_0}}$,

$\displaystyle A_2 = \frac{-n_{010}}{1+e^{-(\beta_0+\beta_T)}} + \frac{n_{011}}{1+e^{\beta_0+\beta_T}}$,

$\displaystyle A_3 = \frac{-n_{100}}{1+e^{-(\beta_0+\beta_Z)}} +\frac{n_{101}}{1+e^{\beta_0+\beta_Z}}$,

$\displaystyle A_4 = \frac{-n_{110}}{1+e^{-(\beta_0+\beta_Z+\beta_T)}} + \frac{n_{111}}{1+e^{\beta_0+\beta_Z+\beta_T}}$,

$\displaystyle \frac{\partial \hat{R}}{\partial \beta_0} = 0 \Rightarrow A_1 + A_2 + A_3 + A_4 = 0$

$\displaystyle \frac{\partial \hat{R}}{\partial \beta_Z} = 0 \Rightarrow A_3 + A_4 = 0$

$\displaystyle \frac{\partial \hat{R}}{\partial \beta_T} = 0 \Rightarrow A_2 + A_4 = 0$

And the system of equations does not seem to have an analytical solution.

---

Edits:

@Noah suggested adding an interaction term to the model, and it works.

Logistic regression model: $\displaystyle \hat{y} = \frac{1}{1+e^{-(\beta_0 + \beta_Z Z + \beta_T T + \beta_{ZT} Z*T)}}$.

The terms $A_1$ to $A_3$ are the same for this model, while $A_4$ becomes $\displaystyle \frac{-n_{110}}{1+e^{-(\beta_0+\beta_Z+\beta_T + \beta_{ZT})}} + \frac{n_{111}}{1+e^{\beta_0+\beta_Z+\beta_T+\beta_{ZT}}}$.

And again, we have:

$\displaystyle \frac{\partial \hat{R}}{\partial \beta_0} = 0 \Rightarrow A_1 + A_2 + A_3 + A_4 = 0$

$\displaystyle \frac{\partial \hat{R}}{\partial \beta_Z} = 0 \Rightarrow A_3 + A_4 = 0$

$\displaystyle \frac{\partial \hat{R}}{\partial \beta_T} = 0 \Rightarrow A_2 + A_4 = 0$

$\displaystyle \frac{\partial \hat{R}}{\partial \beta_{ZT}} = 0 \Rightarrow A_4 = 0$

In this case, the system of equations have analytical solutions:

$\displaystyle e^{\beta_0} = \frac{n_{001}}{n_{000}}, \, e^{\beta_Z} = \frac{n_{101}/n_{100}}{n_{001}/n_{000}}, \, e^{\beta_T} = \frac{n_{011}/n_{010}}{n_{001}/n_{000}}, \, e^{\beta_{ZT}} = \frac{(n_{111}/n_{110}) (n_{001}/n_{000})}{(n_{011}/n_{010}) (n_{101}/n_{100})}$.

Answer 1

All this is explained in details in the book Modern Epidemiology by Sander Greenland. But forgive me for not being able to remember the details, the answers to your questions are:

1. Yes. It's one way to calculate the overall odds ratio. Besides the logistic regression method, there's also the Mantel-Haenszel method, widely used before the age of computers.

2. Logistic regression makes the assumption that the odds ratio in your two $Z$ groups are the same. If so, then it is the most efficient method, in the sense it reaches the Cramer-Rao bound asymptotically and is consistent. If not, then it is "wrong", but still approximates the population causal odds ratio, as you calculated (see some work by Greenland in the book).

3. In the 2x2 case, perhaps an intuitive explanation is that they are both doing maximum likelihood?