I'm trying to understand the behavior of logistic regression in high dimensional problems (i.e. when you are fitting a logistic regression to data with a high number of predictor variables).
Every time I fit a logistic regression with a high number of predictors, I get the following warning in R:
Warning message:
glm.fit: fitted probabilities numerically 0 or 1 occurred
I read somewhere that this is called the Hauck-Donner phenomenon and it is due to the fact that the fitted probabilities are too close to 0/1.
However, my hypothesis is that as you increase the number of predictors, the variances of your fitted probabilities has to increase in the logistic regression. This is because your log odds estimates is essentially a sum of random variables, and a sum of 100 similarly distributed random variables will (likely) have a larger variance than a sum of 10 rvs. Therefore, when you apply logistic regression to very high dimensional problems, your fitted probabilities will be closer to 0/1 (because higher variance), and therefore your coefficient estimates must be biased (incorrect) because of this problem? Is this hypothesis correct?
I've created a simulation with the following code to try to answer this:
genLogit <- function(n,dimens){
dimens <- floor(dimens/2)*2 #make sure dimens is even
xdata <- cbind(replicate(dimens/2,runif(n)),replicate(dimens/2,runif(n,-1,0)))
ydata <- apply(xdata,1,sum)
prob <- exp(ydata)/(1+exp(ydata))
runis <- runif(length(prob))
ydata <- ifelse(runis<prob,1,0)
model <- glm(ydata~.,data = data.frame(cbind(ydata,xdata)),family = binomial(link =
return(summary(model))
}
What the code is doing is basically simulating a logistic regression from the following:
$$\log\left(\frac{p}{1-p}\right) = U_1+U_1+\ldots+U_1 + U_{-1}+U_{-1}+\ldots+U_{-1}$$
where $U_1 = \text{Unif}(0,1)$ and $U_{-1} = \text{Unif}(-1,0)$. You can vary the number of predictors in the model as well as the number of data points generated. Then the function will fit a logistic regression to the simulated data, and you can examine the coefficients, residual deviance, fit, etc.
I understand that I all of my predictors have the same variance (which is not necessarily true when dealing with real data) but is this simulation still sufficient to prove my hypothesis?
