Search code examples
rsjplot

Displaying one response level when using plot_model function of sjplot


I have successfully used the plot_model function of sjplot to plot a multinomial logistic regression model. The regression contains an outcome (Veto Type, with 3 levels) and 2 continuous predictors (Coalition Size and coalition Coalescence). I have also changed the values of coalition Coalescence in the plot_model so as to plot predicted effect outcomes based on the coalition Coalescence mean value and SDs

plot_model(model9, type = "pred", terms = c("CoalitionPercentageAmorim", "CoalescenceAmorim [0,1]"), mdrt.values = "meansd", ci.lvl=0.00, title = "Figure: Predicted Probability of Partial Vetoes of PLs", axis.title = c("Coalition Size", "Predicted Probability"), legend.title = "COALITION COALESCENCE")

enter image description here

However, I want the graph to only display one Veto Type level at a time. I have tried using rm.terms, and subsetting in multiple ways but so far I have not been successful. I would appreciate any help or advice regarding how to display only one level at a time.

Thanks!

addition

I tried using gtable as suggested here but this produced a 1/3 graph rather than one that cover the whole graph. I would still prefer a full graph for each of the levels.

enter image description here


Solution

  • Instead of trying to fiddle around with the plot object via the gtableone option would be to extract the model data using get_model_data and build your desired plot from scratch using ggplot2.

    Using some fake random example data.

    library(sjPlot)
    library(nnet)
    
    # Example Data
    set.seed(123)
    
    dat <- expand.grid(
      VetoType = c("Total Veto", "No Veto", "Partial Veto"),
      "CoalescenceAmorim" = c(0, 1)
    )
    dat <- replicate(10, dat, simplify = FALSE)
    dat <- do.call(rbind.data.frame, dat)
    dat$CoalitionPercentageAmorim[dat$CoalescenceAmorim == 0] <- 
      rnorm(3 * 10, 60, 20)
    
    dat$CoalitionPercentageAmorim[dat$CoalescenceAmorim == 1] <- 
      rnorm(3 * 10, 40, 10)
    
    model9 <- multinom(
      VetoType ~ CoalescenceAmorim + CoalitionPercentageAmorim,
      data = dat
    )
    #> # weights:  12 (6 variable)
    #> initial  value 65.916737 
    #> iter  10 value 65.390983
    #> final  value 65.390969 
    #> converged
    
    plot_model(model9,
      type = "pred",
      terms = c("CoalitionPercentageAmorim", "CoalescenceAmorim [0,1]"),
      mdrt.values = "meansd",
      ci.lvl = 0.00,
      title = "Figure: Predicted Probability of Partial Vetoes of PLs",
      axis.title = c("Coalition Size", "Predicted Probability"),
      legend.title = "COALITION COALESCENCE"
    )
    #> Data were 'prettified'. Consider using `terms="CoalitionPercentageAmorim
    #>   [all]"` to get smooth plots.
    

    
    model_data <- get_model_data(model9,
      type = "pred",
      terms = c("CoalitionPercentageAmorim", "CoalescenceAmorim [0,1]"),
      mdrt.values = "meansd",
      ci.lvl = 0.00
    )
    #> Data were 'prettified'. Consider using `terms="CoalitionPercentageAmorim
    #>   [all]"` to get smooth plots.
    
    library(ggplot2)
    
    model_data |> 
      subset(response.level == "Total Veto") |> 
      ggplot(aes(x , predicted, color = factor(group))) +
      geom_line() +
      scale_color_brewer(type = "qual", palette = 6) +
      scale_y_continuous(
        labels = scales::percent,
        limits = c(0, NA)
      ) +
      facet_wrap(~response.level)