# Plotting lm and glm models with ggplot #rstats

Update
I followed the advice from Tim’s comment and changed the scaling in the sjPlotOdds-function to logarithmic scaling. The screenshots below showing the plotted glm’s have been updated.

Summary
In this posting I will show how to plot results from linear and logistic regression models (`lm` and `glm`) with ggplot. As in my previous postings on ggplot, the main idea is to have a highly customizable function for representing data. You can download all my scripts from my script page.

The inspiration source
My following two functions are based on an idea which I saw at the Sustainable Research Blog. Actually, this was a kind of starting point for me to get started with R and learn more about its data visualization facilities. After playing around some time with ggplot, I built my own function based on the script posted at Sustainable Research.

Plotting odds ratios
Plotting odds ratios gives you mainly two display styles: bars or plots (dots). First, let me show you the dot-style. Assuming you have a glm-object (in my examples, it’s called logreg) and have loaded the function `sjPlotOdds.R` (see my script page for downloads), you can plot the results like this (I have used `oddsLabels=lab` , a vector with label-strings, which are used as axis-labels. If you leave out this parameter, the variable-names from the model will be taken.):

```sjp.glm(logreg,
axisLabels.y=lab,
gridBreaksAt=0.4)```

In the above example, if you do not specifiy axis limits, the boundaries will be calculated according to the lowest and highest confidence interval, thus fitting the diagram to the highest possible „zoom“. The next example demonstrates this with bar charts:

```sjp.glm(logreg,
axisLabels.y=lab,
type="bars",
gridBreaksAt=0.4)```

Both diagrams contain model summaries in the lower right corner. You can change many visual parameters, for instance hiding the summary, changing bar colors, changing border or background colors, line and bar size etc.

If you dislike the grid bars to become narrower with increasing odds ratio values, you can use the `transformTicks` parameter, which uses exponential distances between the tick marks. This results in grid bars with (almost) equal distances. However, the tick values, of course, are accordingly set:

```sjp.glm(logreg,
axisLabels.y=lab,
transformTicks=TRUE,
gridBreaksAt=0.2,
errorBarWidth=0,
errorBarSize=1)```

Plotting betas and standardized betas of linear regressions
Quite similar is my function `sjPlotLinreg.R` which visualizes the results of linear regressions. Thus, it requires a lm-object.

```sjp.lm(linreg,
axisLimits=c(-0.5, 0.9),
axisTitle.x="beta (blue) and std. beta (red)",
sort="std",
axisLabels.y=lab,
axisLabelSize=1,
breakLabelsAt=30)```

As you can see, I have used `predictorLabelSize=1` and `breakLabelsAt=30` due to the long variable labels. By default, each label at the left axis would break into more lines, thus being narrower and worse to read. Then I used `sort="std"` to sort the predictors according to their standardized beta values (default would be ordering according to the beta values).

```sjp.lm(linreg,
axisLabels.y=lab,
axisLabelSize=1,
breakLabelsAt=30,
showStandardBeta=FALSE)```

The `showStandardBeta=FALSE` makes the red dots (standardized beta values) and their connecting line disappear.

```sjp.lm(linreg,
axisLabels.y=lab,
axisLabelSize=1,
breakLabelsAt=30,
showValues=FALSE,
showPValues=FALSE)```

This last example shows how to hide the value labels inside the diagram, so you only have the dots for beta and standardized beta coefficients.

Last remark
In between I have also updated my other scripts. For instance, the sjPlotGroupFrequencies.R function can now also plot box plots or violin plots (see examples at the end of that posting). So make sure you have the latest version from my script page.

## 45 Kommentare zu „Plotting lm and glm models with ggplot #rstats“

1. Tim Churches sagt:

Odds ratios should be plotted on a log scale, not a linear scale as you have used. See the example on the Sustainable Research blog that inspired you – that’s the correct way to do it.

1. Thanks for the hint, Tim! I will change that the next days.

2. For anybody’s interest, I have found two links related to this issue:

Best wishes
Daniel

2. Joseph sagt:

Good looking plots.

3. F sagt:

Daniel, great job. It worked perfectly for me — I’m using the sjPlotOdds function. But I have a doubt though. Is it possible to insert two or three models in the same graph — to compare the results — using the sjPlotOdds function?

1. I already though about that, too. When I have some time after Easter, I will try to figure out how to do this.

4. F sagt:

Cool. I’ll try to figure out a way to do it and will keep you posted. And, again, great job.

5. Daniel sagt:

Hi! Great work! Can you provide your GLM code? I have some troubles with my model. Iam using family = binomial with logit link. Thank you

