################################################################### # Codigo construido para executar o Modelo 4 no WinBUGS a partir # # do R, atraves do pacote R2WinBUGS # # # # GUERINO PIROLLO JUNIOR - 12/2009 # ################################################################### ### Carregamento do pacote R2WinBUGS e entrada de informações require(R2WinBUGS) n.iter=as.numeric(readline("Informe o total de iterações:")) n.burnin=as.numeric(readline("Informe o tamanho do aquecimento (burn in):")) n.thin=as.numeric(readline("informe o tamanho do intervalo (thin):")) a=as.numeric(readline("informe o valor de alpha (>=0):")) b=as.numeric(readline("informe o valor de beta (>0):")) g=as.numeric(readline("informe o valor de gama (>=1):")) ### Especificacao do diretorio de trabalho e variavel monitorada parameters<-c("q.est") model.file<-system.file(package="R2WinBUGS","model","modelo4.txt") WD="c:\\scripts\\graduacao\\bugs4\\" ### Carga e criacao do arquivo de dados data=read.table(file=paste(WD,"dados.txt", sep=""),header=T) n=nrow(data) d=data$mortes e=data$expo data=list("n","d","e") ### Estabelecimento dos valores iniciais dos parametros inits1<-list(alpha=a,beta=b,gama=g) inits<-list(inits1) ### Estabelecimento de parametros para a funcao "bugs" e execucao do ### WinBUGS, gerando traces das cadeias e densidades a posteriori out.sim=bugs(data, inits, parameters, model.file, n.chains=1, n.iter, n.burnin, n.thin, bugs.directory="c:\\Program Files\\WinBUGS14\\", program=c("WinBUGS"), debug=T, working.directory=WD, codaPkg=T) for (i in 1:length(out.sim)){ output=read.coda(out.sim[i],paste(WD,"codaindex.txt", sep=""), quiet=F) plot(output, trace=T, density=T, smooth=T, auto.layout=T, ask=dev.interactive()) } ### Retorna os valores das taxas graduadas, estimativas e estatisticas ### decorrentes da implementacao do modelo e diagnosticos de convergencia logfile=bugs.log(paste(WD,"log.txt", sep="")) print(logfile) print(out.sim) #para execucao desta linha, codaPkg=F na funcao bugs plot(out.sim) #para execucao desta linha, codaPkg=F na funcao bugs