0

I have had difficulties plotting results from GLM logit models to show up as odds ratios on a log scale. Ultimately I want to obtain estimates from different models and plot the results on one graph as shown here (https://www.ctspedia.org/do/view/CTS...ClinAEGraph001). Do you have any insight?

victoria
  • 43
  • 1
  • 7
  • 3
    I think your link is broken... can you try pasting that again? – PaSTE Aug 17 '17 at 22:54
  • 2
    Please add some basic example data and code that shows where you're up to in plotting these. – Marius Aug 17 '17 at 23:01
  • Also have a look at [How to make a reproducible example in R?](https://stackoverflow.com/questions/5963269/how-to-make-a-great-r-reproducible-example) In order of preference, either (a) use built-in data, (b) share code to simulate some sample data, or (c) use `dput()` to share some of your real data in a copy-pasteable way. Also show the code you've tried and where you're stuck - from here it's not clear if you can make a plot of each model's results separately and just don't know how to combine them, if you don't know how to get the results you want out of the model, or somewhere in between. – Gregor Thomas Aug 17 '17 at 23:30
  • Apologies! http://www.jscarlton.net/images/hellbMultipleRegressions.png – victoria Aug 17 '17 at 23:33
  • ggplot(dat, aes(x = X, Y = or, ymin = lcl, ymax = ucl)) + geom_pointrange(aes(col = factor(lag)), position=position_dodge(width=xAxis)) + geom_hline(aes(yintercept = 1)) + scale_color_discrete(name = "Lag") + xlab("") – victoria Aug 17 '17 at 23:34
  • 1
    @victoria - best to edit these things into your question, code doesn't display well in comments. – Marius Aug 17 '17 at 23:43
  • And don't forget to share data! – Gregor Thomas Aug 18 '17 at 00:36

1 Answers1

6

I was making similar plots this week and found that it worked best to produce them with the confidence intervals running vertically at first, then use coord_flip() at the end to convert to horizontal:

library(tidyverse)
library(ggplot2)
library(broom)

set.seed(1234)

# Using the builtin Titanic dataset as example data for a GLM
tdf = as.data.frame(Titanic)

m1 = glm(Survived == "Yes" ~ Class + Sex, data = tdf, family = "binomial", weights = Freq)
m1_preds = tidy(m1, conf.int = TRUE, exponentiate = TRUE) %>%
    mutate(Model = "m1")
# Create modified data by mixing up the frequencies - doesn't do anything meaningful,
#   just a way to get different coefficients
tdf$FreqScrambled = sample(tdf$Freq)
m2 = glm(Survived == "Yes" ~ Class + Sex, data = tdf, 
         family = "binomial", weights = FreqScrambled)
m2_preds = tidy(m2, conf.int = TRUE, exponentiate = TRUE) %>%
    mutate(Model = "m2")

# At this point we have a table of odds ratios and confidence intervals
#   for plotting
ors = bind_rows(m1_preds, m2_preds)
ors


dodger = position_dodge(width = 0.3)
# Elements like pointrange and position_dodge only work when the outcome
#   is mapped to y, need to go through with OR set as y then flip at the
#   end
ggplot(ors, aes(y = estimate, x = term, colour = Model)) +
        geom_pointrange(aes(ymin = conf.low, ymax = conf.high),
                       position = dodger,
                       size = 1.2) +
        geom_hline(yintercept = 1.0, linetype = "dotted", size = 1) +
        scale_y_log10(breaks = c(0.1, 0.2, 0.5, 1.0, 2.0, 5.0, 10),
                      minor_breaks = NULL) +
        labs(y = "Odds ratio", x = "Effect") +
        coord_flip(ylim = c(0.1, 10)) +
        theme_bw()

The result is a forest plot that looks like this:

enter image description here

Marius
  • 58,213
  • 16
  • 107
  • 105