JAGS中持久边界对参数估计的影响

如何解决JAGS中持久边界对参数估计的影响

我正在尝试进行模型拟合，并从Stan开始，但是已经迁移到JAGS尝试进行一些预测，并且遇到了一些奇怪的边界影响，这些影响是我没能进入Stan的参数估计值。

我在Stan中所做的工作与在JAGS中所做的工作之间的主要区别在于：（1）在Stan中，我使用的是RK ODE求解器，而在JAGS中，我必须使用欧拉逼近法；（2）在Stan中，我正在进行全局拟合，而在JAGS中，我正在一步一步拟合。

这是我的JAGS模型脚本：

# Functional response project
## One-step-ahead model with JAGS instead of Stan because I'm stupid

# data {
#   for (i in 2:ts) {
#     New[i] <- y[i] - y[i-1]
#   }
# }

model { # Keep priors simple for Now
  r~dunif(0,5) # Intrinsic growth rate of parasite
  O~dunif(0,1) # Rate at which immune system correctly recognizes parasites
  h~dunif(0,1) # Handelling time
  b~dunif(0,100) # Immigration rate of immune cells to infection site
  c~dunif(0,1) # Growth rate of immune system as a response to infection
  u~dunif(0,1) # Natural mortality rate of immune cells
  
  #r~dunif(2,3) # Intrinsic growth rate of parasite
  #O~dunif(0,0.5) # Rate at which immune system correctly recognizes parasites
  #h~dunif(0,0.5) # Handelling time
  #b~dunif(30,40) # Immigration rate of immune cells to infection site
  #c~dunif(0,1) # Growth rate of immune system as a response to infection
  #u~dunif(0,1) # Natural mortality rate of immune cells
  
  # State vars
  P[1] <- 80
  H[1] <- 200 # Dat
  
  y_hat_P[1] <- 80
  y_hat_H[1] <- 200 # Est

  for (i in 2:L_T){
    P[i] <- P[i-1] + (P[i-1]*r - H[i-1]*(O*P[i-1]/(1 + O*P[i-1]*h)))*dt
    H[i] <- H[i-1] + (b + H[i-1]*(c*(O*P[i-1]/(1 + O*P[i-1]*h)) - u))*dt
    
    y_hat_P[i] <- P[i];
    y_hat_H[i] <- H[i]; # Expectation
    
    C_P[i] ~ dpois(y_hat_P[i]);
    C_H[i] ~ dpois(y_hat_H[i]);
  }
}

这是我的实现：

# Data org.
IC=list(r=2.5,O=0.012,h=0.075,b=35,c=0.3,u=0.41); inits=list(IC,IC,IC) 
data = SmallTS_Data
L_T = length(data$t)
C_P = as.integer(data$P); C_H = as.integer(data$H)
dt <- 0.001 # Same for all subsequnt models; small step size for "better" Euler approx.

# First trial
summary(Fit7) 
plot(Fit7)
# Trials w/out setting seed
# Also,increased iter by 10x
Fit7.2 = run.jags(#model = here('./Current_Code/OSA_JAGS_Para.R'),model = "/Users/mjarviscross/Desktop/GitHub/FcnalRespProj_2/Current_Code/OSA_JAGS_SmallTS.R",data = list('L_T' = L_T,'C_P' = C_P,'C_H' = C_H,'dt' = dt),monitor = c('r','O','h','b','c','u'),n.chains = 4,method = 'rjags',burnin = 100000,adapt = 100000,sample = 10000,inits = inits,summarise = TRUE,plots = TRUE)
summary(Fit7.2) 
plot(Fit7.2)

我尝试混合使用“ inits”，因此每条链都从一个不同的地方开始，没有运气。在上面的块中，init值是参数的真实值。

我还尝试记录/求幂一些参数，以使参数接近0的采样分布更容易，但这也无济于事。

奇怪的是，其中4/6个参数遇到边界效应（r，h，b和u），而O和c没有。

这是输出：

       Lower95       Median     Upper95         Mean           SD         Mode
r 4.977625e+00  4.994850857  4.99999888  4.992436257 0.0077408904  4.997662853
O 2.951738e-02  0.029777153  0.03030740  0.029825597 0.0002001345  0.029726595
h 1.283357e-07  0.002377487  0.01023549  0.003333945 0.0030660994  0.001091999
b 9.604828e+01 98.843977266 99.99991291 98.522219100 1.2420867280 99.477325634
c 7.300956e-02  0.082557460  0.09259498  0.082891674 0.0048552133  0.082253983
u 9.845811e-01  0.996027921  0.99999981  0.994573133 0.0049174393  0.998227142
         MCerr MC%ofSD SSeff     AC.10     psrf
r 3.120977e-04     4.0   615 0.2528493 1.016287
O 2.513034e-05    12.6    63 0.8550039 1.074630
h 4.088013e-04    13.3    56 0.8525586 1.021393
b 8.206921e-02     6.6   229 0.5700555 1.037147
c 3.548487e-04     7.3   187 0.6381600 1.051069
u 2.474836e-04     5.0   395 0.3892634 1.015181

任何帮助将不胜感激！：）

P.S。如果有人对将这种一步一步的模型转换为Stan有任何建议，我很想听听！我已经在Stan论坛上花了很短的时间运气了！