0

我正在尝试参数化 JAGS 中的伽马分布 - 使用分段线性预测器，但我的模型无法运行并出现以下错误消息：

Error: Error in node (ashape/(aexp(mu[59]))) Invalid parent values

当从正态分布中提取 wood.recovery 时，该模型有效，但较低的分位数预测值小于零，这在生物学上是不可能的。我尝试了一些技巧，例如将 0.001 添加到“mu”参数以防绘制零，根据 glm 的输出设置初始值；但两者都不能解决错误消息。任何见解将不胜感激[我正在使用 R2jags]。我的模型：

cat (
"model {

    # UNINFORMATIVE PRIORS
    sd_plot ~ dunif(0, 100)
    tau_plot <- 1/(sd_plot * sd_plot) 
      # precision for plot level variance

    alpha ~ dnorm(0, 1e-06) 
      # normal prior for intercept term 
    shape ~ dunif(0, 100) 
      # shape parameter for gamma     

    log_intensity ~ dnorm(0, 1e-06) 
      # uninformative prior for logging intensity

    beta_1 ~ dnorm (0, 1e-06) 
      # uninformative prior; change in slope for first segment : <=3.6 years
    beta_2 ~ dnorm (0, 1e-06) 
      # uninformative prior; change in slope for first segment : >3.6 years
    InX_1 ~ dnorm (0, 1e-06) 
      # uniformative prior for interaction between tsl and log_intensity : <=3.6 years
    InX_2 ~ dnorm (0, 1e-06) 
      # uniformative prior for interaction between tsl and log_intensity : >3.6 years

    # PLOT LEVEL RANDOM EFFECTS
    for (i in 1:nplots) {

      plot_Eff[i] ~ dnorm(0,tau_plot)

    }


    for (i in 1:Nobs) {

      # PIECEWISE LINEAR PREDICTOR     
      mu[i] <- 
        alpha + 
        beta_1 * (time.since.logged[i] * tsl.DUM1[i]) + 
        log_intensity * log.volume [i] + 
        beta_2 * (time.since.logged[i] * tsl.DUM2[i] - 3.6) + 
        beta_1 * (time.since.logged[i] * tsl.DUM2[i]) + 
        plot_Eff[plot.id[i]] + 
        InX_1 * (time.since.logged[i] * tsl.DUM1[i]) * log.volume [i] + 
        InX_2 * (time.since.logged[i] * tsl.DUM2[i] - 3.6) * log.volume[i] + 
        InX_1 * (time.since.logged[i] * tsl.DUM2[i]) * log.volume[i]


      timber.recovery[i] ~ dgamma(shape,shape/exp(mu[i])) 
        # observed recovery


      pred_timber_recovery[i] ~  dgamma(shape,shape/exp(mu[i])) 
        # posterior predictive distribution

      pearson.residual[i] <- 
        (timber.recovery[i] - pred_timber_recovery[i]) / (sqrt(timber.recovery[i]))
    }


 }",
 fill = TRUE, 
 file = "outputs/piecewise_TIMBER_MODEL_FINAL_GAMMA.txt")
4

0 回答 0