How can I plot the difference between marginal effects for a categorical variable?

I'm trying to plot the difference between marginal effects for a dependent categorical variable. I have tried using emmeans, but I can't get what I want.

I'll try to follow the example in this vignette for emmeans for ordinal models.

I run this model with interaction:

library(ordinal)
library(emmeans)
wine.clm <- clm(rating ~ temp * contact, data = wine)

Then I plot the model:

plot_model(wine.clm, type = "pred",
           terms = c("contact", "temp"))

This isn't bad. But to simplify, I want to plot the difference between the predicted probabilities for warm and cold. This, I hope, should also highlight the interaction.

I've tried with emmeans, which gives me a difference, but only for the odd-logs latent variable.

emmeans(wine.clm, list(pairwise ~ contact|temp))

Instead, I would like to plot the difference in probabilities for each rating category .

One option is to use the marginaleffects package (disclaimer: I am the author).

See the Contrasts vignette for different ways to compare predicted probabilities. Here’s one possibility:

library(ordinal)
library(ggplot2)
library(marginaleffects)

wine.clm <- clm(rating ~ temp * contact, data = wine)

cmp <- comparisons(wine.clm,
    variables = "temp",
    newdata = datagrid(contact = unique))

ggplot(cmp, aes(x = contact, y = comparison, ymin = conf.low, ymax = conf.high)) +
    geom_pointrange() +
    facet_wrap(~group) +
    labs(x = "Contact", y = "P(Y|Warm) - P(Y|Cold)")