I'm interested in understanding how rms::lrm handles fitted probabilities close to 0 or 1, as in these cases it produces model coefficients that are quite different from glm with family = binomial("logit").
df <- data.frame(y = c(0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1),
x1 = c(3, 3, 4, 4, 3, 2, 5, 8, 9, 9, 9, 8, 9, 9, 9),
x2 = c(8, 7, 7, 6, 5, 6, 5, 2, 2, 3, 4, 3, 7, 4, 4))
glm(y ~ x1 + x2, data=df, family=binomial) |> coef()
#> Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
#> (Intercept) x1 x2
#> -75.204641 13.308926 -2.792796
rms::lrm(y ~ x1 + x2, data = df) |> coef()
#> Intercept x1 x2
#> -26.3382603 4.5923387 -0.7945125
Created on 2023-05-08 with reprex v2.0.2