Plot the odds ratios (ORs) for a specific covariate from different logistic regression models

Perhaps you are interested in the effect of a specific covariate from several different fitted models. For example, we might be interested in the effect of gender (sex) from several logistic regression models (all with the same outcome). An illustrative graph would be to plot the OR of gender from each of the different fitted logistic regression models in the same plot.

To do this, assign a similar object for each model (e.g., summfit1, summfit2, etc.) using the summary.Design function of a lrm model fit. The needed information is the Effect, Lower 0.95, and Upper 0.95 columns of the output for the Odds Ratio row (2nd row). Then create a data frame from all three model summaries (plotinfo). Then generate the plot.

For this example, let's use our samplefile.txt data file --- not the best example, but a good illustration of the process.

library(Hmisc)
library(Design)

x<-read.table("samplefile.txt", header=TRUE)
x<-upData(x,
   labels=c(age="Age", race="Race", sex="Sex",
      weight="Weight", visits="No. of Visits",
      tx="Treatment"),
   levels=list(sex=c("Female", "Male"),
      race=c("Black", "Caucasian", "Other"),
      tx=c("Drug X", "Placebo")),
   units=c(age="years", weight="lbs."))
contents(x)

# Generate an outcome variable
x<-upData(x,
   # Specify population model for log odds that Y=1
   L = .4*(sex=='Male') + .045*(age-50) +
      (log(weight - 10)-5.2)*(-2*(sex=='Female') +
      2*(sex=='Male')),
   # Simulate binary y to have Prob(y=1) = 1/[1+exp(-L)]
   y  = ifelse(runif(n) < plogis(L), 1, 0))
ddist <- datadist(x) ; options(datadist='ddist')

fit1 <- lrm(y ~ visits + sex * (age + rcs(weight,4)), data = x, x=TRUE, y=TRUE)
fit2 <- lrm(y ~ visits + tx + sex * (age + rcs(weight,4)), data = x, x=TRUE, y=TRUE)
fit3 <- lrm(y ~ sex * (age + rcs(weight,4)), data = x, x=TRUE, y=TRUE)

summfit1 <- summary(fit1, sex = NA, est.all = FALSE)[2, 
   c("Effect", "Lower 0.95", "Upper 0.95")]
summfit2 <- summary(fit2, sex = NA, est.all = FALSE)[2, 
   c("Effect", "Lower 0.95", "Upper 0.95")]
summfit3 <- summary(fit3, sex = NA, est.all = FALSE)[2, 
   c("Effect", "Lower 0.95", "Upper 0.95")]
plotinfo <- data.frame(rbind(summfit1, summfit2, summfit3), model = 1:3)

with(plotinfo, {
   plot(Effect ~ model, axes = FALSE, xlab = "Model", ylab = "OR", 
      ylim = c(min(Lower.0.95), max(Upper.0.95)))
   axis(side = 2)
   axis(side = 1, at = 1:3, labels = paste("Model", 1:3))
   box()
   arrows(x0 = 1:3, x1 = 1:3, y0 = Lower.0.95, y1 = Upper.0.95, code = 3, angle = 90) 
})
Topic revision: r2 - 14 Nov 2006, TheresaScott
 

This site is powered by FoswikiCopyright © 2013-2017 by the contributing authors. All material on this collaboration platform is the property of the contributing authors.
Ideas, requests, problems regarding Vanderbilt Biostatistics Wiki? Send feedback