Please note this analysis is being done retrospectively (i.e., all data have been collected).
I have a binary response variable $y_i \in \{0, 1\}$ and two covariates $x_{1i}$ and $x_{2i}$. $x_{1i}$ is meant to be a control, and $x_{2i}$ is the covariate with which we'd like to perform inference against the $y_i$. Both $x_{1i}$ and $x_{2i}$ are binary explanatory variables.
With each combination of binary indicators $(x_{1i}, x_{2i}, y_i)$, counts range from 4 to about 15,000. In particular, when $y_i = 1$, the counts get extremely small (all less than 50). When $y_i = 0$, counts are all above 1,000.
My goal is to perform inference to get an idea of the association between $x_{2i}$ and $y_i$ when controlling for $x_{1i}$, and not prediction. How should one handle such a situation? Out-of-the-box logistic regression, I suspect, would not work for this case due to the small group sizes.
I would strongly prefer citations to literature that point to how to deal with such a situation, as this is for a paper.