1. ```# Datensatz einlesen
source("lib/sjImportSPSS.R")
efc <- importSPSS("NWIN-Buch/GER_Services_FU_PV_dt.sav")
efc_vars <- getVariableLabels(efc)
efc_labels <- getValueLabels(efc)

# ---
# Logistische Regression (LogIt-Modell)
# bei binary response (outcome, dependent) Variablen
# ---
# Logistic regression, also called a logit model, is used to model
# dichotomous outcome variables. In the logit model the log odds of the
# outcome is modeled as a linear combination of the predictor variables.
#
# remove missing values for predictors
df <- na.omit(data.frame(tot_sc_e=efc\$tot_sc_e,
e42dep=efc\$e42dep,
neg_c_7=efc\$neg_c_7,
c97qol1=efc\$c97qol1,
c161sex_0=efc\$c161sex_0,
c172high=efc\$c172high,
c172mid=efc\$c172mid))

logreg <- glm(factor(tot_sc_e) ~
e42dep +
neg_c_7 +
c97qol1 +
c161sex_0 +
c172high +
c172mid,

# ------------------------------------------------------
# plot results of log reg
# ------------------------------------------------------
lab <- c(paste(efc_vars['e42dep']),
paste(efc_vars['neg_c_7']),
paste(efc_vars['c97qol1']),
paste(efc_vars['c161sex_0']),
paste(efc_vars['c172high']),
paste(efc_vars['c172mid']))
source("lib/sjPlotOdds.R")
plotModelAssumptions(logreg)
plotOdds(logreg, oddsLabels=lab, type="bars", gridBreaksAt=0.4)
plotOdds(logreg, oddsLabels=lab, gridBreaksAt=0.4, errorBarWidth=0, errorBarSize=1)```
2. The data set is not yet public, but you can look at the questionnaires we used to collect the data:

You can easily find the related questions because the variable names follow the same numbering / order.

6. Owusu Michael sagt:

Hi SJ,
Just saw your code and has been trying to use it on my data. The x-axis limits is however giving me problems. If I use the format : axisLimits=c(0.8, 10.0),
gridBreaksAt=0.2). The x- axis seem not to be equally spaced and shrinks towards the tail end of the x-axis scale. Your graph however looked good and equally spaced.
My odd ratios are 1.15, 2.37, 3.7 and 6.54.
plotOdds(glm1,
oddsLabels=NULL,
axisLimits=c(0.8, 10.0),
gridBreaksAt=0.2)
Kindly advice me on what to do.

7. Aybek sagt:

I’ve tried the script and it seems like not graph the reference (intercept) factor/level. Is it correct?

1. Yes, the intercept is not plotted in both functions. However, in the model summary of the sjPlotLinreg, you find the intercept value as „y=…“. If this is not the correct notation for the intercept, please let me know and I’ll fix that.

1. Aybek sagt:

I think the graph should have all factors/categories, like in „effects“ library.

2. In my experience, the intercept value may deviate a lot from the beta coefficient values, which means you may have a large gap on the x axis between intercept and beta values (for instance, intercept of, let’s say -5, betas ranging from 2 to 3). To avoid this, I decided to not plot the intercept but provide this value in the plot annotation. Perhaps I should make it optional to include the intercept in the graph.

3. Aybek sagt:

Daniel, it would be great to have an option to plot reference factor.

2. I’ve update the script and it’s now possible to plot the intercept for the glm (sjPlotOdds). It’s not implemented for the lm-plotting-function yet, but will come soon.

8. Owusu Michael sagt:

I think your code is good especially relating to studies in epidemiology. I just used the odd ratio plots in a presentation today. Great!! keep it up.

1. Thank you for your comment! Good to know my scripts are useful to someone. 🙂

I’ve just updated some of my scripts and will upload them on my page in the course of this week, hopefully. Some bug fixes for the frequency-plotting-functions, and standardization of parameter names for all functions.

9. Perhaps I’m being dense, but I’m unable to get this script to work.

I try

``````
source('http://www.strengejacke.de/R-Stuff/sjPlotOdds_0_7.R')
plotModelAssumptions.glm(myModel.glm)``````

But I get a few errors. The function „halfnorm“ is missing, which does not appear to come from a CRAN listed package. I am also not getting the package „faraway“. Finally, there doesn’t appear to be a

``plotOdds``

function to call.

I’m probably missing a simple step or there is something slightly different about my environment that wasn’t accounted for…

10. Sorry, I should have updated my blog posts. The function name has changed, so the function call is now: sip.glm(…). See this posting for more examples.

