0

JAGS でガンマ分布をパラメーター化しようとしています - 区分的線形予測子を使用していますが、モデルの実行に失敗し、次のエラー メッセージが表示されます。

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

このモデルは、timber.recovery が正規分布から引き出された場合に機能しますが、下位分位予測がゼロ未満であり、生物学的にあり得ません。glm からの出力に基づいて初期値を設定し、ゼロを描画する場合に備えて「mu」パラメーターに 0.001 を追加するなど、いくつかのトリックを試しました。どちらもエラーメッセージを解決しません。[私は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