0

我正在运行带有 brms 包的序数 logit 模型。我想写出参数后验分布的密度。

<model>

    dfpop$pop = factor(df$pop, levels = c("extinct","<10", "<100", "<1000", ">=1000"), ordered = TRUE) 
dfpop$RDB2000pop = factor(dfpop_chenv$RDB2000pop, levels = c("<10", "<100", "<1000", ">=1000"), ordered = TRUE) 
bmodel<- brms::brm(pop ~  R2000pop + Temperature2003 + Population2003 + Volcanic.area + Agricultural.land2003 + RiverLake2003 + Seashore2003 + Protected.area + 
                         (1+RDB2000pop+Temperature2003+Population2003+ Volcanic.area + Agricultural.land2003 + RiverLake2003 + Seashore2003 + Protected.area |species_id),
             data      = dfpop,
             family    = "cumulative",
             prior   = c(set_prior("normal(0,10)", class = "b")), 
             warmup    = 200,
             iter      = 1000,
             chains    = 4,
             save_all_pars = TRUE)

我能够构建模型，但我一次只能创建一个图。以下代码尝试通过提供变量名称来运行循环函数。但是，我无法在获得的结果中绘制密度。

<plot>

model1tranformed <- ggs(bmodel)

exvars <- rownames(summary(bmodel)$fixed)

plots <- list()
for(i in 1:length(exvars)){
p<-ggplot(filter(model1tranformed,
              Parameter  %in% c("exvars[i]"), 
              Iteration > 0),
       aes(x = value))+
  geom_density(fill  = "darkgreen", 
               alpha = .5)+
  geom_vline(xintercept = 0, 
             col  = "darkred",
             size = 1)+
  scale_x_continuous(name   = "Value",
                     limits = c(-10,10)) + 
  geom_vline(xintercept = summary(bmodel)$fixed[i,3:4],
             col = "blue",
             linetype = 2) +
  theme_light() +
  theme(text = element_text(size = 20))
  labs(title = paste("pd_",exvars[i]))
  plots[[i]] = p
}
library(gridExtra)
pp <- do.call("grid.arrange", c(plots, ncol=4))
ggsave(pp)

在此处输入图像描述

4

0 回答 0