#---------------------------------------------- # Inferenza di p di Bernoulli da una sequenza # + previsione di successi in ulteriori prove # # GdA marzo 2011 #---------------------------------------------- # function pausa pausa <- function() { cat ("\n >> Guarda il plot e dai enter per continuare\n") scan() } library(rjags) modello = "inf_r_bck.bug" # file con il modello dati <- NULL # oggetto con i dati dati$X <- 100 dati$T <- 10 dati$ErB <- 5 dati$sigmaB <- 0.01 jm <- jags.model(modello, dati) # definisce il modello update(jm, 100) # burn in catena <- coda.samples(jm, c("r"), n.iter=10000) # sampling print(summary(catena)) plot(catena) pausa() # trasformiamo il formato della catena catena.df <- as.data.frame( as.mcmc(catena) ) # plot a modo nostro hist(catena.df$r, nc=100, prob=TRUE, col='yellow', xlab='r', ylab='f(r)', main='Inferenza su r, intensity of process')