I have a survey with 25 variables (160 observations for each var). Most are categorical, some are continuous. The variable of interest is categorical (binomial). If I do:
model <- glm(Y~., family=binomial(link='logit'), data=...)
I get a warning message:
1: glm.fit: algorithm did not converge
2: glm.fit: fitted probabilities numerically 0 or 1 occurred
I'm pretty sure there's no perfect separation, so the problem is probably overfitting right?
I tried to do two things. The first thing I did was avoid any variable that gave me that warning message, this left me with 7 variables. I performed an AIC test and got it down to three:
model_noerror <- glm(Y~sex+race+shoes)
Obviously I shouldn't be avoiding variables that give me a separation problem. So I looked it up and people suggested overfitting, which can be handled with LASSO.
So I used LASSO:
set.seed(999) cv.lasso <- cv.glmnet(x, y, family='binomial') penalty <- cv.lasso$lambda.min fit.lasso <- glmnet(x, y, family='binomial', alpha=1, lambda=penalty) coef(fit.lasso)
LASSO told me 10 variables were worth using. But if I include them all, again I get this "0 or 1 occurred" error message.
What could I do?
Here is the summary(model) with the 10 variables LASSO suggested
Coefficients: (1 not defined because of singularities) Estimate Std. Error z value Pr(>|z|) (Intercept) 2.769e+01 1.717e+05 0.000 1.000 raceB -3.756e+01 8.725e+03 -0.004 0.997 raceC -1.809e+01 6.169e+03 -0.003 0.998 genderB -2.131e+00 1.835e+00 -1.161 0.246 genderC 3.969e+01 3.638e+04 0.001 0.999 sexorD -2.608e+01 3.580e+04 -0.001 0.999 sexorE -1.496e+01 7.994e+04 0.000 1.000 collegeB -1.735e+01 7.994e+04 0.000 1.000 collegeC -1.615e+01 7.994e+04 0.000 1.000 ... (Dispersion parameter for binomial family taken to be 1) Null deviance: 94.242 on 119 degrees of freedom Residual deviance: 24.264 on 85 degrees of freedom AIC: 94.264 Number of Fisher Scoring iterations: 22
That's not everything but I think you get the point.
Here's the D matrix from the SVD of my model matrix.
x <- (model.matrix(Y~.,data=data)) svd(x)$d
 1.106585e+02 2.284868e+01 2.004170e+01 1.791028e+01 1.303047e+01 1.134027e+01 1.048113e+01 9.474703e+00 8.521352e+00 8.066918e+00 7.931767e+00  7.606502e+00 7.435396e+00 6.935058e+00 6.687999e+00 6.382451e+00 6.179303e+00 5.904702e+00 5.533401e+00 5.428218e+00 5.335237e+00 5.237369e+00  4.874816e+00 4.612060e+00 4.544942e+00 4.475750e+00 4.349221e+00 4.135874e+00 4.022473e+00 3.842192e+00 3.708770e+00 3.656491e+00 3.518284e+00  3.371328e+00 3.305039e+00 3.087673e+00 2.772643e+00 2.743743e+00 2.678787e+00 2.536824e+00 2.515302e+00 2.347022e+00 2.270346e+00 2.180056e+00  2.162843e+00 2.118737e+00 2.043311e+00 1.900974e+00 1.823667e+00 1.762206e+00 1.686792e+00 1.615979e+00 1.583439e+00 1.430641e+00 1.408800e+00  1.257654e+00 1.162283e+00 1.110813e+00 1.089611e+00 1.023927e+00 9.276928e-01 7.876539e-01 7.600454e-01 6.976948e-01 6.817119e-01 6.322815e-01  5.933869e-01 5.189669e-01 4.339124e-01 3.667706e-01 2.506821e-01 2.274977e-01 1.084797e-14 1.084797e-14
I've tried running a Principal Component Analysis. The results are not promising: https://i.imgur.com/Dt0uBHU.png It seems I can only remove a few variables
new_data <- subset(data, select =-c(Y)) library(dummies) new_data <- dummy.data.frame(data) pca.train <- new_data[1:120,] pca.test <- new_data[121:161,] pca.train <- pca.train[,!apply(pca.train, MARGIN = 2, function(x) max(x, na.rm = TRUE) == min(x, na.rm = TRUE))] prin_comp <- prcomp(pca.train,scale.=T) std_dev <- prin_comp$sdev pr_var <- std_dev^2 prop_varex <- pr_var/sum(pr_var) plot(prop_varex,xlab="Principal Component",ylab="Proportion of Variance Explained",type="b")
I got an "cannot rescale a constant/zero column to unit variance" error without the third line, so I added that in to remove constant/zero columns.
If I run a glm with the first 5 Principal Components from my PCA analysis, the auc is .8864 which is pretty good I think.
But since my overall PCA analysis looks bad, I guess I'm stuck doing univariate analysis for each variable at a time?