ggplot中GLM的多次比较的输出

Question

我想在ggplot中可视化GLM的多个比较的输出。在我的示例中：

#Simulation data set
set.seed(123)
d<-NULL
N<-54
d$Production <- rgamma(N,10)
d$Feature <- ifelse(d$Production >7 & d$Production<10, c("A"),ifelse(d$Production>11, 
c("B"), c("C")))  
#d$Temp<-rnorm(N,20,5)
d$Temp<-rep(1:3,18)
d<-as.data.frame(d)
#

# Gamma glm model
mG<- glm(Production~ Feature + Temp, family= Gamma, data = d)
summary(mG)
anova(mG,test="Chi")

#If for example a have interaction between Temp and Feature
d$BV <- interaction(d$Feature, d$Temp)

# GLM fit
mI  <- glm(Production ~ -1 + BV, data = d, family = "Gamma")
#summary(glht(mI, linfct = mcp(BV = "Tukey")))
cld(glht(mI, linfct = mcp(BV = "Tukey")))


# Prediction values
pred.data = data.frame(
Feature<-d$Feature,
Temp<-d$Temp
)  
pred.data$Production = predict(mG, newdata=pred.data, type="response")

# Select means
library(dplyr)
d2<-d %>%
  group_by(Feature, Temp) %>%
  summarize(Production = mean(Production, na.rm = TRUE))

#Final plot
library("ggplot2")
ggplot(d2, aes(Temp, Production, colour = Feature)) +
  geom_point() +
  stat_smooth(method = "glm", formula = y ~ x, family = Gamma)
#

现在我想在ggplot中可视化cld(glht(mI，linfct = mcp(BV =“Tukey”)的字母作为输出：

​

​

这是可能的吗？谢谢

Answer 1

试试这个：

library(tidyverse)

res <- cld(glht(mI, linfct = mcp(BV = "Tukey")))

d2 <-
  res$mcletters$Letters %>%
  tibble(
    reg = names(.),
    lett = .
  ) %>%
  mutate(
    Feature = substr(reg, 1, 1),
    Temp = substr(reg, 3, 3) %>% as.integer()
  ) %>%
  left_join(d2) %>%
  mutate(lett_pos = if_else(Feature == 'C', Production - .5, Production + .5))

ggplot(d2, aes(Temp, Production, colour = Feature)) +
  geom_point() +
  stat_smooth(method = "glm", formula = y ~ x) +
  geom_text(aes(y = lett_pos, label = lett), color = 'black')

​

​