Как построить график из резюме MuMIn model.avg()

Question

Есть ли способ напрямую отобразить средние итоговые результаты модели из MuMIn model.avg() для разных переменных с доверительными полосами. Ранее я использовал ggplot и ggpredict для построения графиков из реальных моделей, но мне не удалось найти способ отобразить результаты усредненных моделей.

Ясно, что я могу построить наклон и перехватить вручную, но получение точных доверительных интервалов и построение из confint() не является идеальным, и мне еще предстоит получить доверительные интервалы из интервалов, которые выглядят правильно.

library(MuMIn)
#Dummy Data
a <- seq(1:5)
set.seed(1)
b <- sample(1:100,5)
c <- sample(1:100,5)
d <-sample(1:100,5)
df <- data.frame(a,b,c,d)

Dredged <- dredge(lm(a ~ b + c + d, data=df), rank=AIC)
ModelAvg <- model.avg(Dredged, subset=delta<=2)


CI <- confint(ModelAvg, full=T) # get confidence intervals
summary(ModelAvg)


#I want to be able to create a graph for each term from the averaged output with its estimate, SE, and Confidence bands


#Output - I've only left the relevant part of the output, my actual data ends up with 5 component models
#Call:
#model.avg(object = Dredged, subset = delta <= 2)
#Component models: 
#    df logLik   AIC delta weight
#12   4  -1.32 10.63  0.00   0.69
#123  5  -1.10 12.21  1.58   0.31

#Model-averaged coefficients:  
#(full average) 
#             Estimate Std. Error Adjusted SE z value Pr(>|z|)
#(Intercept)  4.933497   1.308953    7.725454   0.639    0.523
#b            0.021946   0.010320    0.048539   0.452    0.651
#c           -0.044848   0.012076    0.067954   0.660    0.509
#d           -0.002275   0.014081    0.088694   0.026    0.980

Answer 1

Я не совсем уверен, что понимаю, почему вы сомневаетесь в выводе "confint()", и достоверность его вывода действительно является отдельным вопросом от вопроса построения графиков.

Чтобы построить график коэффициента +/- SE, прил. SE и 95% CI, попробуйте следующее. Здесь используется полное среднее значение модели, поскольку вы использовалиfull=T аргумент в CI.

График не самый красивый, но он выполняет свою работу - дайте мне знать, если вам нужен более красивый. Я не изобразил точку пересечения, потому что в этом случае оценки намного больше, чем коэффициенты, но все данные представлены в удобном для графического представления формате.

library(MuMIn)
#Dummy Data
a <- seq(1:5)
set.seed(1)
b <- sample(1:100,5)
c <- sample(1:100,5)
d <-sample(1:100,5)
df <- data.frame(a,b,c,d)

options(na.action = "na.fail") # needed for dredge to work
Dredged <- dredge(lm(a ~ b + c + d, data=df), rank=AIC)
ModelAvg <- model.avg(Dredged)
mA<-summary(ModelAvg) #pulling out model averages
df1<-as.data.frame(mA$coefmat.full) #selecting full model coefficient averages

CI <- as.data.frame(confint(ModelAvg, full=T)) # get confidence intervals for full model
df1$CI.min <-CI$`2.5 %` #pulling out CIs and putting into same df as coefficient estimates
df1$CI.max <-CI$`97.5 %`# order of coeffients same in both, so no mixups; but should check anyway
setDT(df1, keep.rownames = "coefficient") #put rownames into column
names(df1) <- gsub(" ", "", names(df1)) # remove spaces from column headers

график со всеми тремя планками погрешностей (SE, прил. SE, 95% ДИ)

ggplot(data=df1[2:4,], aes(x=coefficient, y=Estimate))+ #excluding intercept because estimates so much larger
  geom_point(size=10)+ #points for coefficient estimates
  theme_classic(base_size = 20)+ #clean graph
  geom_errorbar(aes(ymin=Estimate-Std.Error, ymax=Estimate+Std.Error), colour ="red", # SE
             width=.2, lwd=3) +
  geom_errorbar(aes(ymin=Estimate-AdjustedSE, ymax=Estimate+AdjustedSE), colour="blue", #adj SE
              width=.2, lwd=2) +
  geom_errorbar(aes(ymin=CI.min, ymax=CI.max), colour="pink", # CIs
                width=.2,lwd=1) 

Это дает следующий график. Красный - SE, синий - прил. SE и розовый - 95% доверительный интервал.

ИЗМЕНИТЬ с более красивым графиком:

ggplot(data=df1[2:4,], aes(x=coefficient, y=Estimate))+ #again, excluding intercept because estimates so much larger
      geom_hline(yintercept=0, color = "red",linetype="dashed", lwd=1.5)+ #add dashed line at zero
      geom_errorbar(aes(ymin=Estimate-AdjustedSE, ymax=Estimate+AdjustedSE), colour="blue", #adj SE
                  width=0, lwd=1.5) +
      coord_flip()+ # flipping x and y axes
      geom_point(size=8)+theme_classic(base_size = 20)+ ylab("Coefficient")