# [R] Odds ratios in logistic regression models with interaction

peter dalgaard pdalgd at gmail.com
Thu Aug 4 15:12:06 CEST 2016

```I suspect that "you can't get there from here"... a1\$fit and a2\$fit are not independent, so you can't work out the s.e. of their difference using sqrt(a1\$se.fit^2+a2\$se.fit^2).

You need to backtrack a bit and figure out how a1\$fit-a2\$fit relates to coef(fit3.16). I suspect it is actually just the age times the interaction term, but since you give no output and your code uses a bunch of stuff that I haven't got installed, I can't be bothered to check....

Once you have your desired value in the form t(a) %*% coef(...), then use the result that V(t(a) %*% betahat) == t(a) %*% vcov() %*% a  (asymptotically).

-pd

On 03 Aug 2016, at 15:08 , Laviolette, Michael <Michael.Laviolette at dhhs.nh.gov> wrote:

> I'm trying to reproduce some results from Hosmer & Lemeshow's "Applied Logistic Regression" second edition, pp. 74-79. The objective is to estimate odds ratios for low weight births with interaction between mother's age and weight (dichotomized at 110 lb.). I can get the point estimates, but I can't find an interval option. Can anyone provide guidance on computing the confidence intervals? Thanks. -Mike L.
>
> library(dplyr)
> data(birthwt, package = "MASS")
> birthwt <- birthwt %>%
>  mutate(low = factor(low, 0:1, c("Not low", "Low")),
>         lwd = cut(lwt, c(0, 110, Inf), right = FALSE,
>                   labels = c("Less than 110", "At least 110")),
>         lwd = relevel(lwd, "At least 110"))
>
> # p. 77, Table 3.16, Model 3
> fit3.16 <- glm(low ~ lwd * age, binomial, birthwt)
> # p. 78, interaction plot
> visreg::visreg(fit3.16, "age", by = "lwd", xlab = "Age",
>               ylab = "Estimated logit")
> # p. 78, covariance matrix
> vcov(fit3.16)
> # odds ratios for ages 15, 20, 25, 30
> age0 <- seq(15, 30, 5)
> df1 <- data.frame(lwd = "Less than 110", age = age0)
> df2 <- data.frame(lwd = "At least 110", age = age0)
> a1 <- predict(fit3.16, df1, se.fit = TRUE)
> a2 <- predict(fit3.16, df2, se.fit = TRUE)
> # p. 79, point estimates
> exp(a1\$fit - a2\$fit)
>
> # How to get CI's?
> # Age    OR     (95% CI)
> # ----------------------
> # 15   1.04 (0.29, 3.79)
> # 20   2.01 (0.91, 4.44)
> # 25   3.90 (1.71, 8.88)
> # 30   7.55 (1.95, 29.19)
>
>
>
>
> 	[[alternative HTML version deleted]]
>
> ______________________________________________
> R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see
> https://stat.ethz.ch/mailman/listinfo/r-help
> and provide commented, minimal, self-contained, reproducible code.

--
Peter Dalgaard, Professor,
Center for Statistics, Copenhagen Business School
Solbjerg Plads 3, 2000 Frederiksberg, Denmark
Phone: (+45)38153501
Office: A 4.23
Email: pd.mes at cbs.dk  Priv: PDalgd at gmail.com

```