如果这对任何人都有用,我使用下面的代码使用 ggplot 创建了我想要的图形;这可能不是最有效的方法,但它确实有效。我已经包含了另一个术语子模型,它允许您分成两种类型的模型(下图中的圆形和三角形)。
首先创建一个从模型中提取所需信息并将其存储在数据框中的函数:
get_estimates_for_coefplot<- function(mod, time, modname){
test2<-data.frame(summary(mod)$coefficients)
test2['term']<-names(coef(mod))
test2['model']<-modname
test2['estimate']<-test2$Estimate
test2['submodel']<-time
test2['std.error']<-test2$Std..Error
test2['ub']<-test2$estimate+test2$std.error
test2['lb']<-test2$estimate-test2$std.error
Newdata <- test2 %>%
filter(!grepl(".*s.*",term))
return(Newdata)
}
然后将您希望使用的所有模型发送到此函数并将它们绑定在一起:
df<-get_estimates_for_coefplot(sc.mod.env.msrhatmu, 'Time of day \n not included', 'a) Horizontal wind speed, \n RH included')
df<-rbind(df, get_estimates_for_coefplot(sc.mod.env.tkerhatmu, 'Time of day \n not included', 'b) TKE, RH included'))
modcoeff<-df[c('term','estimate','model','std.error','submodel','lb','ub')]
modcoeff <- modcoeff %>%
relabel_predictors(c(MeanSpeed = "Horizontal \n wind speed",
RH = "Relative \n Humidity",
AirTemp = "Temperature",
"I(hour^2)" = "Time of Day^2",
hour = "Time of Day",
MeanUpdraft = "Vertical \n wind speed",
TKE = "TKE"))
modcoeff$term<- factor(modcoeff$term, levels = c("Time of Day^2","Time of Day","Vertical \n wind speed","Temperature", "Relative \n Humidity","TKE", "Horizontal \n wind speed"))
modcoeff$submodel<- factor(modcoeff$submodel, levels = c('Time of day \n not included', 'Time of day \n included'))
此时,您应该有一个数据框,其中包含模型系数和您需要的所有其他内容,包括 lb 和 ub,它们是您用来制作误差线的上限和下限。现在使用 ggplot 进行绘图。
pd <- position_dodge(width=0.5)
ggplot(modcoeff,
aes(x=term,
y=estimate,
color=term,
group=std.error,
ymin=lb,
ymax=ub)) +
geom_hline(yintercept = 0, color='darkgrey') +
geom_point(aes(shape=submodel),size=2.7, position=pd, stat="identity") +
geom_errorbar(aes(width=0.2),position=pd) +
facet_wrap(~model)+
xlab("")+
ylab('Standardised coefficient value')+
theme_bw()+
theme(legend.position = c(0.83, 0.25))+
theme(legend.title = element_blank())+
guides(colour=FALSE)+
coord_flip()
这将产生以下图: