2016-12-15 6 views
1

Я пытаюсь реализовать модель пропорциональных опасностей Вейбулла с фракцией лечения после подхода, описанного Хуэй, Ибрагим и Синха (1999) - Новая байесовская модель для данных о выживании с выживающей фракцией. Однако я не уверен, можно ли определить случайный предел для цикла в JAGS.Можно определить случайный предел для цикла в JAGS?

library(R2OpenBUGS) 
library(rjags) 

set.seed(1234) 
censored <- c(1, 1) 
time_mod <- c(NA, NA) 
time_cens <- c(5, 7) 
tau <- 4 
design_matrix <- rbind(c(1, 0, 0, 0), c(1, 0.2, 0.2, 0.04)) 

jfun <- function() { 

    for(i in 1:nobs) { 
    censored[i] ~ dinterval(time_mod[i], time_cens[i]) 
    time_mod[i] <- ifelse(N[i] == 0, tau, min(Z)) 

    for (k in 1:N[i]){ 
     Z[k] ~ dweib(1, 1) 
    } 

    N[i] ~ dpois(fc[i]) 
    fc[i] <- exp(inprod(design_matrix[i, ], beta)) 

    } 

    beta[1] ~ dnorm(0, 10) 
    beta[2] ~ dnorm(0, 10) 
    beta[3] ~ dnorm(0, 10) 
    beta[4] ~ dnorm(0, 10) 
} 

inits <- function() { 
    time_init <- rep(NA, length(time_mod)) 
    time_init[which(!status)] <- time_cens[which(!status)] + 1 
    out <- list(beta = rnorm(4, 0, 10), 
       time_mod = time_init, 
       N = rpois(length(time_mod), 5)) 
    return(out) 
} 

data_base <- list('time_mod' = time_mod, 'time_cens' = time_cens, 
        'censored' = censored, 'design_matrix' = design_matrix, 
        'tau' = tau, 
        'nobs' = length(time_cens[!is.na(time_cens)])) 


tc1 <- textConnection("jmod", "w") 
write.model(jfun, tc1) 
close(tc1) 

# Calling JAGS 
tc2 <- textConnection(jmod) 
j <- jags.model(tc2, 
       data = data_base, 
       inits = inits(), 
       n.chains = 1, 
       n.adapt = 1000) 

Я наблюдал ошибку ниже:

Error in jags.model(tc2, data = data_base, inits = inits(), n.chains = 1, : 
    RUNTIME ERROR: 
Compilation error on line 6. 
Unknown variable N 
Either supply values for this variable with the data 
or define it on the left hand side of a relation. 
+0

Просто любопытно: BR есть? –

ответ

1

Я не совсем уверен, но я уверен, что вы не можете объявить случайное число узлов в BUGS в целом, так что это не будет специфическая сложность JAGS.

Тем не менее, вы можете обойти это.

С ОШИБКОЙ является декларативным языком вместо процедурного, достаточно объявить произвольный, но детерминированного числа узлов (скажем, «достаточно большой»), а затем связать только случайное числа их с распределение и наблюдаемые данные, оставив остальные узлы детерминированными.

После того, как вы заметили, максимальное значение N[i] (скажем N.max), вы можете передать его в качестве параметра к зазубринам, а затем изменить этот код твои:

for (k in 1:N[i]){ 
    Z[k] ~ dweib(1, 1) 
} 

в это:

for (k in 1:N.max){ 
    if (k <= N[i]){ 
    Z[k] ~ dweib(1, 1) 
    } else { 
    Z[k] <- 0 
    } 
} 

Я надеюсь, что это сделает трюк в вашем случае. Поэтому, пожалуйста, сообщите об этом последним.

Излишне говорить, что если у вас есть некоторая ненулевая, наблюдаемые данные, связанные с детерминированным Z[k], то весь ад потерять внутри Зубцов ...

+1

Спасибо. В конце концов, я нашел статистическое решение, интегрирующее мою скрытую переменную. –

+0

Это хорошо. По крайней мере, у вас есть новый трюк, который можно засунуть в рукав. –