To plot the model assumptions, you need following two packages:
car, faraway

Best wishes
Daniel

11. Rich Tyler sagt:

Thanks for this great script, it is very helpful. I was wondering if there is a way to avoid reordering the variables by the size of the coefficients. I want to compare the same predictors for different outcome variables and so it would be helpful if they were presented in the same order. I have been playing around with the script but unfortunately it has been unsuccessful. Thanks for yours or anyones help!

1. Hi Rich,
I’m planning to enhance the script so you can plot several fitted models in one plot, each outcome variable (or each model) represented by a certain color. This may take some time, so I probably will first supply an update that has a „sort“ parameter. The sjPlotLinreg-function already has such a parameter, but only offers ordering according to beta or standardized beta values, so I would update this function as well.

2. Ok, update is online, sorting should work now. Use „sortOdds“. Please refere to my script page for the updated script and make sure to read the change log first!

1. Rich Tyler sagt:

This is perfect, thanks so much!

12. Pingback: R | Pearltrees
13. Odds-ratio Grapher sagt:

I really like the sjp.glm plots. Is there any chance this can be compatible with the survey library?

I get thrown the following errors when attempting to use:

Error in seq.default(dots[[1L]][[1L]], dots[[2L]][[1L]], length = dots[[3L]][[1L]]) :
‚from‘ cannot be NA, NaN or infinite
1: In sjp.glm(logit1, axisLabels.y = unname(var.names), breakLabelsAt = 30) :
Exp. coefficients and/or exp. confidence intervals may be out of printable bounds. Consider using „axisLimits“ parameter!
2: In logLik.svyglm(fit) : svyglm not fitted by maximum likelihood.
3: In logLik.svyglm(fit) : svyglm not fitted by maximum likelihood.

1. In some cases, the exponentiated Beta-coeffients have such high values (due to exponentiation) that ggplot can’t cope with it. Perhaps it’s something ggplot can actually handle, if so, however, I haven’t found out how…
Does the survey-package in general not work with sjPlot? Or just this specific regression results?

1. exponentiated Beta-coeffients sagt:

Seems like it does have something to do with the exponentiated Beta-coeffients. Some of the other results work with the sjPlot package. Thanks for following up.

14. Hanna sagt:

I’m really glad I happened to find the sjPlot-package, it provides great and easy-to-use tools for visualizing results in epidemiological and social studies. However, I wonder whether there exists a code for removing selected variables from the graphical output, although they are included in the actual logistic regression model? I mean situations where one wants the results to be adjusted for e.g. age and gender, but do not want to show their OR’s and CI’s in the graphic output but rather just insert a footnote below the graph stating that the model has been adjusted for those variables. I have tried to search some hide/show-functions that would do this but I have not succeeded so far… Would love to find a solution.

1. It’s not implemented (hence, not possible) to do this in the current official release 1.8.2. However, I’ve implemented this feature for various plotting functions (`sjp.lm, sjp.glm, sjp.lmm, sjp.glmm, sjp.lmer and sjp.glmer`). You can download and install the current development build from GitHub: https://github.com/sjPlot/devel

1. Hanna sagt:

Thanks for the quick reply! Great to hear you have it in the development build, will have a look into that.

15. daniele Ventura sagt:

Hi Daniel, I tried to dowload your script for the function sjp.glm() but I get this message „Diese Seite konnte leider nicht gefunden werden“. Can yuo provide the correct http address? Thanks so much for your support!

1. The scripts are now all bundled in a package called sjPlot. You can download the package from CRAN and find more online manuals here.

16. Sherri sagt:

Hi Daniel, just wanted to say that you are amazing and thank you so much for sharing your lovely script. I am a public health student new to learning R and ggplot… and good data analysis in general… and I thank the stars every day that there are people like you willing to help us newbies! THANK YOU!

1. Hi Sherri, I’m glad you like my package. 🙂 I hope you are referring to the latest package-updates (sjPlot version 2.1.0). The up-to-date tutorials are on this page.

17. Herlon sagt:

Dear Daniel, i’m a beginner at R but doing with help of some friends many complex analysis. But now want to go alone. I want to do some graphics like in this paper „Earthworms increase plant production: a meta-analysis“ from – Jan Willem van Groenigen et al.- . Can you see and give me the way to do this. I saw that you help some people to write the scripts. Thank you very much and sorry for mistakes in my english i m from Brazil.

1. For meta analyses, including plots, I suggest looking at the metafor-package and homepage, where you can see examples on how to create such plots. The focus of my package lies on descriptives or estimates of regression models (excluding meta analysis).