Download PDF
ads:
ANÁLISE BAYESIANA DE INTERFERÊNCIA E
TOXIDEZ EM ENSAIOS DE DILUIÇÃO SERIAL
ANDRÉA CRISTIANE DOS SANTOS
2008
ads:
Livros Grátis
http://www.livrosgratis.com.br
Milhares de livros grátis para download.
ANDRÉA CRISTIANE DOS SANTOS
ANÁLISE BAYESIANA DE INTERFERÊNCIA E TOXIDEZ EM
ENSAIOS DE DILUIÇÃO SERIAL
Tese apresentada à Universidade Federal de
Lavras, como parte das exigências do Programa
de Pós-Graduação em Estatística e Experimen-
tação Agropecuária, para a obtenção do título de
"Doutor".
Orientador
Prof. Dr. Júlio Sílvio de Sousa Bueno Filho
LAVRAS
MINAS GERAIS-BRASIL
2008
ads:
Ficha Catalográfica Preparada pela Divisão de Processos Técnicos da
Biblioteca Central da UFLA
Santos, Andréa Cristiane dos.
Análise bayesiana de interferência e toxidez em ensaios de diluição se-
rial/ Andréa Cristiane dos Santos. - - Lavras: UFLA, 2008.
97p. : il.
Tese(Doutorado) - Universidade Federal de Lavras, 2008.
Orientador: Júlio Sílvio de Sousa Bueno Filho.
Bibliografia.
1. Algoritmo Metropolis-Hastings. 2. Ensaio de diluição. 3. Inferência
bayesiana. 4. Número mais prováve. I. Universidade Federal de Lavras. II.
Título.
CDD-519.542
ANDRÉA CRISTIANE DOS SANTOS
ANÁLISE BAYESIANA DE INTERFERÊNCIA E TOXIDEZ EM
ENSAIOS DE DILUIÇÃO SERIAL
Tese apresentada à Universidade Federal de
Lavras, como parte das exigências do Programa
de Pós-Graduação em Estatística e Experimen-
tação Agropecuária, para a obtenção do título de
"Doutor".
APROVADA em 06 de março de 2008
Prof.
a
Dr.
a
Roseli Aparecida Leandro USP
Prof.
a
Dr.
a
Roberta Hilsdorf Piccoli UFLA
Prof. Dr. Daniel Furtado Ferreira UFLA
Prof.
a
Dr.
a
Thelma Sáfadi UFLA
Prof. Dr. Júlio Sílvio de Sousa Bueno Filho
UFLA
(Orientador)
LAVRAS
MINAS GERAIS-BRASIL
A Deus,
Porque dEle e por Ele, e para Ele, são todas as coisas;
glória pois a Ele eternamente. Rm 11.36
Ofereço
A minha mãe, Maria Aparecida, por sua incansável dedicação
Ao meu marido, Luiz Cláudio, por seu amor
Ao meu filho, Pedro Henrique, razão do meu viver
Dedico
AGRADECIMENTOS
À Universidade Federal de Lavras e ao Departamento de Ciências Exatas,
pela oportunidade de realizar este curso.
À CAPES, pela concessão da bolsa de estudos.
Ao professor Júlio, pela orientação e paciência.
À professora Roberta, do Departamento de Ciências dos Alimentos pela
co-orientação.
Aos professores participantes da banca de defesa, pelas contribuições dadas
a este trabalho.
Aos professores do Departamento de Ciências Exatas da UFLA, pelos en-
sinamentos transmitidos.
Aos funcionários do Departamento de Ciências Exatas da UFLA, pela
prestabilidade
Aos colegas de curso, pelo companherismo, em especial: Claudiney, Denise,
Imaculada, Lourdinha, Rosi, Renato e Waldemar .
A Imaculada, pela grande ajuda no desenvolvimento do meu trabalho.
A Nádia, pelas indispensáveis ajudas no Latex.
A Jennifer e sua família, por dedicarem carinho e cuidados ao meu filho
Pedro Henquire.
A toda minha família, pelo apoio que sempre me dedicam, em especial,
meus avós Leontina e Lázaro, por suas preciosas orações.
Enfim, a todos que contribuíram para a realização deste trabalho.
SUMÁRIO
LISTA DE TABELAS . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . i
LISTA DE FIGURAS . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . iii
RESUMO . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . vi
ABSTRACT . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . vii
1 INTRODUÇÃO . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1
2 REFERENCIAL TEÓRICO . . . . . . . . . . . . . . . . . . . . . . . . 3
2.1 Modelos probabilísticos . . . . . . . . . . . . . . . . . . . . . . . . . 3
2.1.1 Distribuição Binomial . . . . . . . . . . . . . . . . . . . . . . . . . 3
2.1.2 Distribuição de Poisson . . . . . . . . . . . . . . . . . . . . . . . . 4
2.2 Ensaios de diluição . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6
2.3 Número mais provável . . . . . . . . . . . . . . . . . . . . . . . . . . 8
2.4 Modelos amostrais para ensaios de diluição . . . . . . . . . . . . . . . 9
2.4.1 Modelo Usual . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10
2.4.2 Modelo de Toxidez . . . . . . . . . . . . . . . . . . . . . . . . . . . 12
2.4.3 Modelo de Interferência . . . . . . . . . . . . . . . . . . . . . . . . 15
2.5 Inferência bayesiana . . . . . . . . . . . . . . . . . . . . . . . . . . . 16
2.5.1 Método de Monte Carlo via cadeias de Markov . . . . . . . . . . . . 19
2.5.2 Diagnóstico de convergência . . . . . . . . . . . . . . . . . . . . . . 22
2.5.3 Fator de Bayes . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23
3 MATERIAL E MÉTODOS . . . . . . . . . . . . . . . . . . . . . . . . . 24
3.1 Material experimental . . . . . . . . . . . . . . . . . . . . . . . . . . 24
3.2 Metodologia . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25
3.2.1 Modelo Completo . . . . . . . . . . . . . . . . . . . . . . . . . . . 25
3.3 Distribuições a priori . . . . . . . . . . . . . . . . . . . . . . . . . . . 26
3.3.1 Processo de geração das cadeias de Markov . . . . . . . . . . . . . . 27
4 RESULTADOS E DISCUSSÃO . . . . . . . . . . . . . . . . . . . . . . 29
4.1 Resultados metodológicos: obtenção das distribuições a posterioris con-
junta para os modelos . . . . . . . . . . . . . . . . . . . . . . . . . . . 29
4.1.1 Distribuição a posteriori conjunta para o modelo Usual . . . . . . . . 29
4.1.2 Distribuição a posteriori conjunta para o modelo de Toxidez . . . . . 30
4.1.3 Distribuição a posteriori conjunta para o modelo de Interferência . . 30
4.1.4 Distribuição a posteriori conjunta para o modelo Completo . . . . . 31
4.2 Distribuições condicionais completas a posteriori . . . . . . . . . . . . 31
4.2.1 Obtenção da distribuição condicional completa a posteriori para o pa-
râmetro do modelo Usual . . . . . . . . . . . . . . . . . . . . . . . . 32
4.2.2 Obtenção da distribuição condicional completa a posteriori para os
parâmetros do modelo de Toxidez . . . . . . . . . . . . . . . . . . . 32
4.2.3 Obtenção da distribuição condicional completa a posteriori para os
parâmetros do modelo de Interferência . . . . . . . . . . . . . . . . . 33
4.2.4 Obtenção da distribuição condicional completa a posteriori para os
parâmetros do modelo Completo . . . . . . . . . . . . . . . . . . . . 34
4.3 Probabilidade de aceitação no algoritmo de Metropolis-Hastings . . . . 35
4.3.1 Probabilidade de aceitação no algoritmo de Metropolis-Hastings para
o modelo Usual . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 35
4.3.2 Probabilidade de aceitação no algoritmo de Metropolis-Hastings para
o modelo de Toxidez . . . . . . . . . . . . . . . . . . . . . . . . . . 36
4.3.3 Probabilidade de aceitação no algoritmo de Metropolis-Hastings para
o modelo de Interferência . . . . . . . . . . . . . . . . . . . . . . . 37
4.3.4 Probabilidade de aceitação no algoritmo de Metropolis-Hastings para
o modelo Completo . . . . . . . . . . . . . . . . . . . . . . . . . . . 38
4.4 Exemplos de aplicação, considerando cinco diluições . . . . . . . . . . 39
4.4.1 Convergência do processo de Markov e estimação do número mais
próvavel para resultado experimental (01300) . . . . . . . . . . . . . 40
4.4.2 Convergência do processo de Markov e estimação do número mais
provável para o resultado experimental (33300) . . . . . . . . . . . . 48
4.5 Exemplos de aplicação considerando três diluições . . . . . . . . . . . 56
4.5.1 Convergência do processo de Markov e estimação do número mais
provável para o resultado experimental (030) . . . . . . . . . . . . . 56
4.5.2 Convergência do processo de Markov e estimação do número mais
provável para o resultado experimental (330) . . . . . . . . . . . . . 63
4.6 Considerações gerais . . . . . . . . . . . . . . . . . . . . . . . . . . . 69
5 CONCLUSÕES . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 70
6 REFERÊNCIAS BIBLIOGRÁFICAS . . . . . . . . . . . . . . . . . . . 71
ANEXOS A . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 74
ANEXOS B . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 82
ANEXOS C . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 95
LISTA DE TABELAS
1 Interpretação do Fator de Bayes . . . . . . . . . . . . . . . . . . 24
2 Convergência dos parâmetros dos modelos para o resultado expe-
rimental (01300). . . . . . . . . . . . . . . . . . . . . . . . . . . 41
3 Estimativas das médias a posteriori e intervalos de credibilidade
(I.C.95%), dos parâmetros dos modelos e erro de Monte Carlo
(Erro M.C.) para o resultado (01300). . . . . . . . . . . . . . . . 47
4 Fator de Bayes aproximado e grau de evidência em favor do mo-
delo da margem esquerda, em comparação ao modelo da margem
superior, considerando o resultado experimental (01300). . . . . . 48
5 Convergência dos parâmetros dos modelos para o resultado expe-
rimental (33300). . . . . . . . . . . . . . . . . . . . . . . . . . . 49
6 Estimativas das médias a posteriori e intervalos de credibilidade
(I.C.95%), dos parâmetros dos modelos e erro de Monte Carlo
(Erro M.C.) para o resultado (33300). . . . . . . . . . . . . . . . 54
7 Fator de Bayes aproximado e grau de evidência em favor do mo-
delo da margem esquerda, em comparação ao modelo da margem
superior, considerando o resultado experimental (33300). . . . . . 55
8 Convergência dos parâmetros dos modelos para o resultado expe-
rimental (030). . . . . . . . . . . . . . . . . . . . . . . . . . . . 57
9 Estimativas das médias a posteriori e intervalos de credibilidade
(I.C.95%), dos parâmetros dos modelos e erro de Monte Carlo
(Erro M.C.) para o resultado (030). . . . . . . . . . . . . . . . . . 61
10 Fator de Bayes aproximado e grau de evidência em favor do mo-
delo da margem esquerda, em comparação ao modelo da margem
superior, considerando o resultado experimental (030). . . . . . . 62
i
11 Convergência dos parâmetros dos modelos para o resultado expe-
rimental (330). . . . . . . . . . . . . . . . . . . . . . . . . . . . 63
12 Estimativas das médias a posteriori e intervalos de credibilidade
(I.C.95%), dos parâmetros dos modelos e erro de Monte Carlo
(Erro M.C.) para o resultado (330). . . . . . . . . . . . . . . . . . 67
13 Fator de Bayes aproximado e grau de evidência em favor do mo-
delo da margem esquerda, em comparação ao modelo da margem
superior, considerando o resultado experimental (330). . . . . . . 68
ii
LISTA DE FIGURAS
1 Traço da cadeia gerada (a) e densidade a posteriori (b) do parâme-
tro λ, do modelo Usual, para o resultado experimental (01300). . . 42
2 Traço das cadeias geradas (a) e (c) e densidades a posteriori (b)
e (d) dos parâmetros λ, e λ
i
, do modelo de Interferência, para o
resultado experimental (01300). . . . . . . . . . . . . . . . . . . 42
3 Traço da cadeia gerada (a) e densidade a posteriori (b) do parâme-
tro λ, do modelo de Toxidez, para o resultado experimental (01300). 42
4 Traço das cadeias geradas (a), (c), (e), (g) e (i) e densidades a
posteriori (b), (d), (f), (h) e (j) dos parâmetros T
1
, T
2
, T
3
, T
4
e T
5
,
do modelo de Toxidez, para o resultado experimental (01300). . . 43
5 Traço das cadeias geradas (a) e (c) e densidades a posteriori (b),
(d) do parâmetro λ e λ
i
, respectivamente, para o modelo Com-
pleto, para o resultado experimental (01300). . . . . . . . . . . . 44
6 Traço da cadeia gerada (a) e densidade a posteriori (b) do parâ-
metros T
1
, do modelo Completo, para o resultado experimental
(01300). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 44
7 Traço das cadeias geradas (a), (c), (e) e (g) e densidades a poste-
riori (b), (d), (f), (h) e (j) dos parâmetros T
2
, T
3
, T
4
e T
5
do modelo
Completo, para o resultado experimental (01300). . . . . . . . . . 45
8 Traço da cadeia gerada (a) e densidade a posteriori (b) do parâme-
tro λ, do modelo Usual, para o resultado experimental (33300). . . 50
9 Traço das cadeias geradas (a) e (c) e densidades a posteriori (b)
e (d) dos parâmetros λ, e λ
i
, do modelo de Interferência, para o
resultado experimental (33300). . . . . . . . . . . . . . . . . . . 50
iii
10 Traço da cadeia gerada (a) e densidade a posteriori (b) do parâme-
tro λ, para o modelo de Toxidez, segundo o resultado experimental
(33300). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 50
11 Traço das cadeias geradas (a), (c), (e), (g) e (i) e densidades a
posteriori (b), (d), (f), (h) e (j) dos parâmetros T
1
, T
2
, T
3
, T
4
e T
5
,
do modelo de Toxidez, para o resultado experimental (33300). . . 51
12 Traço das cadeias geradas (a) e (c) e densidades a posteriori (b),
(d) do parâmetro λ e λ
i
, respectivamente, do modelo Completo,
para o resultado experimental (33300). . . . . . . . . . . . . . . . 52
13 Traço da cadeia gerada (a) e densidade a posteriori (b) do parâme-
tro T 1, do modelo Completo, para o resultado experimental (33300). 52
14 Traço das cadeias geradas (a), (c), (c) e (f) e densidades a poste-
riori (b), (d), (f) e (h) dos parâmetros T
2
, T
3
, T
4
e T
5
do modelo
Completo, para o resultado experimental (33300). . . . . . . . . . 53
15 Traço da cadeia gerada (a) e densidade a posteriori (b) do parâme-
tro λ, do o modelo Usual, para o resultado experimental (030). . . 57
16 Traço das cadeias geradas (a) e (c) e densidades a posteriori (b)
e (d) dos parâmetros λ, e λ
i
, do modelo de Interferência, para o
resultado experimental (030). . . . . . . . . . . . . . . . . . . . . 58
17 Traço da cadeia gerada (a) e densidade a posteriori (b) do parâme-
tro λ, do modelo de Toxidez, para o resultado experimental (030). 58
18 Traço das cadeias geradas (a), (c) e (e) e densidades a posteriori
(b), (d) e (f) dos parâmetros T
1
, T
2
e T
3
, do modelo de Toxidez,
para o resultado experimental (030). . . . . . . . . . . . . . . . . 59
iv
19 Traço das cadeias geradas (a) e (c) e densidades a posteriori (b)
e (d) do parâmetro λ e λ
i
, respectivamente, do modelo Completo,
para o resultado experimental (030). . . . . . . . . . . . . . . . . 60
20 Traço das cadeias geradas (a), (c) e (e) e densidades a posteriori
(b), (d) e (f) dos parâmetros T
1
, T
2
e T
3
do modelo Completo, para
o resultado experimental (030). . . . . . . . . . . . . . . . . . . . 60
21 Traço da cadeia gerada (a) e densidade a posteriori (b) do parâme-
tro λ, do modelo Usual, para o resultado experimental (330). . . . 64
22 Traço das cadeias geradas (a) e (c) e densidades a posteriori (b)
e (d) dos parâmetros λ, e λ
i
, do modelo de Interferência, para o
resultado experimental (330). . . . . . . . . . . . . . . . . . . . . 64
23 Traço da cadeia gerada (a) e densidade a posteriori (b) do parâme-
tro λ, do modelo de Toxidez, para o resultado experimental (330). 65
24 Traço das cadeias geradas (a), (c) e (d) e densidades a posteriori
(b), (d) e (f) dos parâmetros T
1
, T
2
e T
3
, do modelo de Toxidez,
para o resultado experimental (330). . . . . . . . . . . . . . . . . 65
25 Traço das cadeias geradas (a) e (c) e densidades a posteriori (b),
(d) do parâmetro λ e λ
i
, respectivamente, do modelo Completo,
para o resultado experimental (330). . . . . . . . . . . . . . . . . 66
26 Traço das cadeias geradas (a) e (c) e densidades a posteriori (b)
e (d) dos parâmetros T
1
, T
2
e T
3
, do modelo Completo, para o
resultado experimental (330). . . . . . . . . . . . . . . . . . . . . 66
v
RESUMO
SANTOS, Andréa Cristiane. Análise bayesiana de interferência e toxidez em
ensaios de diluição serial. 2008. 97p. Tese (Doutorado em Estatística e Experi-
mentação Agropecuária) - Universidade Federal de Lavras, Lavras, MG.
*
O método do número mais provável (NMP) foi introduzido por McCrady
(1915) como resultado da estimativa da densidade de microrganismos provenientes
de diluições seriadas. O modelo padrão utilizado para obter tabelas de estimati-
vas do NMP não contempla situações tais como Toxidez e interferência. A análise
bayesiana foi utilizada para realizar inferências sobre estes modelos de acordo com
alguns resultados experimentais. Por meio do algoritmo de Metropolis-Hastings,
foi possível gerar cadeias para cada parâmetro dos modelos. A convergência
dessas cadeias foi monitorada por meio de analises gráficas e pelos critérios de
Raftery & Lewis (1992) e Heidelberg & Welch (1983). A comparação do ajuste
dos modelos foi feita por meio do Fator de Bayes. Foi possível discriminar, entre
modelos com parâmetros de interferência e Toxidez mais precisos que o modelo
Usual sob certos resultados experimentais. Em situações experimentais em que o
crescimento do microrganismo de interesse é inibido em baixas diluições, o mo-
delo Usual não é adequado e tem valores subestimados para o NMP. No caso em
que crescimento de microrganismo em baixas diluições, o modelo Usual tem
estimativa para o NMP próxima das apresentadas nos modelos mais complexos,
embora os intervalos de credibilidade sejam maiores. A rotina implementada em
R pode ser extendida para uma ampla gama de planos de diluição (número de
diluições e número de tubos por diluição) e pode ser usada em substituição das
tabelas nos laboratórios.
Palavras Chave: Algoritmo Metropolis-Hastings; Ensaio de diluição; Inferência
bayesiana; Número Mais Provável.
*
Comitê Orientador: Prof.
o
Dr.
o
Júlio Sílvio de Sousa Bueno Filho - UFLA (Orientador) e
Prof.
a
Dr.
a
Roberta Hilsdorf Piccoli - UFLA (co-orientadora)
vi
ABSTRACT
SANTOS, Andréa Cristiane. Bayesian Analysis of interfering and toxicity ef-
fects in dilution assays. 2008. 97p. Thesis (Doctorate in Statistics and Agri-
cultural Experimentation) - Federal University of Lavras, Lavras, Minas Gerais,
Brazil.
*
Estimation of microorganism densities by means of the Most Probable
Number (MPN) is a technique introduced by McCrady (1915) to analyse serial di-
lution assays. The standard model used to generate MPN tables does not consider
medium toxicity nor interference due to competitor microorganisms. In this work
we aim to develop a Bayesian framework to analyze these phenomena. MCMC
methods using Metropolis-Hastings algorithm were used to get posterior distri-
butions given some experimental results. Convergence was monitored using gra-
phical display and both Raftery & Lewis (1992) and Heidelberg & Welsh (1983)
criteria. Model comparison was done using Bayes Factors. It was possible to sort
out models with interfering and toxicant parameters that were more probable than
standard model for some experimental results. When microorganism do not grow
in initial dilutions, the standard model underestimates MPN. In the situations in
which standard model is the most probable, MPN estimates from any model are
similar, although standart model is the best with smaller credibility interval. A
very flexible R routine was implemented. It can manage a wide range of dilution
designs with more dilutions and more tubes per dilution and is a suitable tool for
replacing standard tables in laboratory.
Key-words: MCMC; Bayesian Inference; Dilution Assay; Metropolis-Hastings
Algorithm; Microbiology; Most Probable Number.
*
Guidance Committee: Prof. DSc. Júlio Sílvio de Sousa Bueno Filho - UFLA (Adviser) e
Prof. DSc. Roberta Hilsdorf Piccoli - UFLA (Co-Adviser)
vii
1 INTRODUÇÃO
A enumeração de microrganismos é uma técnica de grande importância e muito
utilizada em procedimentos laboratoriais. Na análise de alimentos processados, a
presença de um determinado microrganismo é considerada um indicativo de con-
taminação pós-sanificação ou pós-processamento, e de que as práticas de higieni-
zação estão aquém dos padrões estabelecidos (Silva & Junqueira, 1997).
Um delineamento experimental muito útil para estimar a concentração de
microrganismos em alimentos é o ensaio de diluições seriadas. Este processo con-
siste em realizar sucessivas diluições, a partir de uma amostra inicial, em que se
mede a concentração de determinados indicadores associados a microrganismos.
Na versão mais simples e mais utilizada deste tipo de ensaio, medem-se apenas a
presença ou a ausência de microrganismos em um ou mais tubos de uma diluiçao.
Ao interpretar os resultados obtidos nestes ensaios, usando um modelo
probabilístico hierárquico, calcula-se o estimador de máxima verossimilhança do
número mais provável de microganismos. O método de obtenção de tais estima-
tivas a partir desses delineamentos ficou conhecido simplesmente como número
mais provável. Este conceito foi introduzido por McCrady (1915) e, desde então,
tem sido amplamente utilizado.
As estimativas de número mais provável para três diluições, considerando-
se 3, 5, 8 e 10 tubos, podem ser encontradas em tabelas que são usadas de forma
rotineira em laboratórios de microbiologia. Tais tabelas são adequadas para um
número limitado de resultados experimentais, os quais apresentem crescimento de
microrganismos em baixas diluições e diminuição gradativa de tubos que mostrem
crescimento, à medida que aumenta a diluição das amostras.
O modelo probabilístico mais direto utilizado nestas situações não se ajusta
a resultados experimentais que fogem do padrão da tabela, como, por exemplo para
1
resultados em que não crescimento de microrganismos em baixas diluições. A
inibição do crescimento de microrganismos pode estar associada a fatores, como
toxidez no meio de cultura e interferência por competição entre espécies, além de
outros de difícil explicação. Blodgett & Garthright (1998) sugerem modelos al-
ternativos que seriam mais prováveis para alguns desses resultados experimentais,
sem, no entanto, apresentar métodos integrados de estimação e avaliação do ajuste
em tais modelos.
Não na literatura, nenhum enfoque bayesiano para o cálculo do número
mais provável, que envolva os modelos descritos por Blodgett & Garthright (1998).
No entanto, é possível encontrar estudos sobre estimativas de microrganismos,
utilizando-se inferência bayesiana, em técnicas de plaqueamento. Artigos como
os de Gelman at al. (2004) e Clough (2005), tratam deste assunto.
Além de apresentar um novo enfoque para o cálculo do número mais
provável, neste trabalho também é proposto um modelo para obter estimativas de
microrganismos, apresentado na seção 3.2.1.
O objetivo deste trabalho foi desenvolver uma abordagem bayesiana para
alguns modelos utilizados para determinar o número mais provável, que consi-
derem a possibilidade de interferência de outros microrganismos e de toxidez do
meio de cultura. Espera-se estabelecer uma rotina flexível que permita avaliar ex-
perimentos de uma ampla gama de planos de diluição e torne dispensável o uso de
tabelas na avaliação de ensaios posteriores.
2
2 REFERENCIAL TEÓRICO
2.1 Modelos probabilísticos
É comum o uso de variáveis na forma de contagem em pesquisas, nas mais
diversas áreas do conhecimento, em particular, na área de ciências de alimentos
e, mais especificamente, em microbiologia de alimentos, é rotineira nos proces-
sos laboratoriais a enumeração de microrganismos em uma dada cultura. Nesta
situação, o uso de modelos discretos é mais indicado para descrever o fenômeno,
no entanto, o uso de aproximações normais em modelos contínuos é bastante co-
mum. Nas seções 2.1.1 e 2.1.2 deste trabalho são abordadas algumas distribuições
discretas de probabilidade que foram utilizadas neste trabalho.
2.1.1 Distribuição Binomial
Quando uma variável é classificada com base em apenas dois resultados
possíveis, tais como presença ou ausência de crescimento de microrganismo em
tubo ou placa, planta atacada ou não por fungos, define-se a chamada distribuição
de Bernoulli. Em tais situações, pode-se definir uma variável aleatória que assume
valores iguais a 0 ou 1, associados aos dois resultados possíveis. O evento a partir
do qual a variável assume valor 1 é chamado de sucesso e a probabilidade de este
evento ocorrer é chamada de probabilidade de sucesso.
Diferentes ensaios de Bernoulli envolvem, cada um, uma variável aleatória.
Por exemplo, seja a probabilidade π de um tubo apresentar crescimento de mi-
crorganismos. Neste caso, a observação de um tubo corresponde a um ensaio de
Bernoulli.
Supondo que n tubos são observados com relação ao crescimento do mi-
crorganismo e sendo a probabilidade de um tubo apresentá-lo igual a π, de tal
3
modo que esta seja a mesma para cada tubo, então, pode-se definir uma variável
aleatória associada ao número de tubos que apresentam crescimento (Y=0, 1, 2, 3,
. . . , n). Esta tal variável aleatória apresenta distribuição de probabilidade chamada
binomial, com parâmetro n e π.
Como cada tubo está associado a uma distribuição de Bernoulli, então,
a soma das realizações dessas variáveis corresponde exatamente ao número de
sucessos. Por isso, a distribuição binomial também é definida como sendo a soma
de n variáveis independentes com distribuição de Bernoulli.
Uma variável Y é definida como tendo uma distribuição binomial, se a
densidade de Y é dada por:
P (Y = y) =
n
y
π
y
(1 π)
(ny)
(2.1)
em que:
n é o tamanho da amostra;
π é a probabilidade de sucesso;
1 π é a probabilidade de fracasso;
n
y
é o número combinatório.
A média e a variância da distribuição binomial são dadas, respectivamente
por: E(Y)= e V(Y)=(1 π). Ver por exemplo, Mood et al. (1974)
2.1.2 Distribuição de Poisson
Supondo que, em um laboratório de microbiologia, microrganismos são
cultivados em tubos e deseja-se quantificar a ocorrência de um microrganismo es-
pecífico, sendo essa ocorrência um evento raro, a probilidade π para sua determi-
nação é muito baixa. Além disso, o número de microrganismos por tubo é em geral
4
elevado. Em casos como este, é indicada a utilização do modelo (distribuição) de
Poisson.
O modelo Poisson é o caso limite da distribuição binomial quando o nú-
mero de tentativas n tende para o infinito e a probabilidade π do evento em cada
tentativa é próxima de zero. Em essência, a distribuição de Poisson é a distribuição
binomial adequada para eventos independentes e raros, ocorrendo em um período
praticamente infinito de intervalos. Destaca-se que a unidade de medida é contínua
(em geral, o tempo, a distância, a área, o volume ou outra unidade), mas a variável
aleatória (número de ocorrência) é discreta.
Uma variável Y é definida como tendo uma distribuição de Poisson, se a
densidade de Y é dada por
P (Y = y) =
e
λ
λ
y
y!
(2.2)
em que:
λ (λ > 0) é o número esperado e e = 2, 71183 . . . é o número neperiano
A média e a variância da distribuição Poisson são dadas, respectivamente,
por: E(Y)= λ e V(Y)=λ. Ver por exemplo, Mood et al.(1974)
Seja X(X = 1, 2, . . . , n) o número de microrganismos em um tubo, com
volume z, então, o número esperado de microrganismos neste volume é E(X) =
λz. Sendo assim, a probabilidade de se encontrar x microrganismos no tubo é
dada por:
P (X = x) =
e
λz
(λz)
x
x!
(2.3)
5
2.2 Ensaios de diluição
Ensaios de diluição são utilizados para estimar a concentração de micror-
ganismos (número por unidade de volume, de área, de peso, etc.) em uma amostra
quando a contagem direta não é possível, mas a presença ou a ausência do mi-
crorganismo em subamostras podem ser detectadas (Cochran, 1950). Não são
necessárias, portanto, as contagens dos organismos e, em geral, registra-se apenas
a presença ou a ausência destes. Pode-se detectar se um determinado microrga-
nismo está presente, ou não, em um líquido por um teste de cor, ou se um fungo
está presente, ou não, em uma amostra de solo, plantando-se uma planta suscep-
tível nesse solo e verificando se a planta apresenta sintomas da doença (Demétrio,
2002).
O método consiste em fazer sucessivas diluições a partir de uma amostra
inicial, observando-se a ausência ou a presença de crescimento de microrganismo.
De modo geral, nos procedimentos de diluição seriada, faz-se a homo-
geneização do material que será analisado em um diluente adequado. No caso
de líquidos, eles são colocados em frascos com espaço suficiente para a agitação.
Embora existam várias exceções, que dependem do material a ser analisado e da
quantidade de microrganismos existentes, normalmente, utilizam-se 25 ml do ma-
terial a ser analisado e 225 ml de diluente na homogeneização. A escolha do
diluente depende do material de análise. Depois da homogeneização, transfere-
se uma porção de 1 ml da amostra homogeneizada para um tubo contendo 9 ml
de diluente. Esta primeira diluição é denotada por 10
1
. Como é um processo
de diluições sucessivas, para a diluição 10
2
, é retirado 1 ml da diluição 10
1
, o
qual é colocado em um tubo com 9 ml do diluente; para a diluição 10
3
, o procedi-
mento é o mesmo e assim sucessivamente, sempre tirando 1 ml da diluição anterior
e inoculando-se em 9 ml de diluente. Em outras palavras, pode-se dizer que, na
6
primeira diluição, a quantidade inoculada é 10
1
g ou ml; na segunda diluição, a
quantidade inoculada é 10
2
g ou ml e, na terceira diluição, a quantidade inocu-
lada é 10
3
g ou ml. A mesma forma é utilizada para realizar as demais diluições
após a terceira. A quantidade inoculada da amostra determina a qual diluição ela
pertence, de acordo com o números de diluições utilizadas.
Após o processo de diluições sucessivas, realiza-se o teste presuntivo, que
consiste em inocular uma série de três tubos de Durhan por diluição no meio de
cultura adequado, como por exemplo, caldo lauril sulfato triptose (LST). Este meio
de cultura é utilizado para favorecer o crescimento do microrganismo que será
analisado. Em cada um desses três tubos será adicionado 1 ml da solução diluída a
10 ml de LST. Sendo assim, tem-se uma série nove tubos, ou seja, de três diluições
e, para cada diluição, três tubos. Nesta fase, conta-se a quantidade de tubos que
foram positivos, ou seja, que apresentaram crescimento (Silva & Junqueira, 1997).
Supondo um ensaio experimental com 3 diluições e 3 tubos para cada di-
luição, um possível resultado seria (321). Este resultado significa que, na primeira
diluição, houve crescimento de microrganismo em 3 tubos; na segunda diluição,
apenas 2 tubos mostraram crescimento e, na terceira somente 1 tubo mostrou
crescimento.
De modo geral, em ensaios laboratoriais, utilizam-se três diluições e, para
cada diluição, três tubos. A utilização deste número de diluição e tubos é favore-
cida, pois, para estes, existem estimativas de número mais provável tabeladas.
Como ocorrem sucessivas diluições, de modo geral, uma diminuição
na quantidade de tubos com evidência de crescimento de microrganismo. Por esta
razão, algumas vezes, ensaios de diluição seriada são chamados de método de
diluição de extinção (Loyer & Hamilton, 1984).
7
2.3 Número mais provável
O conceito de número mais provável (NMP), foi introduzido por McCrady
(1915) para estimar a densidade de microrganismos resultantes de ensaios de di-
luição. Desde então, tem sido amplamente utilizado.
Fisher (1922), citado por McCullagh & Nelder (1989), discute estimação
da densidade de microrganismos utilizando-se diluições seriadas, baseada no mo-
delo de Poisson. Segundo Loyer & Hamilton (1984), o método utilizado para
estimar a quantidade de microrganismos é o de máxima verossimilhança e o valor
estimado é dito número mais provável (NMP).
Os resultados obtidos no teste presuntivo (fase na qual conta-se a quanti-
dade de tubos que apresentaram crescimento) ou no teste confirmativo são com-
parados com uma tabela de Número Mais Provável, que fornecerá a estimativa do
número de microrganismos da substância não diluída. No Manual de Bacterio-
logia Analítica (Bacteriological Analytical Manual),(AOAC International, 2001)
são encontradas tabelas para resultados experimentais com 3 diluições, para 3, 5,
8 e 10 tubos; as quantidades inoculadas da amostra são 10
1
, 10
2
e 10
3
. Nes-
tas tabelas são apresentadas as estimativa do número mais provável por ponto e
por intervalo. As tabelas de NMP, neste trabalho, serão denominadas de tabelas
padrão, as quais são o resultado do método de máxima verossimilhança e estão no
anexo C.
Quando o pesquisador utiliza em seu experimento quantidades inoculadas
que não estão tabeladas, lança-se mão de artifícios para obter as estimativas do
NMP, fazendo-se associações às tabelas padrão. Por exemplo, se as quantidades
inoculadas utilizadas são 10
2
, 10
3
e 10
4
, para 3 tubos em cada diluição, então,
multiplica-se a estimativa de NMP tabelada por 10, para que o resultado obtido
se equipare ao da tabela padrão utilizada. Por meio desta multiplicação, uma
8
compensação para a diluição não utilizada (AOAC International, 2001).
Para o exemplo em questão, a tabela padrão que seria consultada é a que
combina resultados de 3 tubos e 3 diluições, a mesma pode ser consultada no anexo
C. Supondo que o resultado de uma diluição seriada é (310) e, consultando-se a
tabela de NMP para o referido resultado, tem-se que o valor estimado de micror-
ganismo é 43 por grama ou ml. Para se obter o valor estimado de microrganismos
para as quantidades inoculadas 10
2
, 10
3
e 10
4
, utilizando-se a tabela padrão,
multiplica-se 43 por 10. Sendo assim, para este resultado da diluição seriada, a
estimativa do NMP é 430.
As estimativas de NMP são utilizadas para verificar se um alimento pode
ou não ser consumido. Para ter uma conclusão e interpretação dos resultados das
análises microbiológicas de alimentos, a Agência Nacional de Vigilância Sanitária
(Anvisa), dispõe de tabelas implementadas do Regulamento da Anvisa (2001), nas
quais são encontradas valores de tolerância máxima da ocorrência de um determi-
nado microrganismo em alimentos destinados ao consumo humano.
Até este item apresentou-se o NMP e o processo experimental pelo qual é
obtido. Nos próximos formaliza-se matematicamente o processo de estimação do
NMP.
2.4 Modelos amostrais para ensaios de diluição
De modo geral, para se obter estimativas de NMP, determina-se a probabi-
lidade de crescimento nos tubos, em seguida, a função de verossimilhança e, por
último, o valor das estimativas dos parâmetros.
Para especificar modelos hierárquicos que representem ensaios de diluição,
devem-se ser levadas em conta duas probabilidades, que são: a probabilidade de
crescimento de microrganismo no tubo e a probabilidade do número de tubos que
9
apresentam crescimento de microrganismo na diluição.
A probabilidade de tubos que apresentam crescimento por diluição é dada
pela distribuição binomial e a probabilidade de ocorrência de microrganismos no
tubo é dada por meio da distribuição Poisson, como descrito nas seções 2.1.1 e
2.1.2, respectivamente.
2.4.1 Modelo Usual
Segundo Blodgett & Garthright (1998), as seguintes suposições são neces-
sárias para a definição do modelo:
1. os tubos são independentes;
2. os microrganismos devem estar aleatoriamente distribuídos nos tubos. O
número de microrganismo no tubo segue uma distribuição de Poisson;
3. o tubo mostrará crescimento de microrganismos se estes estiverem inocula-
dos no tubo;
Supondo um ensaio de diluição em que k denota o número de diluições
seriadas; z
j
, j = 1, . . . k, é a quantidade inoculada da amostra na j-ésima diluição
e λ a concentração de microrganismos na substância original.
Sendo X o número de microrganismos dentro do tubo, determina-se, por
meio da distribuição Poisson, a probabilidade de não ocorrer crescimento, ou seja,
P (X = 0) =
e
λz
j
(λz
j
)
0
0!
= e
λz
j
(2.4)
Seja y o número de tubos que apresentam crescimento na j-ésima diluição
Y = (y
1
, y
2
. . . y
k
), então Y é binomialmente distribuído, com probabilidade de
sucesso dada por:
10
π = 1 - (probabilidade de não ocorrer crescimento) = 1 e
λz
j
Considerando n
j
tubos na diluição z
j
, de acordo com a distribuição bino-
mial a função de verossimilhança é definida por:
L(Y |λ) =
k
j=1
n
j
y
j
(1 e
λz
j
)
y
j
e
λz
j
(n
j
y
j
)
(2.5)
Utilizando-se o método da máxima verossimilhança encontra-se o valor
de λ, que maximiza a equação (2.5), que é o método que determina as estimativas
pontuais as quais podem ser encontradas nas tabelas padrão.
Blodgett & Garthright (1998) chamam a expressão (2.5) de "modelo Usual",
que é a denominação usada neste trabalho.
O modelo Usual é adequado a situações nas quais o crescimento do mi-
crorganismo é favorecido em baixas diluições.
Segundo Blodgett & Garthright (1998), situações em que o ensaio tem
um resultado não esperado. Isso ocorre, por exemplo, quando, em baixas diluições,
não crescimento de microrganismos, ocorrendo crescimento apenas em altas
diluições. Quando esta situação ocorre, pode ser uma indicação de que o meio
de cultura, por algum motivo, não está favorecendo o crescimento dos micror-
ganismos de interesse. Sendo assim, a utilização do modelo Usual para estimar
a quantidade de microrganismos não seria adequada. Os autores sugerem dois
modelos alternativos para ensaios de diluição com tubos, considerando estas situ-
ações não contemplada pelo modelo Usual. Estes modelos são o de Toxidez e o de
Interferência.
Segundo Cochran (1950), em baixas diluições, não crescimento de mi-
crorganismos e o NMP teria resultados subestimados, pois uma das condições para
que a utilização do processo de estimação fosse eficaz era o crescimento de mi-
11
crorganismos nos tubos.
2.4.2 Modelo de Toxidez
Este modelo assume que o caldo inicial contém uma substância tóxica para
o microrganismo de interesse. Como ocorrem sucessivas diluições a concentração
desta substância diminuirá com o aumento das diluições, criando uma situação
propícia para o crescimento do microrganismo, o que se dará apenas em altas
diluições.
As duas primeiras suposições feitas no modelo Usual também são válidas
para este modelo.
Além disso duas novas suposições são acrescentadas;
1. a substância tóxica é uniformemente dispersa no caldo;
2. o efeito da substância tóxica em cada microrganismo no tubo é indepen-
dente.
Neste modelo, a terceira suposição é a de que o tubo mostrará crescimento
do microrganismo se estes não forem inibidos pela substância tóxica. A probabili-
dade de que a substância tóxica iniba o crescimento do microrganismo na diluição
z
j
é indicada pro T
j
Da mesma maneira que no modelo Usual, determina-se, primeiramente, a
probabilidade de não haver crescimento no tubo. Segundo Blodgett & Garthright
(1998), não haverá crescimento de microrganismo no tubo da diluição z
j
se a
substância tóxica inibir o crescimento de cada microrganismo no tubo.
12
Esta probabilidade é dada por:
P (X = 0) =
i=0
P r(haver i microrganismos)T
i
j
=
=
i=0
e
λz
j
(λz
j
)
i
i!
T
i
j
=
i=0
e
λz
j
(λz
j
T
j
)
i
i!
=
= e
λz
j
i=0
(λz
j
T
j
)
i
i!
(2.6)
Desenvolvendo o último somatório da equação (2.6), tem-se o seguinte
resultado:
i=0
(λz
j
T
j
)
i
i!
= 1 +
λz
j
T
j
1
+
(λz
j
T
j
)
2
2
+ . . . +
(λz
j
T
j
)
i
i!
+ . . . (2.7)
Antes de prosseguir com os cálculos, considera-se, a expansão da série de
Taylor para a função f(x)= e
x
em torno de c = 0. Para isto, admite-se a notação
f
k
(c) para a k-ésima derivada de f em c. Sendo assim tem-se
e
x
= f(c) +
f
(c)(x c)
1!
+
f

(c)(x c)
2
2!
+
f

(c)(x c)
3
3!
+ . . . (2.8)
Considerando-se que em c = 0 f
(0) = f

(0) = f

(0) = . . . = e
0
= 1 e
substituindo-se na equação (2.8) tem-se
e
x
= e
0
+
1(x 0)
1!
+
1(x 0)
2
2!
+
1(x 0)
3
3!
+ . . . (2.9)
13
Simplificando o resultado da equação (2.9) obtém-se:
e
x
= 1 +
x
1!
+
x
2
2!
+
x
3
3!
+ . . . (2.10)
Fazendo x = λz
j
T
j
e substituindo em (2.10) o resultado será:
e
λz
j
T
j
= 1 +
λz
j
T
j
1!
+
(λz
j
T
j
)
2
2!
+
(λz
j
T
j
)
3
3!
+ . . . (2.11)
Assim, verifica-se que a equação (2.11) é igual ao resultado obtido por
meio do desenvolvimento do somatório da equação (2.7), e portanto pode-se rees-
crever a equação (2.6) da seguinte maneira:
P (X = 0) = e
λz
j
i=0
(λz
j
T
j
)
i
i!
= e
λz
j
e
λz
j
T
j
= e
λz
j
(1T
j
)
(2.12)
A equação (2.12) representa a probabilidade de não haver crescimento no
tubo.
Sendo Y a variável aleatória que denota o número de tubos que apresentam
crescimento de microrganismos, e considerando o modelo binomial, a função de
verossimilhança deste modelo é dada por:
L(Y |λ, T
j
) =
k
j=1
n
j
y
j
(1 e
λz
j
(1T
j
)
)
y
j
× (e
λz
j
(1T
j
)
)
(n
j
y
j
)
(2.13)
14
Pode-se notar que o que diferem a expressão (2.5) da expressão (2.13) é
a grandeza (1 T
j
) multiplicando o expoente. Se T
j
for igual a zero, a equação
(2.13) se reduz a equação (2.5).
Em seu artigo, Blodgett & Garthright (1998) consideram a substância tó-
xica como um valor fixo.
2.4.3 Modelo de Interferência
Para este modelo, assume-se que, no mesmo tubo, estejam crescendo duas
espécies de microrganismos. Seja λ a concentração de microrganismos da espé-
cie de interesse e λ
i
a concentração da espécie que interfere no crescimento do
microrganismo de interesse. Esta interferência pode ser predação, competição por
nutrientes ou outro mecanismo.
A primeira suposições feita no modelo Usual não se altera para este mo-
delo. A segunda suposição de que o número de microrganismo segue uma dis-
tribuição de Poisson é válida para ambos os microrganismos. Supõe-se que as
duas espécies sejam distribuídas independentemente.
A terceira suposição é a de que um tubo mostrará crescimento se não hou-
ver nenhum outro microrganismo que possa estar interferindo em seu crescimento.
A probabilidade de crescimento do microrganismo de interesse é a mesma
do modelo Usual e é dada por : 1 e
λz
j
.
Por outro lado a probabilidade de não haver uma espécie interferindo é
dada por:
P (X = 0) =
e
λ
i
z
j
(λ
i
z
j
)
0
0!
= e
λ
i
z
j
(2.14)
Supondo que as duas espécies sejam distribuídas independentemente, a
probabilidade de crescimento do microrganismo de interesse em um tubo na
15
j-ésima diluição é dada pela multiplicação da probabilidade de não haver cresci-
mento de uma espécie que produz algum tipo de interferência pela probabilidade
de crescimento do microrganismo de interesse. Sendo assim, a probabilidade de
observar crescimento do microrganismo de interesse é:
e
λ
i
z
j
(1 e
λz
j
)
A função de verossimilhança do modelo é dada por:
L(Y |λ, λ
i
) =
k
j=1
n
j
y
j
[e
λ
i
z
j
(1 e
λz
j
)]
y
j
×[1 e
λ
i
z
j
(1 e
λz
j
)]
(n
j
y
j
)
(2.15)
Se λ
i
for igual a zero, a equação (2.15) se reduz à expresão (2.5).
2.5 Inferência bayesiana
Na inferência bayesiana modela-se a presença de observações y, por uma
densidade f(y|θ). A quantidade θ serve como indexador da família de distribuições
das observações, representando uma característica de interesse (Gamerman, 1997).
O verdadeiro valor do parâmetro θ é desconhecido, porém, a medida da
incerteza a respeito de θ pode variar, dependendo do contexto ou do pesquisador,
sendo de interesse, portanto, tentar reduzir essa incerteza. Para a inferência baye-
siana, o parâmetro θ também é considerado uma variável aleatória, possuindo dis-
tribuição de probabilidade associada, enquanto que, para a inferência frequentista,
os parâmetros são valores fixos ou constantes. A inferência bayesiana usa toda a
informação disponível a priori, associando à sua análise a densidade a priori de θ,
p(θ), a qual é especificada antes da observação do valor de y (Box & Tiao, 1992).
A abordagem clássica permite que as informações ou inferências sobre estas quan-
16
tidades sejam feitas a partir de informações fornecidas por uma amostra aleatória
selecionada da população.
As distribuições a priori podem ser informativas ou não informativas. As
prioris informativas são usadas quando se tem algum conhecimento prévio sobre
o parâmetro, e as não informativas são usadas quando se tem pouca ou nenhuma
informação sobre o parâmetro (Paulino et al., 2003).
A informação a priori de que dispomos sobre θ, resumida probabilisti-
camente por p(θ), pode ser aumentada observando-se uma quantidade aleatória
relacionada com θ. Esta quantidade aleatória, representada pela função de veros-
similhança, é combinada com a distribuição a priori, resultando na distribuição de
probabilidade a posteriori.
O Teorema de Bayes é a regra de atualização da informação que se tem
sobre o parâmetro e é utilizado para quantificar esse aumento de informação, sendo
um elemento essencial na análise bayesiana, pois toda inferência é feita a partir da
distribuição a posteriori (Gamerman, 1997). Se p(θ) é a densidade a priori de θ,
então a densidade a posteriori de θ, p(θ|Y), é dada pelo Teorema de Bayes:
p(θ|Y ) =
p(Y |θ)p(θ)
p(Y |θ)p(θ)
(2.16)
em que p(Y |θ) é a função de verossimilhança e Y = (y
1
, y
2
, . . . , y
n
) é o vetor de
dados. O denominador funciona como uma constante normalizadora, que não
depende de θ. Assim, o teorema pode ser reescrito como:
p(θ|Y ) p(y|θ)p(θ)
em que representa proporcionalidade.
Em outras palavras, p(θ|Y ) é:
Posteriori Verossimilhança × Priori
17
No caso de θ ser multivariado, θ = (θ
1
, θ
2
, . . . , θ
p
), as distribuições mar-
ginais de cada componente θ
i
, a partir das quais as inferências para cada parâme-
tro serão feitas, devem ser obtidas por meio de integração da densidade conjunta
a posteriori p(θ
1
, θ
2
, . . . , θ
p
|y
1
, . . . , y
n
), com relação a cada parâmetro. No en-
tanto, a resolução desta integral, analiticamente, é, em muitos casos, impraticável.
Uma das alternativas neste caso é a utilização dos métodos aproximados de infe-
rência, baseados em simulação estocástica, os quais ganharam destaque entre os
pesquisadores a partir do trabalho de Gelfand & Smith (1990), utilizando algorit-
mos dos métodos Monte Carlo via cadeias de Markov (MCMC, do inglês Markov
Chain Monte Carlo).
Os métodos de simulação estocástica consideram as distribuições condi-
cionais completas a posteriori de cada parâmetro para gerar amostras que con-
vergem para a densidade marginal, com o aumento do tamanho dessa amostra
(Gelfand, 2000).
A distribuição condicional completa do parâmetro θ
i
, denotada por
p(θ
i
|θ
i
, Y ), é obtida considerando que, na densidade conjunta, os demais pa-
râmetros θ
i
são conhecidos e, assim, a expressão se torna menos complexa,
que as constantes podem ser desconsideradas.
Quando a expressão da condicional completa tem a forma de uma densi-
dade conhecida, um método de simulação indicado é o Amostrador de Gibbs, um
processo iterativo que gera valores que convergem para a distribuição a posteriori
conjunta, sem que se conheça a sua expressão (Gelfand, 2000). Se a distribuição
condicional completa a posteriori não tem a forma de uma densidade conhecida,
um dos métodos de simulação indicado é o algoritmo de Metropolis-Hastings
(Chib & Greenberg, 1995; Hastings, 1970;). Os algoritmos do amostrador de
Gibbs e do Metropolis-Hastings estão descritos a seguir.
18
2.5.1 Método de Monte Carlo via cadeias de Markov
Os métodos de Monte Carlo via cadeias de Markov (MCMC) são uma
alternativa aos métodos não iterativos para a resolução de integrais complexas
(Ehlers, 2007). Uma cadeia de Markov é um processo estocástico. Os dados
são gerados a partir de uma distribuição de probabilidade em que o valor gerado
em um estado atual depende do valor gerado no passo anterior, caracterizando,
assim, dependência entre os valores gerados.
Sendo assim, pode-se dizer que a essência do método MCMC é simular um
passeio aleatório no espaço do parâmetro θ, o qual converge para uma distribuição
estacionária, que é a distribuição a posteriori de interesse p(θ|Y ) (Gelman et al.,
2000). Os métodos MCMC mais utilizados são o amostrador de Gibbs e o algo-
ritmo de Metropolis-Hastings.
2.5.1.1 Amostrador de Gibbs
O amostrador de Gibbs é um esquema iterativo de amostragem de uma
cadeia de Markov, cujo núcleo de transição é formado pelas distribuições condi-
cionais completas (Gamerman, 1997). Esse método, inicialmente utilizado por
Geman & Geman (1984), passou a ser muito discutido e utilizado desde então
(Casella & George, 1992; Gelfand, 2000; Gelfand & Smith, 1990).
Para descrever o algoritmo, suponha-se que a distribuição de interesse
p(θ|Y), em que θ = (θ
1
, θ
2
, . . . , θ
p
|y
1
, . . . , y
n
), corresponde à distribuição a poste-
riori. Considere-se, ainda, que as densidades condicionais completas a posteriori
p(θ
i
|θ
i
, Y), i = 1, ..., p estão disponíveis. O algoritmo do amostrador de Gibbs é
descrito da seguinte forma:
1 - iniciar o contador de iterações da cadeia t = 1 e estabelecer valores iniciais
θ
0
= (θ
0
1
, θ
0
2
, . . . , θ
0
p
) ;
19
2 - obter um novo valor θ
t
= (θ
t
1
, θ
t
2
, . . . , θ
t
p
) a partir de θ
t1
por meio de sucessi-
vas gerações de valores:
θ
(t)
1
p(θ
1
|θ
(t1)
2
, θ
(t1)
3
. . . , θ
(t1)
p
)
θ
(t)
2
p(θ
2
|θ
(t)
1
, θ
(t1)
3
, . . . , θ
(t1)
p
)
.
.
.
θ
(t)
p
p(θ
p
|θ
(t)
1
, θ
(t)
2
. . . , θ
(t)
p1
);
3 - muda-se o contador t para t + 1 e volta-se ao passo 2 até a convergência ser
alcançada.
À medida que o número t de iterações aumenta, a seqüência de valores
gerados se aproxima da distribuição de equilíbrio, ou seja, da densidade marginal
desejada para cada parâmetro, quando se assume que a convergência foi atingida
(Gamerman, 1997).
2.5.1.2 Algoritmo de Metropolis-Hastings
A partir dos trabalhos de Metropolis et al. (1953) e Hastings (1970),
definiu-se um método MCMC para simular distribuições multivariadas, denomi-
nado Metropolis-Hastings (Chib & Greenberg, 1995). O algoritmo, e justificativas
teóricas para seu uso, podem ser encontradas com maiores detalhes em Gamerman,
(1997), e é apresentado a seguir.
Suponha-se que o objetivo seja gerar valores de uma distribuição de inte-
resse π(θ), cuja densidade é desconhecida. O método consiste em gerar valores
candidatos θ
de uma densidade auxiliar q(θ, θ
) que possa ser amostrada e rejeitar
ou aceitar esses valores com probabilidade α(θ, θ
). A expressão mais comumente
20
usada para a probabilidade de aceitação é:
α(θ, θ
) = min
1,
π(θ
)q(θ|θ
)
π(θ)q(θ
|θ)
(2.17)
Em termos práticos, o algoritmo de Metropolis-Hastings pode ser descrito
da seguinte forma:
1 - iniciar o contador de iterações t = 0 e especifique um valor inicial θ
t
;
2 - gerar um valor θ
da distribuição q(θ
t
, .); e um valor u de uma uniforme (0,1);
3 - calcular α(θ
t
, θ
);
4 - se u α aceitar o novo valor e fazer θ
t+1
= θ
;
5 - incrementar o contador t para t + 1 e voltar ao passo 2
O método define uma cadeia de Markov, pois as transições dependem ape-
nas das posições no estágio anterior. O núcleo de transição q define apenas uma
proposta de movimento que pode ou não ser confirmada pela probabilidade de
aceitação α.
Segundo Chib & Greenberg (1995), várias são as famílias de densidades
que podem ser escolhidas. Para densidade candidata q(θ, θ
), algumas das quais
foram citadas por Hastings (1970). Uma escolha bastante eficiente, quando
disponível, é explorar a forma de π(θ) para especificar a densidade candidata. Por
exemplo, se π(θ) pode ser escrita como π(θ) = ψ(θ)h(θ), sendo h(θ) possível de
ser amostrada e ψ(θ) é uniformemente limitada, pode-se fazer q(θ, θ
)=h(θ) para
gerar as amostras candidatas. Neste caso, a probabilidade de movimento se reduz
à seguinte equação:
α(θ, θ
) = min
1,
ψ
(θ)
ψ(θ)
(2.18)
21
2.5.2 Diagnóstico de convergência
À medida que o número de iterações aumenta, a cadeia gerada se apro-
xima da distribuição de equilíbrio, ou seja, da densidade marginal desejada de
cada parâmetro. Assim, assume-se que a convergência é atingida em uma iteração
cuja distribuição esteja arbitrariamente próxima da distribuição de equilíbrio. A
grande dificuldade é determinar quantas iterações são necessárias para se atingir
esta condição. No intuito de alcançar a convergência minimizando o tempo e o
esforço computacional, vários métodos estão disponíveis para monitorar a con-
vergência das seqüências geradas pelas cadeias de Markov, entre eles o critério
proposto por Raftery & Lewis, (1992) e o critério proposto por Heidelberger &
Welch (1983), citados por Nogueira et al. (2004).
2.5.2.1 Raftery & Lewis (1992)
O diagnóstico de Raftery & Lewis (1992) baseia-se em obter uma amostra
piloto estimativas do número de iterações necessárias para se obter a convergência,
do número de iterações iniciais que devem ser descartadas (burn-in) e da distância
mínima (k); k é o fator de dependência de uma iteração à outra para se obter uma
amostra independente. Segundo os autores, se o fator de dependência, for maior
que 5, pode-se concluir que a cadeia apresenta problemas de convergência.
2.5.2.2 Heidelberger & Welch (1983)
Para testar a hipótese nula da estacionaridade de uma cadeia gerada se-
gundo o critério de Heidelberger & Welch (1983), dois testes simultâneos são
utilizados. Se a hipótese nula é rejeitada para um dado valor, o teste é repetido
depois de descartados os 10% valores iniciais da seqüência. Se a hipótese é nova-
mente rejeitada, mais 10% dos valores iniciais são descartados e assim sucessiva-
22
mente até serem descartados os 50% valores iniciais. Se a hipótese for novamente
rejeitada, isso indica que é necessário um número maior de iterações. Caso con-
trário, o número inicial de iterações descartadas é indicado como o tamanho do
"burn-in". Considerando-se que uma sequência de valores amostrados passou no
teste de estacionaridade, aplica-se, nesta parte, o teste de Half-Width para verificar
se a média estimada está sendo calculada com uma acurácia pré-especificada. A
convergência fica garantida quando se aceita a hipótese nula nos dois testes.
2.5.3 Fator de Bayes
Sejam M
i
o modelo 1 e M
j
o modelo 2. Um método que pode ser utilizado
para fazer comparações entre modelos é o Fator de Bayes (FB) e pode ser definido,
segundo Kass & Raftery (1995), como:
F B
ij
=
P (Y |M
i
)
P (Y |M
j
)
=
L(Y |θ
i
, M
i
)P (θ
i
|M
i
)
i
L(Y |θ
j
, M
j
)P (θ
j
|M
j
)
j
(2.19)
em que: L(Y |θ
i
, M
i
) e L(Y |θ
j
, M
i
) são a função de verossimilhança dos mo-
delos M
i
e M
j
respectivamente e P (θ
i
|M
i
) e P (θ
j
|M
j
) são, respectivamente, a
distribuição a priori correspondente ao modelo M
i
e M
j
.
Se o fator de Bayes (F B
ij
) for maior que 1, pode-se dizer que o modelo
do numerador (M
i
) é mais plausível. Caso contrário, escolhe-se o modelo do
denominador (M
j
).
A solução da equação (2.19), geralmente não é analítica. Um método uti-
lizado para obter uma aproximação da marginal dos dados, é a amostragem por
importância.
Este processo permite utilizar como aproximação da marginal dos dados, a
estimativa da média (Raftery, 1995) ou da média harmôrnica (Raftery & Newton,
2007), quando a função de importância são respectivamente a distribuição a priori
23
ou a distribuição a posteriori dos parâmetros dos modelos. Segundo Raftery &
Newton (2007), se a função de importância for a distribuição a posteriori, uma
aproximação para P (Y |M
i
) é dada por:
P (Y |M
i
) =
1
B
B
t=1
1
p(Y |θ
t
)
1
(2.20)
em que t= 1, 2, . . . , B representa cada iteração dos algoritmos MCMC.
O estimador da distribuição na equação (2.20) é a média harmônica dos
valores da pseudo-verossimilhança calculados para a amostra da distribuição a
posteriori obtida pelo algoritmo de Metropolis-Hastings
Uma interpretação do fator de Bayes sugerido por Kass & Raftery (1995),
é apresentada na Tabela 1.
TABELA 1 Interpretação do Fator de Bayes
FB 2log(FB) Evidência a favor de M
i
< 1 < 0 Negativa (apoia M
j
)
1 3 0 2 Fraca
3 20 2 6 Positiva
20 150 6 10 Forte
> 150 > 10 Muito Forte
3 MATERIAL E MÉTODOS
3.1 Material experimental
Foram analisados neste trabalho, vários possíveis resultados experimentais
de diluições seriadas. Utilizaram-se experimentos com cinco e três diluições, e,
para cada uma dessas diluições, foram considerados três tubos. As quantidades
24
inoculadas das amostras para o caso de três diluições foram 0,1; 0,01 e 0,001 e no
caso de 5 diluições foram 0,1; 0,01; 0,001; 0,0001 e 0,00001; .
Os resultados experimentais utilizadas no ajuste dos modelos foram (00330),
(01300), (03320), (32300) e (33300) para 5 diluições e (030),(032) e (330) para 3
diluições.
3.2 Metodologia
3.2.1 Modelo Completo
Nesta seção, apresenta-se o modelo Completo como proposta de modela-
gem do número mais provável deste trabalho. Este modelo é uma combinação dos
modelos Usual, de Toxidez e Interferência apresentados na seção 2.4. Os parâ-
metros deste modelo são λ, λ
i
e T
j
. Os termos utilizados neste modelo foram
descritos nos demais modelos. As suposições associadas a este modelo são:
1. os tubos são independentes;
2. os microrganismos devem estar aleatoriamente distribuídos nos tubos e su-
põem-se que duas espécies sejam distribuídas independentemente. O número
de microrganismo no tubo para as duas espécies segue uma distribuição de
Poisson;
3. a substância tóxica é uniformemente dispersa no caldo;
4. o efeito da substância tóxica em cada microrganismo no tubo é indepen-
dente;
5. o tubo mostrará crescimento de microrganismos se não ocorrer espécies que
interferem nem substância tóxica presente no tubo;
25
A função de verossimilhança do modelo é dada por:
L(Y |λ, λ
i
, T
j
) =
k
j=1
n
j
y
j
[e
λ
i
z
j
(1T
j
)
(1 e
λz
j
(1T
j
)
]
y
j
× [1 e
λ
i
z
j
(1T
j
)
(1 e
λz
j
(1T
j
)]
(n
j
y
j
)
(3.1)
Verifica-se que, quando λ
i
=0, a equação (3.1) se reduz ao modelo de To-
xidez dado pela equação (2.13). Se T
j
= 0, a equação (3.1) se reduz ao modelo
de Interferência dado pela equação (2.15). Se λ
i
=0 e T
j
= 0, a equação (3.1) se
reduz ao modelo Usual dado pela equação (2.5).
Considerando-se os modelos apresentados em na seção 2.4 e os resulta-
dos experimentais, as inferências sobre os parâmetros dos modelos foram feitas
utilizando-se de técnicas de inferência bayesiana.
3.3 Distre1nicas d1sicas(tˇa955 Tf 0 -36.876 Td[e88 6i-.2f 0 -36.876 Td[(3.3)-1000(Dis1.636 Td[(=)]TJ/F157Td[(i)-250(modelo)-250(U2-20.324 Td[(u20.324 Td[(reduz)-250(ando-s.TJ/3)]TJ/F15 10.90uufu82o.(3.3s)-346(=U2-20.324 Td[(20 Td[-Q4 Tn875dosd[(j)]2 i26)2.5)8324020se)-3 Tf 3.903 1.215a(V)02i26)2.5)8324020i7li(bayesi.v39 583.229 Td[(%e)-25e)-250(de)0i5 3.229 Td[(%e)-25e)-250(de) se)-3 Tf 0(d4ni26)2.5;o)-346(2.4)tcz.-341(in -20.323 Td[(dos93 -1.215 Td[(j)].5;o)-346(.03niaDistre1nical02.15).as232;o)-346(2073o)-346(.03niaDis2e1nicas)-25e)-250(d1sicas(tˇa955 Tf 0 -36.876 Td[e88 6i-.2f 0 -36.876 Td[(3.3)-1000(Dis1.636 Td[(=)]TJ/F157Td[(i)-250(modelo)-250(U2-231sica.97 Tf 37.694 13.636 d[(%e)-25e)-250(de6.(u20.324 Td[(r Tf 4.18.d[(3.3)-1007 .j481p88)65020.324 Td[(r Tf 3)-10005n7o25e)-23x4)-295(se)]TJ -270.57 -20.324 T93m tTf 6.364 -1.637 se)-3 Tf 0 Tf 0 -36.87Tf 0=3.3)-141.(u)86istr[(24is9501)]TJ/F3aDis2e1nicas)-25e2icas(t(equação)-256us.tr[(24is9501)]i
em que α
j
= 0 e β
j
= 1
3.3.1 Processo de geração das cadeias de Markov
As cadeias dos parâmetros de cada modelo foram geradas considerando
cinco e três diluições e três tubos para cada diluição.
Para determinar a quantidade de iterações iniciais descartadas (descarte) e
o tamanho do salto entre cada iteração (salto), foi feita uma amostra piloto e calcu-
lado o diagnóstico de Raftery & Lewis (1992). As configurações apresentadas têm
valores de descarte e salto muito além dos sugeridos pelo diagnóstico de Raftery
& Lewis (1992). Isso foi feito para garantir que todos os parâmetros dos modelos
atingissem a convergência.
No processo de geração das cadeias foi obtida uma amostra final de tama-
nho 5000 considerando-se um descarte inicial de 50000 valores e um salto a cada
200 iterações.
Atribuiu-se, para os parâmetros λ e λ
i
, distribuições a priori gama com
parâmetros (2,100).
Para o parâmetro T
j
, utilizou-se como distribuição a priori, distribuições
uniformes. A quantidade de parâmetros T, está associada a quantidade de diluições
existentes no experimento. Sendo assim estimou-se uma probabilidade de Toxidez
para cada diluição do experimento estudado.
A distribuição utilizada como geradora de valores candidatos para os parâ-
metros λ e λ
i
é a gama com parâmetros (α
c
, β
c
) e (α
c
i
, β
c
i
) para os respectivos pa-
râmetros. Os parâmetros desta distribuição foram continuamente atualizados com
uma amostra dos últimos 1000 valores da cadeia da distribuição condicional com-
pleta a posteriori. Este processo consistiu em calcular valores de α e β, segundo
as definições de média e variância da distribuição gama e são respectivamente:
27
E(Y ) = αβ e V (Y ) = αβ
2
.
A distribuição geradora de valores candidatos para T foi a uniforme, e o
parâmetro desta distribuição para a primeira diluição é (0,1). Os demais parâme-
tros estarão entre zero e o valor da probabilidade de Toxidez da diluição anterior.
Estas distribuições podem ser definidas da seguinte maneira:
T
c
1
U[0, 1]
T
c
2
U[0, T
c
1
]
.
.
.
T
c
j
U[0, T
c
j1
]
A probabilidade conjunta destas probabilidades é dada por:
P (T
c
1
, T 2
c
, . . . , T j
c
) = 1 ×
1
T
c
1
×
1
T
c
2
× . . . ×
1
T
c
j1
=
k
j=1
1
T
c
j1
É importante ressaltar que, embora tenham sido considerados, apenas re-
sultados experimentais envolvendo cinco e três diluições com três tubos para cada
diluição, é possível obter estimativas para outras combinações de diluições e tubos.
Para cada situação experimental considerada neste trabalho, foram ajusta-
dos os quatro modelos estudados.
A implementação do algoritmo de Metropolis-Hastings considerando as
configurações apresentadas, foi feita no pacote estatístico R (R Development Core
Team, 2006) e podem ser encontradas no Anexo B. A probabilidade de aceitação
utilizada em cada caso é definida pela equação (2.17) e estão devidamente especi-
ficadas na seção 4.3.
A avaliação de convergência das cadeias geradas para cada parâmetro foi
28
feita por meio dos diagnósticos de Raftery & Lewis (1992) e Heidelberger &
Welch (1983) para tanto foram usados os pacotes BOA ("Bayesyan Output Analy-
sis") e CODA ("Convergence Diagnostics and Output Analysis") do software R. A
comparação entre os modelos foi feito por meio do Fator de Bayes.
4 RESULTADOS E DISCUSSÃO
Nesta seção foram apresentadas as obtenções das distribuições a posteriori con-
junta dos modelos amostrais, suas respectivas condicionais completas, a proba-
bilidade de aceitação no algoritmo de Metropolis-Hastings e alguns exemplos de
aplicação.
4.1 Resultados metodológicos: obtenção das distribuições a posterio-
ris conjunta para os modelos
Nesta seção é apresentada a obtenção da distribuição conjunta a posteriori
para cada modelo.
As distribuições conjuntas a posteriori de todos os parâmetros dos mode-
los, dado que uma amostra aleatória foi realizada, foram obtidas a partir do Teo-
rema de Bayes, multiplicando-se a verossimilhança e todas as prioris correspon-
dentes de cada modelo. Portanto, as distribuições conjuntas a posteriori, supondo
independência entre os parâmetros, foram apresentadas nesta seção.
4.1.1 Distribuição a posteriori conjunta para o modelo Usual
A distribuição conjunta a posteriori do modelo Usual pode ser obtida
por meio da multiplicação da função de verossimilhança do modelo, definida na
29
equação (2.5), pela distribuição a priori definida na equação (3.2).
Sendo assim, pode-se escrever distribuição a posteriori conjunta do mo-
delo Usual como:
P (λ|Y ) L(Y |λ)P (λ)
P (λ|Y )
k
j=1
n
j
y
j
(1 e
λz
j
)
y
j
e
λz
j
(n
j
y
j
)
× λ
α1
e
λ
β
(4.1)
4.1.2 Distribuição a posteriori conjunta para o modelo de Toxidez
Neste modelo, têm-se dois parâmetros a serem estimados, o parâmetro λ e
o parâmetro T
j
, especificados na seção 2.4.2. A distribuição a priori para λ está
definido na equação (3.2) e do parâmetro T
j
, esta definida na equação (3.4).
Sendo assim, a distribuição posteriori conjunta para este modelo é definida
como:
P (λ, T
j
|Y ) L(Y |λ, T
j
)P (λ)P (T
j
)
P (λ, T
j
|Y )
k
j=1
n
j
y
j
(1 e
λz
j
(1T
j
)
)
y
j
× e
λz
j
(1T
j
)
(n
j
y
j
)
× λ
α1
e
λ
β
(4.2)
4.1.3 Distribuição a posteriori conjunta para o modelo de Interferência
A obtenção da distribuição a posteriori conjunta do modelo de Interferên-
cia é análoga aos modelos anteriores. Neste caso, também existem dois parâmetros
a serem estimados, λ e λ
i
.
A distribuição a priori para λ esta definida na equação (3.2), a dis-
tribução a priori para λ
i
esta definida na equação (3.3)
30
A distribuição a posteriori conjunta é, então, dada pela equação:
P (λ, λ
i
|Y ) L(Y |λ, λ
i
)P (λ)P (λ
i
)
P (λ, λ
i
|Y )
k
j=1
n
j
y
j
[e
λ
i
z
j
(1 e
λz
j
)]
y
j
× [1 e
λ
i
z
j
(1 e
λz
j
)]
(n
j
y
j
)
× λ
α1
e
λ
β
× λ
α
i
1
i
e
β
i
λ
i
(4.3)
4.1.4 Distribuição a posteriori conjunta para o modelo Completo
Da mesma forma que nos modelos anteriores a distribuição a posteriori
conjunta para o modelo Completo é definida. As distribuições a priori dos parâ-
metros estão definidas nas equações (3.2), (3.3) e (3.4).
A distribuição a posteriori conjunta para este modelo é:
P (λ, λ
i
, T
j
|Y ) L(Y |λ, λ
i
, T
j
)P (λ)P (λ
i
)P (T
j
)
P (λ, λ
i
, T
j
|Y )
k
j=1
n
j
y
j
[e
λ
i
z
j
(1T
j
)
(1 e
λz
j
(1T
j
)
)]
y
j
× [1 e
λ
i
z
j
(1T
j
)
(1 e
λz
j
(1T
j
)
)]
(n
j
y
j
)
× λ
α1
e
β
λ
× λ
α
i
1
i
e
β
i
λ
i
(4.4)
4.2 Distribuições condicionais completas a posteriori
Nesta seção é apresentada a obtenção da distribuição condicional completa
a posteriori para cada parâmetro dos modelos.
31
P (T
j
|λ, Y ) L(Y |λ, T
j
)P (λ)P (T
j
)
P (T
j
|λ, Y )
k
j=1
n
j
y
j
(1 e
λz
j
(1T
j
)
)
y
j
× e
λz
j
(1T
j
)
(n
j
y
j
)
(4.7)
4.2.3 Obtenção da distribuição condicional completa a posteriori para os
parâmetros do modelo de Interferência
Por meio da equação (4.3), derivam-se as distribuições condicionais com-
pletas a posteriori para os parâmetros do modelo de interferência, λ e λ
i
.
A distribuição condicional completa a posteriori para o parâmetro λ, é
dada por:
P (λ|λ
i
, Y ) L(Y |λ, λ
i
)P (λ)P (λ
i
)
P (λ|λ
i
, Y )
k
j=1
n
j
y
j
[e
λ
i
z
j
(1 e
λz
j
)]
y
j
× [1 e
λ
i
z
j
(1 e
λz
j
)]
(n
j
y
j
)
× λ
α1
e
β
λ
(4.8)
A distribuição condicional completa a posteriori para o parâmetro λ
i
é
definida por:
P (λ
i
|λ, Y ) L(Y |λ, λ
i
)P (λ)P (λ
i
)
P (λ
i
|λ, Y )
k
j=1
n
j
y
j
[e
λ
i
z
j
(1 e
λz
j
)]
y
j
× [1 e
λ
i
z
j
(1 e
λz
j
)]
(n
j
y
j
)
× λ
α
i
1
i
e
β
i
λ
i
(4.9)
33
4.2.4 Obtenção da distribuição condicional completa a posteriori para os pa-
râmetros do modelo Completo
As distribuições condicionais completas a posteriori dos parâmetros do
modelo proposto são obtidas por meio da equação (4.4).
A distribuição condicional completa a posteriori para o parâmetro λ, é
dada por:
P (λ|λ
i
, T
j
, Y ) L(Y |λ, λ
i
, T
j
)P (λ)P (λ
i
)P (T
j
)
P (λ|λ
i
, T
j
, Y )
k
j=1
n
j
y
j
[e
λ
i
z
j
(1T
j
)
(1 e
λz
j
(1T
j
)
)]
y
j
× [1 e
λ
i
z
j
(1T
j
)
(1 e
λz
j
(1T
j
)
)]
(n
j
y
j
)
× λ
α1
e
β
λ
(4.10)
Define-se a distribuição condicional completa a posteriori do parâmetro
λ
i
, por:
P (λ
i
|λ, T
j
, Y ) L(Y |λ, λ
i
, T
j
)P (λ)P (λ
i
)P (T
j
)
P (λ
i
|λ, T
j
, Y )
k
j=1
n
j
y
j
[e
λ
i
z
j
(1T
j
)
(1 e
λz
j
(1T
j
)
)]
y
j
× [1 e
λ
i
z
j
(1T
j
)
(1 e
λz
j
(1T
j
)
)]
(n
j
y
j
)
× λ
α
i
1
i
e
β
i
λ
i
(4.11)
A equação a seguir é a distribuição condicional completa a posteriori para
34
o parâmetro T
j
:
P (T
j
|λ, λ
i
, Y ) L(Y |λ, λ
i
, T
j
)P (λ)P (λ
i
)P (T
j
)
P (T
j
|λ, λ
i
, Y )
k
j=1
n
j
y
j
[e
λ
i
z
j
(1T
j
)
(1 e
λz
j
(1T
j
)
)]
y
j
× [1 e
λ
i
z
j
(1T
j
)
(1 e
λz
j
(1T
j
)
)]
(n
j
y
j
)
(4.12)
As distribuições condicionais completas a posteriori para todos os parâ-
metros dos modelos não possuem formas de densidades conhecidas. Sendo assim,
é necessário utilizar o algoritmo de Metropolis-Hastings para se amostrar destas
distribuições.
4.3 Probabilidade de aceitação no algoritmo de Metropolis-Hastings
Nesta seção, será apresentada a razão R da probabilidade de aceitação do
algortimo de Metropolis-Hastings, para os parâmetros dos modelos.
Para cada parâmetro dos modelos estudados, será feita esta razão con-
siderando a distribuição gama para o parâmetro λ, gama para o parâmetro λ
i
e
uniforme para o parâmetro T
j
. A distribuição geradora de valores candidatos uti-
lizada para os parâmetros λ, λ
i
e T
j
é dinâmica e a distribuição a priori destes
parâmetros é estática.
Em todos os modelos, a obtenção do valor R será feita de forma análoga.
4.3.1 Probabilidade de aceitação no algoritmo de Metropolis-Hastings para
o modelo Usual
No caso do modelo Usual apenas um parâmetro a ser estimado, sendo
assim tem-se apenas R
λ
U(razão na probabilidade de aceitação para o parâmetro λ
no modelo Usual). Para se obter a razão da probalilidade de aceitação utiliza-se a
35
distribuição condicional completa para o parâmetro λ do modelo Usual dada pela
equação (4.5) multiplicado pela densidade candidata gama (α, β). Sendo assim
R
λ
U pode ser definida por:
R
λ
U =
k
j=1
n
j
y
j
(1 e
λ
z
j
)
y
j
e
λ
z
j
(n
j
y
j
)
k
j=1
n
j
y
j
(1 e
λz
j
)
y
j
e
λz
j
(n
j
y
j
)
×
(λ
)
α1
e
β
λ
× λ
α
c
1
e
β
c
λ
λ
α1
e
βλ
× (λ
)
α
c
1
e
β
c
λ
(4.13)
4.3.2 Probabilidade de aceitação no algoritmo de Metropolis-Hastings para
o modelo de Toxidez
De forma análoga, define-se R
λ
T (razão na probabilidade de aceitação para
o parâmetro λ no modelo de Toxidez) e R
T
T (razão na probabilidade de aceitação
para o parâmetro T
j
no modelo de Toxidez). Neste caso, utiliza-se a distribuição
condicional completa para os parâmetros λ e T
j
, dados pela equação (4.6) e (4.7),
respectivamente.
A razão na probabilidade de aceitação para o parâmetro λ pode ser descrita
da seguinte forma:
R
λ
T
j
=
k
j=1
n
j
y
j
(1 e
λ
z
j
(1T
j
)
)
y
j
e
λ
z
j
(1T
j
)
(n
j
y
j
)
k
j=1
n
j
y
j
(1 e
λz
j
(1T
j
)
)
y
j
e
λz
j
(1T
j
)
(n
j
y
j
)
×
(λ
)
α1
e
β
λ
× λ
α
c
1
e
β
c
λ
λ
α1
e
βλ
× (λ
)
α
c
1
e
β
c
λ
(4.14)
36
A razão para o parâmetro T
j
é feita da mesma maneira e é definida por:
R
T
T =
k
j=1
n
j
y
j
(1 e
λz
j
(1T
j
)
)
y
j
e
λz
j
(1T
j
)
(n
j
y
j
)
k
j=1
n
j
y
j
(1 e
λz
j
(1T
j
)
)
y
j
e
λz
j
(1T
j
)
(n
j
y
j
)
×
k
j=1
1
T
c
j1
k
j=1
1
(T
c
j1
)
(4.15)
4.3.3 Probabilidade de aceitação no algoritmo de Metropolis-Hastings para
o modelo de Interferência
A obtenção de R no modelo de Interferência é feita de forma análoga aos
modelos Usual e de Toxidez. Neste caso, determina-se R
λ
I(razão da probabili-
dade de aceitação para o parâmetro λ no modelo de Interferência) e R
λ
i
I (razão
da probabilidade de aceitação para o parâmetro λ
i
no modelo de Interferência).
As distribuições condicionais completas para os parâmetros λ e λ
i
estão represen-
tadas nas equações (4.8) e (4.9), respectivamente. As disribuições candidatas para
os parâmetros λ e λ
i
, são respectivamente gama (α,β) e gama (α
i
,β
i
)
Define-se a razão de aceitação para o parâmetro λ deste modelo como:
R
λ
I =
k
j=1
n
j
y
j
[e
λ
i
z
j
(1 e
λ
z
j
)]
y
j
[1 e
λ
i
z
j
(1 e
λ
z
j
)]
(n
j
y
j
)
k
j=1
[e
λ
i
z
j
(1 e
λz
j
)]
y
j
[1 e
λ
i
z
j
(1 e
λz
j
)]
(n
j
y
j
)
(4.16)
×
(λ
)
α1
e
β
λ
× λ
α
c
1
e
β
c
λ
λ
α1
e
β
λ
× (λ
)
α
c
1
e
β
c
λ
A razão na probabilidade de aceitação para o parâmetro λ
i
, ou seja, R
λ
i
I,
37
é definida da seguinte forma:
R
λ
i
I =
k
j=1
n
j
y
j
[e
(λ
i
)
z
j
(1 e
λz
j
)]
y
j
[1 e
(λ
i
)
z
j
(1 e
λz
j
)]
(n
j
y
j
)
k
j=1
n
j
y
j
[e
λ
i
z
j
(1 e
λz
j
)]
y
j
[1 e
λ
i
z
j
(1 e
λz
j
)]
(n
j
y
j
)
×
(λ
i
)
α
i
1
e
βλ
i
× λ
α
c
i
1
i
e
β
c
λ
i
λ
α
i
1
i
e
β
λ
i
× (λ
i
)
α
c
i
1
e
β
c
λ
i
(4.17)
4.3.4 Probabilidade de aceitação no algoritmo de Metropolis-Hastings para
o modelo Completo
Finalmente, será definida a razão da probabilidade de aceitação para os
parâmetros do modelo Completo. Estas razões são definidas da mesma maneira
que as dos modelos anteriores. Neste caso, serão definidas as razões das probabili-
dades de aceitação para os parâmetros λ, λ
i
e T
j
, dadas respectivamente por: R
λ
C,
R
λ
i
C e R
T
C. A definição da nomenclatura das razões foi apresentada nos mo-
delos anteriores. Ressalta-se apenas que a diferença está simplesmente relacionada
ao modelo o qual a razão se refere, e cada uma delas diferencia-se pela utilização
da letra inicial maiúscula de cada modelo.
Para o parâmetro λ, a razão da probabilidade de aceitação é dada por:
R
λ
T
j
=
k
j=1
[e
λ
i
z
j
(1T
j
)
(1 e
λ
z
j
(1T
j
)
)]
y
j
k
j=1
[e
λ
i
z
j
(1T
j
)
(1 e
λz
j
(1T
j
)
)]
y
j
×
[1 e
λ
i
z
j
(1T
j
)
(1 e
λ
z
j
(1T
j
)
)]
(n
j
y
j
)
[1 e
λ
i
z
j
(1T
j
)
(1 e
λz
j
(1T
j
)
)]
(n
j
y
j
)
×
(λ
)
α1
e
β
λ
× λ
α
c
1
e
β
c
λ
λ
α1
e
β
λ
× (λ
)
α
c
1
e
β
c
λ
(4.18)
38
A razão da probabilidade de aceitação para o parâmetro λ
i
é definida por:
R
λ
i
C =
k
j=1
[e
(λ
i
)
z
j
(1T
j
)
(1 e
λz
j
(1T
j
)
)]
y
j
k
j=1
[e
λ
i
z
j
(1T
j
)
(1 e
λz
j
(1T
j
)
)]
y
j
×
[1 e
(λ
i
)
z
j
(1T
j
)
(1 e
λz
j
(1T
j
)
)]
(n
j
y
j
)
[1 e
λ
i
z
j
(1T
j
)
(1 e
λz
j
(1T
j
)
)]
(n
j
y
j
)
×
(λ
i
)
α
i
1
e
βλ
i
× λ
α
c
i
1
i
e
β
c
λ
i
λ
α
i
1
i
e
β
λ
i
× (λ
i
)
α
c
i
1
e
β
c
λ
i
(4.19)
Para o parâmetro T
j
, a razão R
T
C é descrita da seguinte forma:
R
T
C =
k
j=1
[e
λ
i
z
j
(1T
j
)
(1 e
λz
j
(1T
j
)
)]
y
j
k
j=1
[e
λ
i
z
j
(1T
j
)
(1 e
λz
j
(1T
j
)
)]
y
j
×
[1 e
λ
i
z
j
(1T
j
)
(1 e
λz
j
(1T
j
)
)]
(n
j
y
j
)
[1 e
λ
i
z
j
(1T
j
)
(1 e
λz
j
(1T
j
)
)]
(n
j
y
j
)
×
k
j=1
1
T
c
j1
k
j=1
1
(T
c
j1
)
(4.20)
4.4 Exemplos de aplicação, considerando cinco diluições
São apresentadas, nesta seção, as comparações entre os ajustes os modelos
para alguns resultados experimentais selecionados. Para esta seção, foram ilustra-
dos exemplos de aplicações do método utilizando-se cinco diluições. Em seções
subsequentes são apresentados resultados com três diluições. Outros resultados
estão no Anexo A e seguem o mesmo padrão de interpretação do que está aqui
apresentado.
As tabelas foram organizadas em sequência, de forma que, primeiro, são
39
apresentados resultados referentes à convergência das cadeias dos parâmetros, de-
pois os resultados referentes às médias a posteriori, seus respectivos intervalos de
credibilidade (HPD; do inglês highest probability density) e erros de Monte Carlo.
Por último, é apresentada a comparação entre os modelos, por meio do Fator de
Bayes.
4.4.1 Convergência do processo de Markov e estimação do número mais
próvavel para resultado experimental (01300)
Em ensaios de diluição podem ocorrer resultados nos quais, em baixas
diluições, o crescimento de um microrganismo de interesse seja inibido. Fatores,
como Toxidez no meio de cultura, duas espécies que competem por nutrientes no
mesmo tubo ou outro mecanismo, podem ser as causas deste tipo de resultado ex-
perimental. Para ilustrar esta situação, foi utilizado o exemplo (01300). Neste caso
observa-se que não crescimento de microrganismo em nenhum tubo da primeira
diluição e, na segunda diluição, não crescimento em dois tubos. Este tipo de
resultado é pouco provável sob o modelo Usual. Outros resultados experimentais,
semelhantes a este estão apresentados no Anexo A e seguem o mesmo padrão de
interpretação do resultado (01300).
Os resultados obtidos pelos critérios de convergência das cadeias geradas
pelo método MCMC encontram-se nas Tabelas 2. Por meio dos resultados desta
tabela 2, verifica-se que a convergência de todos os parâmetros dos modelos foi
assegurada, segundo o critério de Heidelberg & Welch (1983). Verifica-se também
que apenas os parâmetros T
3
, T
4
e T
5
do modelo de Toxidez não atingiram con-
vergência por meio do critério de Raftery & Lewis (1992), pois os valores do fator
de dependência (FD) para estes parâmetros foram maiores que 5.
A monitoração informal da convergência das cadeias de cada parâmetro
40
TABELA 2 Convergência dos parâmetros dos modelos para o resultado experi-
mental (01300).
Raftery & Lewis Heidelberg & Welch
Modelo Parâmetros FD estacionária Half-width
Usual λ 0,9506 sim sim
Interferência λ 0,9666 sim sim
λ
i
0,9989 sim sim
Tóxico λ 1,0152 sim sim
T
1
2,5974 sim sim
T
2
3,0048 sim sim
T
3
18,5496 sim sim
T
4
30,4871 sim sim
T
5
12,8051 sim sim
Completo λ 0,9663 sim sim
λ
i
1,0320 sim sim
T
1
1,0152 sim sim
T
2
1,0250 sim sim
T
3
0,9914 sim sim
T
4
1,0082 sim sim
T
5
0,9506 sim sim
dos modelos Usual, de Interferência, de Toxidez e Completo foi realizada por
meio de técnicas gráficas e está apresentada nas Figuras 1 a 7.
O traço e a densidade a posteriori da cadeia do parâmetro λ para o modelo
Usual estão representados na Figura de 1. Verifica-se, nesta figura, que a cadeia
deste parâmetro convergiu. O mesmo se verifica para as cadeias dos parâmetros
do modelo de Interferência, representadas na da Figura 2. Ambos os parâmetros
têm convergência garantida, pois as cadeias estão num estado de equilíbrio.
Por meio da Figura 3, verifica-se que a cadeia do parâmetro λ no modelo
de Toxidez atingiu convergência. Verifica-se, na Figura 4, que convergência nas
cadeias dos parâmetros T
1
, T
2
, T
3
, T
4
e T
5
deste modelo.
41
(a) (b)
FIGURA 1 Traço da cadeia gerada (a) e densidade a posteriori (b) do parâmetro
λ, do modelo Usual, para o resultado experimental (01300).
(a) (b)
0 1000 2000 3000 4000 5000
0 100 200 300 400 500
Iterações
λλ
I
0 500 1000 1500
0.0000 0.0005 0.0010 0.0015 0.0020
Densidade
λλ
I
(c) (d)
FIGURA 2 Traço das cadeias geradas (a) e (c) e densidades a posteriori (b) e (d)
dos parâmetros λ, e λ
i
, do modelo de Interferência, para o resultado
experimental (01300).
(a) (b)
FIGURA 3 Traço da cadeia gerada (a) e densidade a posteriori (b) do parâmetro
λ, do modelo de Toxidez, para o resultado experimental (01300).
42
0 1000 2000 3000 4000 5000
0.80 0.85 0.90 0.95 1.00
Iterações
T
1
0.80 0.85 0.90 0.95 1.00
0 20 40 60 80
Densidade
T
1
(a) (b)
(c) (d)
(e) (f)
0 1000 2000 3000 4000 5000
0.0 0.2 0.4 0.6 0.8
Iterações
T
4
0.0 0.2 0.4 0.6 0.8
0.0 0.5 1.0 1.5 2.0 2.5 3.0 3.5
Densidade
T
4
(g) (h)
0 1000 2000 3000 4000 5000
0.0 0.1 0.2 0.3 0.4 0.5 0.6
Iterações
T
5
0.0 0.1 0.2 0.3 0.4 0.5 0.6
0 2 4 6 8
Densidade
T
5
(i) (j)
FIGURA 4 Traço das cadeias geradas (a), (c), (e), (g) e (i) e densidades a pos-
teriori (b), (d), (f), (h) e (j) dos parâmetros T
1
, T
2
, T
3
, T
4
e T
5
, do
modelo de Toxidez, para o resultado experimental (01300).
43
Observando-se a Figura 5, nota-se que as cadeia dos parâmetros λ e λ
i
do
modelo Completo atingiram convergência. Verifica-se, nas Figuras 6 e 7, que
convergência nas cadeias dos parâmetros T
1
, T
2
, T
3
, T
4
e T
5
deste modelo.
0 1000 2000 3000 4000 5000
0 500 1000 1500
Iterações
λλ
0 500 1000 1500
0.0000 0.0005 0.0010 0.0015 0.0020
Densidade
λλ
(a) (b)
(c) (d)
FIGURA 5 Traço das cadeias geradas (a) e (c) e densidades a posteriori (b), (d)
do parâmetro λ e λ
i
, respectivamente, para o modelo Completo, para
o resultado experimental (01300).
0 1000 2000 3000 4000 5000
0.0 0.2 0.4 0.6 0.8 1.0
Iterações
T
1
(a) (b)
FIGURA 6 Traço da cadeia gerada (a) e densidade a posteriori (b) do parâmetros
T
1
, do modelo Completo, para o resultado experimental (01300).
44
0 1000 2000 3000 4000 5000
0.0 0.2 0.4 0.6 0.8 1.0
Iterações
T
2
0.0 0.2 0.4 0.6 0.8 1.0
0.0 0.5 1.0 1.5 2.0 2.5
Densidade
T
2
(a) (b)
(c) (d)
(e) (f)
(g) (h)
FIGURA 7 Traço das cadeias geradas (a), (c), (e) e (g) e densidades a posteriori
(b), (d), (f), (h) e (j) dos parâmetros T
2
, T
3
, T
4
e T
5
do modelo Com-
pleto, para o resultado experimental (01300).
45
Os resultados das Tabela 3 referem-se às estimativas das médias a poste-
riori de cada parâmetro dos modelos, com seus respectivos intervalos de credibili-
dade e erro de Monte Carlo.
Verifica-se, por meio da Tabela 3, que a estimativa do NMP (média a pos-
teriori para o parâmetro (λ), do modelo Usual, está subestimada em relação ao
mesmo parâmetro nos demais modelos. É provável que, em resultados onde
inibição de crescimento de microrganismos em baixas diluições, exista algum tipo
de fator que possa estar influenciando a falta de crescimento. Dessa forma, o
modelo Usual se torna pouco provável sob este tipo de resultado experimental,
subestimando λ. Este fato fora constatado por Cochram (1950). Os valores mé-
dios para o parâmetro λ nos diferentes modelos foram diferentes. Destaca-se que
a maior estimativa foi obtida por meio do modelo de interferência. Foi possível
estimar, por meio dos modelos de Interferência e Completo, microrganismos que
interferem, (λ
i
), no crescimento do microrganismo em estudo (λ). De modo geral,
os intervalos de credibilidade dos parâmetros dos modelos são amplos e o erros de
Monte Carlo associados, pequenos.
46
TABELA 3 Estimativas das médias a posteriori e intervalos de credibilidade
(I.C.95%), dos parâmetros dos modelos e erro de Monte Carlo (Erro
M.C.) para o resultado (01300).
I.C. 95%
Modelos Parâmetros Média LI LS Erro M.C.
Usual λ 196,4000 7,1770 472,3158 2,0875
Interferência λ 405,4000 82,3324 770,3252 3,2027
λ
i
137,2000 16,073 277,9222 1,0195
Tóxico λ 290,7000 58,4778 571,6920 2,9097
T
1
0,9902 0,9700 0,9999 0,0003
T
2
0,7432 0,4540 0,9789 0,0059
T
3
0,3011 0,0355 0,5949 0,0065
T
4
0,1854 0,0108 0,4493 0,0060
T
5
0,0889 0,0001 0,2597 0,0036
Completo λ 371,1000 64,0547 741,6136 3,8180
λ
i
153,2000 17,1127 318,0987 1,9881
T
1
0,5033 0,0092 0,9529 0,0034
T
2
0,2493 0,0000 0,6919 0,0028
T
3
0,1251 0,0000 0,4524 0,0022
T
4
0,0636 0,0000 0,2581 0,0013
T
5
0,0317 0,0000 0,1459 0,0008
Os resultados das comparações dos modelos por meio do fator de Bayes
encontram-se na Tabela 4.
Os resultados do fator de Bayes são interpretados segundo a sugestão de
Kass & Raftery (1995), contida na Tabela 1. Verifica-se na Tabela 4, que a evidên-
cia a favor do modelo Completo, quando comparado com o de Toxidez, é forte.
Este mesmo modelo comparado com o de Interferência tem evidência positiva a
favor do modelo Completo. A evidência a favor do modelo de Toxidez, quando
comparado com o de Interferência, é fraca. Em relação ao modelo Usual, forte
evidência a favor dos modelos que estão no numerador.
O modelo Completo, quando comparado com os demais modelos, tem in-
dicação favorável, para o ajuste deste resultado experimental. Um segundo modelo
indicado seria o de Toxidez. Nota-se que o modelo Usual, que comumente é uti-
lizado, é o menos provável para o tipo do resultado analisado.
47
TABELA 4 Fator de Bayes aproximado e grau de evidência em favor do modelo
da margem esquerda, em comparação ao modelo da margem superior,
considerando o resultado experimental (01300).
Denominador
Numerador Toxidez Interferência Usual
Completo 2, 6209
F
7, 4376
P
5, 0344 × 10
153
MF
Toxidez 2, 8378
F r
1, 9209 × 10
153
MF
Interferência 6, 7689 × 10
153
MF
MF: Muito forte; F: Forte; P: Positiva; Fr:Fraca
4.4.2 Convergência do processo de Markov e estimação do número mais
provável para o resultado experimental (33300)
O resultado experimental (33300) representa a situação na qual cresci-
mento de microrganismos em baixas diluições. Observa-se que, nas três primeiras
diluições, os três tubos apresentaram crescimento de microrganismo.
Observa-se, por meio da Tabela 5, que a convergência das cadeias dos parâ-
metros dos modelos foi satisfeita por meio do critério de Raftery & Lewis (1992).
Com exceções do parâmetro λ
i
no modelo de Interferência e dos parâmetros T
3
,
T
4
e T
5
no modelo de Toxidez, todos os demais tiveram suas médias estimadas
& .emBT/F52 10.909 Tf 130.581 582u- 130.581mJ9u- -f2-25ros
TABELA 5 Convergência dos parâmetros dos modelos para o resultado experi-
mental (33300).
Raftery & Lewis Heidelberg & Welch
Modelo Parâmetros FD estacionária Half-width
Usual λ 1,0491 sim sim
Interferência λ 0,9826 sim sim
λ
i
1,0152 não sim
Tóxico λ 1,0491 sim sim
T
1
1,0320 sim sim
T
2
0,9823 sim sim
T
3
1,0843 não não
T
4
1,0664 sim não
T
5
1,1022 sim não
Completo λ 0,9823 sim sim
λ
i
0,9663 sim sim
T
1
0,9823 sim sim
T
2
0,9823 sim sim
T
3
1,0491 sim sim
T
4
1,0664 sim sim
T
5
1,0758 sim sim
Verifica-se que o traço da cadeia do parâmetro λ no modelo Usual, na
Figura 8, se encontra num estado de equilíbrio, ou seja, esta cadeia atingiu a con-
vergência.
Na Figura 9, observa-se que tanto a cadeia do parâmetro λ como a do
parâmetro λ
i
, convergiram no modelo de Interferência, ao contrário do que indica
o critério de Heidelberger & Welch (1983), para o parâmetro λ
i
.
Na Figura 4.10, verifica-se que a cadeia do parâmetro λ no modelo de Toxi-
dez atingiu convergência. Também observa-se, na Figura 11, que convergência
nas cadeias dos parâmetros T
1
, T
2
, T
3
, T
4
e T
5
deste modelo.
49
0 1000 2000 3000 4000 5000
0 200 400 600 800 1000
Iterações
λλ
0 200 400 600 800 1000 1200
0.000 0.001 0.002 0.003
Densidade
λλ
(a) (b)
FIGURA 8 Traço da cadeia gerada (a) e densidade a posteriori (b) do parâmetro
λ, do modelo Usual, para o resultado experimental (33300).
(a) (b)
0 1000 2000 3000 4000 5000
0 10 20 30 40
Iterações
λλ
I
0 500 1000 1500
0.0000 0.0005 0.0010 0.0015 0.0020 0.0025
Densidade
λλ
I
(c) (d)
FIGURA 9 Traço das cadeias geradas (a) e (c) e densidades a posteriori (b) e (d)
dos parâmetros λ, e λ
i
, do modelo de Interferência, para o resultado
experimental (33300).
(a) (b)
FIGURA 10 Traço da cadeia gerada (a) e densidade a posteriori (b) do parâme-
tro λ, para o modelo de Toxidez, segundo o resultado experimental
(33300).
50
(a) (b)
(c) (d)
(e) (f)
0 1000 2000 3000 4000 5000
0.0 0.2 0.4 0.6
Iterações
T
4
0.0 0.2 0.4 0.6 0.8
0 1 2 3 4
Densidade
T
4
(g) (h)
(i) (j)
FIGURA 11 Traço das cadeias geradas (a), (c), (e), (g) e (i) e densidades a pos-
teriori (b), (d), (f), (h) e (j) dos parâmetros T
1
, T
2
, T
3
, T
4
e T
5
, do
modelo de Toxidez, para o resultado experimental (33300).
51
Observando-se, na Figura 12, nota-se que as cadeias dos parâmetros λ e
λ
i
do modelo Completo atingiram convergência. Verifica-se, nas Figuras 13 e 14,
que convergência nas cadeias dos parâmetros T
1
, T
2
, T
3
, T
4
e T
5
deste modelo.
0 1000 2000 3000 4000 5000
0 500 1000 1500
Iterações
λλ
0 500 1000 1500
0.0000 0.0005 0.0010 0.0015 0.0020
Densidade
λλ
(a) (b)
(c) (d)
FIGURA 12 Traço das cadeias geradas (a) e (c) e densidades a posteriori (b), (d)
do parâmetro λ e λ
i
, respectivamente, do modelo Completo, para o
resultado experimental (33300).
0 1000 2000 3000 4000 5000
0.0 0.2 0.4 0.6 0.8 1.0
Iterações
T
1
0.0 0.2 0.4 0.6 0.8 1.0
0.0 0.2 0.4 0.6 0.8 1.0
Densidade
T
1
(a) (b)
FIGURA 13 Traço da cadeia gerada (a) e densidade a posteriori (b) do parâmetro
T 1, do modelo Completo, para o resultado experimental (33300).
52
(a) (b)
(c) (d)
0 1000 2000 3000 4000 5000
0.0 0.1 0.2 0.3 0.4 0.5 0.6
Iterações
T
4
0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7
0 5 10 15
Densidade
T
4
(e) (f)
0 1000 2000 3000 4000 5000
0.0 0.1 0.2 0.3 0.4 0.5
Iterações
T
5
0.0 0.1 0.2 0.3 0.4 0.5
0 10 20 30 40
Densidade
T
5
(g) (h)
FIGURA 14 Traço das cadeias geradas (a), (c), (c) e (f) e densidades a posteriori
(b), (d), (f) e (h) dos parâmetros T
2
, T
3
, T
4
e T
5
do modelo Completo,
para o resultado experimental (33300).
53
A média a posteriori, o intervalo de credibilidade e o erro de Monte Carlo
para o resultado experimental analisado encontram-se na Tabela 6. As estimati-
vas do NMP para os diferentes modelos são ligeiramente diferentes, no entanto
os intervalos de credibilidade se sobrepõem. Esta diferença é menos acentuada na
comparação entre os modelos de Toxidez e Completo. Foi possível obter estima-
tivas para o parâmetro (λ
i
), tanto no modelo de Interferência quanto no modelo
Completo. De modo geral, o erro de Monte Carlo foi pequeno e os intervalos de
credibilidade amplos.
TABELA 6 Estimativas das médias a posteriori e intervalos de credibilidade
(I.C.95%), dos parâmetros dos modelos e erro de Monte Carlo (Erro
M.C.) para o resultado (33300).
I.C. 95%
Modelos Parâmetros Média LI LS Erro M.C.
Usual λ 201,1867 3,8736 466,1230 1,7793
Interferência λ 377,0000 67,0797 716,6476 2,7262
λ
i
7,0500 0,1897 16,9798 0,0844
Tóxico λ 442,9499 112,2653 849,1340 2,9598
T
1
0,6557 0,3065 0,9540 0,0024
T
2
0,3884 0,1060 0,6790 0,0020
T
3
0,2321 0,2885 0,4609 0,0017
T
4
0,1521 0,0025 0,3410 0,0016
T
5
0,0758 0,0000 0,2150 0,0012
Completo λ 412,3000 95,2652 797,3210 2,7989
λ
i
11,9600 0,0223 30,6229 0,1552
T
1
0,5045 0,0160 0,9637 0,0044
T
2
0,2470 0,0000 0,6858 0,0037
T
3
0,1224 0,0000 0,4279 0,0024
T
4
0,0599 0,0000 0,2471 0,0013
T
5
0,0301 0,0000 0,1357 0,0007
54
TABELA 7 Fator de Bayes aproximado e grau de evidência em favor do modelo
da margem esquerda, em comparação ao modelo da margem superior,
considerando o resultado experimental (33300).
Denominador
Numerador Toxidez Interferência Usual
Completo 7, 7915 × 10
7
N
1, 1230 × 10
4
N
2, 9131 × 10
1
F
Toxidez 1, 4414 × 10
2
F
3, 7388 × 10
7
MF
Interferência 2, 5940 × 10
5
MF
N: Negativa; MF: Muito forte; F: Forte; Fr:Fraca
Os resultados das comparações dos modelos por meio do fator de Bayes
para o resultado experimental (33300) encontram-se na Tabela 7. Segundo a sug-
estão contida na Tabela 1, verifica-se, na Tabela 7, que a comparação do modelo
Completo, entre os modelos de Toxidez e de Interferência, leva a uma evidência
negativa a favor do modelo Completo, ou seja, o fator de Bayes indica o modelo
do denominador como mais plausível. Quando este mesmo modelo é comparado
com o modelo Usual, verifica-se que a evidência a seu favor é fraca, ou seja, seu
valores são aproximados. Na comparação do modelo de Toxidez com o modelo de
Interferência, observa-se que forte evidência a favor do primeiro modelo. Esta
evidência se torna muito forte quando este modelo é comparado com o Usual. A
evidência a favor do modelo de Interferência é muito forte quando comparado com
o Usual.
O modelo Toxidez, quando comparado com os demais modelos, tem indi-
cação favorável para o ajuste deste resultado experimental. Este resultado não era
esperado, uma vez que o resultado experimental (33300) sugere o modelo Usual
como o mais plausível. Uma possível explicação é a quantidade de parâmetros
do modelo de Toxidez (6) em relação a do modelo Usual (1) e a priori com alta
probabilidade de valores intermediários do efeito tóxico (U[0,1]).
55
4.5 Exemplos de aplicação considerando três diluições
Nesta seção foram ilustrados exemplos de aplicações do método utilizando-
se três diluições. Outros resultados estão no Anexo e seguem o mesmo padrão de
interpretação. Os resultados foram organizadas numa sequência em que, primeiro,
é apresentada a monitoração da convergência das cadeias dos parâmetros. Depois,
são apresentados os resultados descritivos dos modelos e, por último, as compara-
ções dos modelos por meio do Fator de Bayes. As estimativas de NMP, para os
casos de três diluições e três tubos, são tabeladas e seus resultados podem ser
comparados.
4.5.1 Convergência do processo de Markov e estimação do número mais
provável para o resultado experimental (030)
O resultado experimental (030) representa a situação em que o crescimento
do microrganismo é inibido em baixas diluições. Espera-se que o modelo Usual
não seja o mais adequado para esta situação.
A monitoração da convergência das cadeias dos parâmetros é apresentada
na Tabela 8. Por meio desta tabela, observa-se que todos os parâmetros dos mo-
delos convergiram, pelo critério de Raftery & Lewis. Com exceção do parâmetro
λ no modelo de Completo, os demais parâmetros tiveram convergência satisfeita
por meio do critério de Heidelberger & Welch.
56
TABELA 8 Convergência dos parâmetros dos modelos para o resultado experi-
mental (030).
Raftery & Lewis Heidelberg & Welch
Modelo Parâmetros FD estacionária Half-width
Usual λ 0,9823 sim sim
Interferência λ 0,9666 sim sim
λ
i
1,0152 sim sim
Tóxico λ 1,1022 sim sim
T
1
1,1022 sim sim
T
2
2,0448 sim sim
T
3
2,1276 sim sim
Completo λ 0,9986 não sim
λ
i
1,0152 sim sim
T
1
0,9823 sim sim
T
2
0,9986 sim sim
T
3
1,0152 sim sim
Verifica-se que o traço da cadeia do parâmetro λ no modelo Usual, na
Figura 15, se encontra num estado de equilíbrio, ou seja, atingiu convergência.
(a) (b)
FIGURA 15 Traço da cadeia gerada (a) e densidade a posteriori (b) do parâmetro
λ, do o modelo Usual, para o resultado experimental (030).
57
Por meio da Figura 16, observa-se que tanto a cadeia do parâmetro λ como
a cadeia do parâmetro λ
i
convergiram no modelo de Interferência.
(a) (b)
0 1000 2000 3000 4000 5000
0 100 200 300
Iterações
λλ
I
0 200 400 600 800 1000
0.000 0.001 0.002 0.003 0.004
Densidade
λλ
I
(c) (d)
FIGURA 16 Traço das cadeias geradas (a) e (c) e densidades a posteriori (b) e (d)
dos parâmetros λ, e λ
i
, do modelo de Interferência, para o resultado
experimental (030).
Por meio da Figura 17, verifica-se que a cadeia do parâmetro λ no modelo
de Toxidez atingiu convergência. Observa-se, na Figura 4.18, que convergência
nas cadeias dos parâmetros T
1
, T
2
e T
3
deste modelo.
(a) (b)
FIGURA 17 Traço da cadeia gerada (a) e densidade a posteriori (b) do parâmetro
λ, do modelo de Toxidez, para o resultado experimental (030).
58
(a) (b)
0 1000 2000 3000 4000 5000
0.0 0.2 0.4 0.6 0.8 1.0
Iterações
T
2
0.0 0.2 0.4 0.6 0.8 1.0
0.0 0.5 1.0 1.5
Densidade
T
2
(c) (d)
0 1000 2000 3000 4000 5000
0.0 0.2 0.4 0.6 0.8
Iterações
T
3
0.0 0.2 0.4 0.6 0.8 1.0
0 1 2 3 4
Densidade
T
3
(e) (f)
FIGURA 18 Traço das cadeias geradas (a), (c) e (e) e densidades a posteriori (b),
(d) e (f) dos parâmetros T
1
, T
2
e T
3
, do modelo de Toxidez, para o
resultado experimental (030).
Observando-se a Figura 19, nota-se que as cadeias dos parâmetros λ e λ
i
do modelo Completo atingiram convergência. Verifica-se, na Figura 4.20 que
convergência nas cadeias dos parâmetros T
1
, T
2
e T
3
deste modelo.
59
0 1000 2000 3000 4000 5000
0 200 400 600 800 1000
Iterações
λλ
0 200 400 600 800 1000 1200
0.000 0.001 0.002 0.003
Densidade
λλ
(a) (b)
(c) (d)
FIGURA 19 Traço das cadeias geradas (a) e (c) e densidades a posteriori (b) e
(d) do parâmetro λ e λ
i
, respectivamente, do modelo Completo, para
o resultado experimental (030).
(a) (b)
(c) (d)
(e) (f)
FIGURA 20 Traço das cadeias geradas (a), (c) e (e) e densidades a posteriori (b),
(d) e (f) dos parâmetros T
1
, T
2
e T
3
do modelo Completo, para o
resultado experimental (030).
60
A média a posteriori, o intervalo de credibilidade e o erro de Monte Carlo
para o resultado experimental analisado (030) encontram-se na Tabelas 9.
Foi possível obter estimativas para o parâmetro (λ
i
), tanto no modelo de
Interferência quanto no modelo Completo. Os intervalos de credibilidade se so-
brepõem em todos os modelos, mas as estimativas pontuais são diferentes. A in-
clusão de novos parâmetro reduziu os intervalos de credibilidade de λ. Ressalta-se
que, neste caso, o parâmetro λ não foi subestimado pelo modelo Usual. De modo
geral, o erro de Monte Carlo foi pequeno e os intervalos de credibilidade amplos.
A estimativa do NMP tabelado para o modelo Usual encontrado na litera-
tura, Anexo C, é 9,4 e seu intervalo de confiança tem como limite inferior o valor
3,6 e, como superior, o valor 38. Observa-se que tanto a estimativa pontual do
NMP como o intervalo de credibilidade encontrados neste trabalho são maiores
que os da literatura.
TABELA 9 Estimativas das médias a posteriori e intervalos de credibilidade
(I.C.95%), dos parâmetros dos modelos e erro de Monte Carlo (Erro
M.C.) para o resultado (030).
I.C. 95%
Modelos Parâmetros Média LI LS Erro M.C.
Usual λ 208,1000 4,6926 499,0277 2,6701
Interferência λ 196,2938 21,5420 431,4838 1,8544
λ
i
66,2979 3,7448 159,2506 0,7220
Tóxico λ 178,4559 15,5273 384,0858 1,6059
T
1
0,9735 0,9091 0,9999 0,0006
T
2
0,3387 0,0141 0,7297 0,0033
T
3
0,1778 0,0000 0,5030 0,0027
Completo λ 216,0392 15,8690 486,6142 2,3138
λ
i
85,6435 0,0413 0,9897 0,9384
T
1
0,5037 0,0413 0,9897 0,0039
T
2
0,2514 0,0000 0,7071 0,0030
T
3
0,1211 0,0000 0,4327 0,0021
61
TABELA 10 Fator de Bayes aproximado e grau de evidência em favor do modelo
da margem esquerda, em comparação ao modelo da margem supe-
rior, considerando o resultado experimental (030).
Denominador
Numerador Toxidez Interferência Usual
Completo 6, 24467 × 10
1
N
9, 14778 × 10
1
N
4, 4520 × 10
159
F
Toxidez 1, 4649 × 10
0
F
7, 1293 × 10
159
MF
Interferência 4, 8668 × 10
159
MF
N: Negativa; MF: Muito forte; F: Forte
Os resultados das comparações dos modelos por meio do fator de Bayes
encontram-se na Tabela 4.9.
Os resultados do fator de Bayes são interpretados segundo a sugestão de
Kass & Raftery (1995), contida na Tabela 1. Verifica-se, pela Tabela 10, que a
evidência a favor do modelo Completo, quando comparado com o de Toxidez e o
de Interferência é negativa, ou seja, o modelo do denominador é mais provável.
Por outro lado, forte evidência a favor deste modelo, quando é comparado com
o Usual. A evidência a favor do modelo de Toxidez, quando comparado com o de
Interferência, é forte. Este evidência é muito forte quando o modelo de Toxidez é
comparado com o Usual. Esta interpretação também é válida na comparação entre
o modelo de Interferência e o Usual.
Pode-se notar que o modelo que é indicado, por meio do fator de Bayes,
como o mais provável para o resultado experimental analisado, é o de Toxidez.
62
4.5.2 Convergência do processo de Markov e estimação do número mais
provável para o resultado experimental (330)
Nesta seção, são apresentados os resultados dos ajustes dos modelos uti-
lizando o resultado (330). Neste exemplo, não inibição do crescimento de mi-
crorganismos em baixas diluições e é esperado que o modelo Usual seja adequado.
A monitoração da convergência das cadeias dos parâmetros para o resul-
tado (330) encontra-se na Tabela 11. Por meio do critério de Raftery & Lewis,
todos os parâmetros dos modelos atingiram convergência. Apenas o parâmetro T
3
do modelo de Toxidez não convergiu, segundo o critério de Heidelberger & Welch.
TABELA 11 Convergência dos parâmetros dos modelos para o resultado experi-
mental (330).
Raftery & Lewis Heidelberg & Welch
Modelo Parâmetros FD estacionária Half-width
Usual λ 0,9663 sim sim
Interferência λ 0,9911 sim sim
λ
i
0,9826 sim sim
Tóxico λ 1,0152 sim sim
T
1
1,0320 sim sim
T
2
0,9986 sim sim
T
3
1,0491 não sim
Completo λ 1,0491 sim sim
λ
i
0,9823 sim sim
T
1
1,0320 sim sim
T
2
0,9823 sim sim
T
3
0,9823 sim sim
A monitoração informal da convergência das cadeias de cada parâmetro
dos modelos pode ser observada por meio das Figuras 21 a 26.
Verifica-se que o traço da cadeia do parâmetro λ no modelo Usual, na
Figura 21, se encontra num estado de equilíbrio, ou seja, atingiu convergência.
63
(a) (b)
FIGURA 21 Traço da cadeia gerada (a) e densidade a posteriori (b) do parâmetro
λ, do modelo Usual, para o resultado experimental (330).
Por meio da Figura 22, observa-se que tanto a cadeia do parâmetro λ como
a do parâmetro λ
i
convergiram no modelo de Interferência.
(a) (b)
(c) (d)
FIGURA 22 Traço das cadeias geradas (a) e (c) e densidades a posteriori (b) e (d)
dos parâmetros λ, e λ
i
, do modelo de Interferência, para o resultado
experimental (330).
Por meio da Figura 23, verifica-se que a cadeia do parâmetro λ no modelo
de Toxidez atingiu convergência. Observa-se, na Figura 24, que convergência
nas cadeias dos parâmetros T
1
, T
2
e T
3
deste modelo.
Observando-se a Figura 25, nota-se que as cadeias dos parâmetros λ e λ
i
64
0 1000 2000 3000 4000 5000
0 200 400 600 800 1000 1200
Iterações
λλ
0 200 400 600 800 1000
0.0000 0.0010 0.0020
Densidade
λλ
(a) (b)
FIGURA 23 Traço da cadeia gerada (a) e densidade a posteriori (b) do parâmetro
λ, do modelo de Toxidez, para o resultado experimental (330).
0 1000 2000 3000 4000 5000
0.0 0.2 0.4 0.6 0.8 1.0
Iterações
T
1
(a) (b)
(c) (d)
(e) (f)
FIGURA 24 Traço das cadeias geradas (a), (c) e (d) e densidades a posteriori (b),
(d) e (f) dos parâmetros T
1
, T
2
e T
3
, do modelo de Toxidez, para o
resultado experimental (330).
do modelo Completo atingiram convergência. Verifica-se, na Figura 26, que
convergência nas cadeias dos parâmetros T
1
, T
2
e T
3
deste modelo.
65
0 1000 2000 3000 4000 5000
0 200 400 600 800 1000
Iterações
λλ
0 200 400 600 800 1000
0.0000 0.0010 0.0020 0.0030
Densidade
λλ
(a) (b)
(c) (d)
FIGURA 25 Traço das cadeias geradas (a) e (c) e densidades a posteriori (b), (d)
do parâmetro λ e λ
i
, respectivamente, do modelo Completo, para o
resultado experimental (330).
(a) (b)
(a) (b)
0 1000 2000 3000 4000 5000
0.0 0.2 0.4 0.6 0.8 1.0
Iterações
T
3
0.0 0.2 0.4 0.6 0.8 1.0
0 1 2 3 4 5 6 7
Densidade
T
3
(a) (b)
FIGURA 26 Traço das cadeias geradas (a) e (c) e densidades a posteriori (b) e (d)
dos parâmetros T
1
, T
2
e T
3
, do modelo Completo, para o resultado
experimental (330).
66
Os resultados referentes à média a posteriori, seus respectivos intervalos
de credibilidade e erros de Monte Carlo para os parâmetros dos modelos encontram-
se na Tabela 12. Observa-se, nesta tabela, que, embora as estimativas do NMP para
os diferentes modelos sejam diferentes, estas não são muito discrepantes. Foi pos-
sível obter estimativas para o parâmetro (λ
i
), tanto no modelo de Interferência
quanto no modelo Completo. De modo geral, o erro de Monte Carlo foi pequeno
e os intervalos de credibilidade amplos. No Anexo C, encontra-se o valor tabelado
deste resultado experimental que é, 240 e seu intervalo de confiança tem como
limite inferior 42 e, como limite superior, 1000. Pode-se notar que o intervalo de
credibilidade apresentado na Tabela 12 é menor.
TABELA 12 Estimativas das médias a posteriori e intervalos de credibilidade
(I.C.95%), dos parâmetros dos modelos e erro de Monte Carlo (Erro
M.C.) para o resultado (330).
I.C. 95%
Modelos Parâmetros Média LI LS Erro M.C.
Usual λ 216,4344 33,2168 429,6894 1,7664
Interferência λ 234,3187 49,8072 470,5193 1,3434
λ
i
7,3756 0,0391 18,2483 0,0774
Tóxico λ 300,7133 49,5183 615,8017 2,0278
T
1
0,6076 0,2330 0,9540 0,0032
T
2
0,3277 0,0259 0,6453 0,0024
T
3
0,1803 0,0000 0,4615 0,0020
Completo λ 253,9760 49,0583 522,2672 3,6565
λ
i
12,5890 0,0793 33,9589 0,1299
T
1
0,4920 0,0003 0,9455 0,0046
T
2
0,2490 0,0000 0,6954 0,0032
T
3
0,1244 0,0000 0,4427 0,0020
Os resultados das comparações dos modelos por meio do fator de Bayes
para o resultado experimental (330) encontram-se na Tabela 13.
67
TABELA 13 Fator de Bayes aproximado e grau de evidência em favor do modelo
da margem esquerda, em comparação ao modelo da margem supe-
rior, considerando o resultado experimental (330).
Denominador
Numerador Toxidez Interferência Usual
Completo 7, 6115 × 10
9
N
6, 9353 × 10
6
N
2, 7215 × 10
9
N
Toxidez 9, 1116 × 10
2
MF
3, 5755 × 10
1
N
Interferência 3, 9242 × 10
4
N
N: Negativa; MF: Muito forte
Segundo a sugestão contida na Tabela 1, verifica-se na Tabela 13, que a
comparação do modelo Completo com os modelos de Toxidez, de Interferência
e Usual tem evidência negativa, o que significa dizer que o modelo que está no
denominador é mais provável que o do numerador.
O modelo de Toxidez, quando comparado com o de Interferência, tem
evidência muito forte, a favor do primeiro modelo. Por outro lado a evidência a
favor do modelo de Toxidez quando comparado com o Usual é negativa, ou seja,
o modelo Usual é mais plausível que o de Toxidez. Para o resultado experimental
(330), em questão, conclui-se que o modelo Usual é o mais provável para ajustar
os dados analisados, como era de se esperar.
68
4.6 Considerações gerais
O parâmetro T apresenta estimativas muito altas em resultados experimen-
tais nos quais isto não seria esperado, por exemplo (33300). Isso pode ser devido à
distribuição a priori uniforme para T. Seria, provavelmente, mais realista assumir
distribuições a priori do tipo Beta (a, b), com a + b = 1. Sendo assim, não seria
muito provável encontrar valores intermediários de toxidez.
Em situações em que, provavelmente, não efeito tóxico, a inclusão deste
parâmetro aumenta a amplitude do intervalo de credibilidade das estimativas de λ
e seria melhor não utilizar tal modelo.
Em geral, o modelo Completo foi melhor, mas deve-se ter cuidado com a
sensibilidade deste modelo à especificação da distribuição a priori utilizada. Sabe-
se, no entanto, que, para a maioria dos problemas de microbiologia de alimentos,
podem-se especificar distribuições a priori mais informativas do que as usadas
para λ e λ
i
(λ
i
deve ser por volta de 10 vezes menor que λ) e distribuições a priori
mais pesadas nos extremos para T . Isto é, claramente, uma boa perspectiva para o
uso do modelo Completo.
O programa implementado não se limita apenas a três e cinco diluições
com três tubos cada e sua utilização pode ser estendida para outros possíveis re-
sultados experimentais. Isso abre possibilidades para a implementação de outros
delineamentos mais adequados para estimar T , por exemplo.
Para estudos futuros, modelagens conjuntas de experimentos de diluição
seriada com auxílio de programas e modelos mais flexíveis, como os apresentados,
podem ser utilizadas.
69
5 CONCLUSÕES
Foi estabelecida análise bayesiana flexível para diversos números de dilui-
ções e tubos.
Foi possível comparar modelos com e sem parâmetros de interferência e
Toxidez.
Em situações experimentais em que o crescimento do microrganismo de
interesse é inibido em baixas diluições, o modelo Usual não é adequado e tem
valores subestimados para o NMP. Nas situações em que o modelo Usual é o mais
adequado, o modelo Completo tem estimativas que não são muito diferentes das
deste.
O modelo Completo pode ser utilizado em programas mais gerais para a
análise de experimentos de diluição seriada.
70
6 REFERÊNCIAS BIBLIOGRÁFICAS
AGÊNCIA NACIONAL DE VIGILÂNCIA SANITÁRIA. Resolução-RDC n. 12,
de 2 de jan. de 2001. Disponível em:<www.abic.com.br/arquivos/leg_resolucao12
_01_anvisa.pdf>. Acesso em: 2007.
AOAC International. Bacteriological Analytical Manual. Gaithersburg, Mary-
land, 2001. (Most Probable Number from Dilution, Appendix, 2).
BLODGETT, R. J.; GARTHRIGHT, W. E. Several MPM models for serial dilution
with suppressed growth at low dilutions. The Food Microbiology, v.15, p.91-99,
1998.
BOX, G.E.; TIAO, G.C. Bayesian inference in statistical analysis. New York:
Wiley, 1992. 588p.
CASELA, G.; GEORGE, G. I. Explaining the gibbs sampler. The American Sta-
tistician, v.46, n.3, 167-174, 1992.
CHIB, S.; GREENBERG, E. Understanding the Metropolis-Hasting Algorithm.
The American Statistician, Alexandria, v.49, n.4, p.327-335, Nov. 1995.
CLOUGH, H. E.; CLANCY, D; O’NEILL, P.D.; ROBINSON, S.E.; FRENCH,
N.P. Quantifying uncertainty associated with microbial count data: a Bayesian ap-
proach. Biometrics, v.61, p.610-616, June 2005
COCHRAN, W. G. Estimating of bacteriological densities by means of the "Most
Probable Number". The Biometrics, v.6, p. 105-116, 1950.
DEMÉTRIO, C. G. B. Modelos lineares generalizados em experimentação agro-
nômica, Esalq/Usp. Piracaba, 2002. 113p.
EHLERS, R. S. Introducão a inferência bayesiana Departamento de Estatística,
UFPR. Disponível em: <http:://www.est.ufrp.br/ ehlers/notas>. Acesso em 20 de
jul. 2007.
FISHER, R. A. On the matematical foundation of theorical statistics.Philosophical
Transaction of the Royal Society, 222, p.309-368, 1922.
GEMAM, S.; GEMAM, D.Stochastic relaxation, Gibbs distributions and the baye-
sian restoration of images. IEEE Transaction on pattern analysis and machine
intellingence, Los Alamitos, v. 6, n. 6, p. 721-741,1984.
71
GELFAND, A.E. Gibbs Sampling. Journal of the American Statistical Associ-
ation, v.95, p.1300-1304, 2000.
GELFAND A.E. and A.F.M. Smith Sampling based approaches to calculating
marginal densities Journal of the American Statistical AND ginv074 0 T,RE-0tl-18(ginv.)3l4cDngthe theAme-1304a4Fa5662d5663r9 689.348 ci48(of)-247of of9E23r.28(ginv.)18(.020n3AD(the)B49(Gibbs)-4u6.3890B20n3AD(the)BA9v2FRE-0tl-18(gin5Rd0CUn5Rd147(the)3)3)3.0n5Rd0CUnA147(theoJ)18(U2(ginv.9 68.35Fo9 689e4lh.L1)Fa5662o9 689eEc40tl-689eEc40tl-p8spTM9 T4doW1-HM2mGFRE-5663r)1A785I 10.90989eEc44s2o92dtbs
Lisboa: Fundação Calouste Gulbenkian, 2003. 446p.
RAFTERY, A. E.; LEWIS, S. How many iterations in the Gibbs sampler? In:
BERNARDO, J. M.; BERGER, J. O.; DAWID, A. P.; SMITH, A. F. M. (Ed.).
Bayesian Statistics 4. Oxford: University Press, 1992. p.763-773b.
RAFTERY, A. E.; NEWTON, A. M. Estimating the integrated likelihood via pos-
terior simulation using the harmonic mean identity. Bayesian Estatistics. Oxford,
USA: University Press, v.8, p.1-45, 2007.
RAFTERY, A. Hypotesis testings and model selection via posterior simulation.
In: GILKS, W. R.; SPIEGELHALTER, D. J.; RICHARDSON, S. (Ed). Pratical
Markov chain Monte Carlo. Oxford, USA: University Press, 1995.
R DEVELOPMENT CORE TEAM The R Manuals: current version: 2.6.1.
2007. Disponível em:<http:www.R-project.org>. Acesso em: 2007.
SILVA, N.; JUNQUEIRA, N. F. A. The Manual de métodos de análise micro-
biológica de alimentos, Editora Varela: São Paulo, 1997. 295p
73
ANEXOS A
TABELAS Páginas
TABELA 1A Convergência dos parâmetros dos modelos para o resultado experi
mental (00330). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 76
TABELA 2A Estimativas das médias a posteriori e intervalos de credibilidade
(I.C.95%), dos parâmetros dos modelos e erro de Monte Carlo
(Erro M.C.) para o resultado experimental (00330). . . . . . . . . . . . . 76
TABELA 3A Fator de Bayes aproximado e grau de evidência em favor do modelo
da margem esquerda, em comparação ao modelo da margem
superior, considerando o resultado experimental (00330). . . . . . . . 77
TABELA 4A Convergência dos parâmetros dos modelos para o resultado experi
mental (03320). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 77
TABELA 5A Estimativas das médias a posteriori e intervalos de credibilidade
(I.C.95%), dos parâmetros dos modelos e erro de Monte Carlo
(Erro M.C.) para o resultado (03320). . . . . . . . . . . . . . . . . . . . . . . . . 78
TABELA 6A Fator de Bayes aproximado e grau de evidência em favor do modelo
da margem esquerda, em comparação ao modelo da margem
superior, considerando o resultado experimental (03320). . . . . . . . 78
74
TABELA 7A Convergência dos parâmetros dos modelos para o resultado experi
mental (32300). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 79
TABELA 8A Estimativas das médias a posteriori e intervalos de credibilidade
(I.C.95%), dos parâmetros dos modelos e erro de Monte Carlo
(Erro M.C.) para o resultado experimental (32300). . . . . . . . . . . . . 79
TABELA 9A Fator de Bayes aproximado e grau de evidência em favor do modelo
da margem esquerda, em comparação ao modelo da margem
superior, considerando o resultado experimental (32300). . . . . . . . 80
TABELA 10A Convergência dos parâmetros dos modelos para o resultado experi
mental (032). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 80
TABELA 11A Estimativas das médias a posteriori e intervalos de credibilidade
(I.C.95%), dos parâmetros dos modelos e erro de Monte Carlo
(Erro M.C.) para o resultado experimental (032). . . . . . . . . . . . . . . 81
TABELA 12A Fator de Bayes aproximado e grau de evidência em favor do modelo
da margem esquerda, em comparação ao modelo da margem
superior, considerando o resultado experimental (032). . . . . . . . . . 81
75
TABELA 1A Convergência dos parâmetros dos modelos para o resultado experimental
(00330).
Raftery & Lewis Heidelberg & Welch
Modelo Parâmetros FD estacionária Half-width
Usual λ 1,0320 sim sim
Interferência λ 0,9826 sim sim
λ
i
1,0152 não sim
Tóxico λ 0,9823 sim sim
T
1
4,2819 sim sim
T
2
3,1425 sim sim
T
3
35,0560 não sim
T
4
121,5928 sim não
T
5
28,1230 sim não
Completo λ 1,0491 sim sim
λ
i
1,0320 sim sim
T
1
0,9663 sim sim
T
2
0,9823 sim sim
T
3
0,9986 sim sim
T
4
0,9823 sim sim
T
5
0,9986 sim sim
TABELA 2A Estimativas das médias a posteriori e intervalos de credibilidade (I.C.95%),
dos parâmetros dos modelos e erro de Monte Carlo (Erro M.C.) para o re-
sultado experimental (00330).
I.C. 95%
Modelos Parâmetros Média LI LS Erro M.C.
Usual λ 218,3000 5,6004 511,1276 2,2727
Interferência λ 644,7000 219,0647 1.138,1640 3,2997
λ
i
228,8000 48,1664 447,4568 1,7686
Tóxico λ 519,2000 174,9129 973,4095 5,0606
T
1
0,9932 0,9784 0,9999 0,0005
T
2
0,9840 0,7594 0,9989 0,0066
T
3
0,2561 0,0330 0,5971 0,0111
T
4
0,1177 0,0020 0,3228 0,0064
T
5
0,0647 0,0000 0,2016 0,0050
Completo λ 652,5000 206,9912 1.147,4950 4,4804
λ
i
238,3000 35,3088 494,8204 3,3278
T
1
0,5056 0,0084 0,9575 0,0033
T
2
0,2506 0,0000 0,7086 0,0028
T
3
0,1238 0,0000 0,4315 0,0018
T
4
0,0615 0,0000 0,2493 0,0012
T
5
0,0309 0,0000 0,1327 0,0007
76
TABELA 3A Fator de Bayes aproximado e grau de evidência em favor do modelo da
margem esquerda, em comparação ao modelo da margem superior, consi-
derando o resultado experimental (00330).
Denominador
Numerador Toxidez Interferência Usual
Completo 3, 7132 × 10
1
N
1, 0770 × 10
1
N
3, 0120 × 10
188
N
Toxidez 2, 9006 × 10
3
N
8, 1117 × 10
186
MF
Interferência 2, 7965 × 10
189
MF
N: Negativa; MF: Muito forte
TABELA 4A Convergência dos parâmetros dos modelos para o resultado experimental
(03320).
Raftery & Lewis Heidelberg & Welch
Modelo Parâmetros FD estacionária Half-width
Usual λ 0,9986 sim sim
Interferência λ 0,9826 sim sim
λ
i
0,9666 não sim
Tóxico λ 1,0491 sim sim
T
1
7,7976 sim sim
T
2
113,0165 sim sim
T
3
95,62258 não sim
T
4
105,62305 sim não
T
5
63,84735 sim não
Completo λ 1,0152 sim sim
λ
i
1,0152 não sim
T
1
1,0491 sim sim
T
2
1,0664 sim sim
T
3
1,0320 sim sim
T
4
1,0320 sim sim
T
5
1,0152 sim sim
77
TABELA 5A Estimativas das médias a posteriori e intervalos de credibilidade (I.C.95%),
dos parâmetros dos modelos e erro de Monte Carlo (Erro M.C.) para o re-
sultado experimental (03320).
I.C. 95%
Modelos Parâmetros Média LI LS Erro M.C.
Usual λ 208,2000 5,0886 486,0129 2,5468
Interferência λ 544,9000 172,2502 1.020,2540 3,5270
λ
i
67,3500 3,1913 151,0249 0,7421
Tóxico λ 605,7000 200,9531 1.068,4620 4,1670
T
1
0,9967 0,9908 0,9999 0,0002
T
2
0,4286 0,1715 0,7464 0,0194
T
3
0,2161 0,0302 0,3941 0,0084
T
4
0,1149 0,0072 0,2496 0,0072
T
5
0,0607 0,0005 0,1873 0,0064
Completo λ 582,8000 172,2280 1.083,8610 3,8878
λ
i
81,1000 7,1491 175,6933 0,8247
T
1
0,4954 0,0041 0,9457 0,0039
T
2
0,2488 0,0000 0,6830 0,0027
T
3
0,1241 0,0000 0,4350 0,0017
T
4
0,0607 0,0000 0,2455 0,0010
T
5
0,0304 0,0000 0,1360 0,0007
TABELA 6A Fator de Bayes aproximado e grau de evidência em favor do modelo da
margem esquerda, em comparação ao modelo da margem superior, consi-
derando o resultado experimental (03320).
Denominador
Numerador Toxidez Interferência Usual
Completo 1, 5927 × 10
2
N
5, 5373 × 10
1
N
2, 6500 × 10
159
MF
Toxidez 3, 4766 × 10
1
F
1, 6638 × 10
161
MF
Interferência 4, 7857 × 10
159
MF
N: Negativa; Forte; MF: Muito forte
78
TABELA 7A Convergência dos parâmetros dos modelos para o resultado experimental
(32300).
Raftery & Lewis Heidelberg & Welch
Modelo Parâmetros FD estacionária Half-width
Usual λ 1,0843 não sim
Interferência λ 1,0152 sim sim
λ
i
0,9826 não sim
Tóxico λ 0,9823 sim sim
T
1
1,0491 sim sim
T
2
0,9986 sim sim
T
3
1,1390 não sim
T
4
1,1022 sim sim
T
5
1,1390 sim sim
Completo λ 1,0320 sim sim
λ
i
1,0320 não sim
T
1
1,1022 sim sim
T
2
1,0157 sim sim
T
3
0,9986 sim sim
T
4
1,0157 sim sim
T
5
0,9823 sim sim
TABELA 8A Estimativas das médias a posteriori e intervalos de credibilidade (I.C.95%),
dos parâmetros dos modelos e erro de Monte Carlo (Erro M.C.) para o re-
sultado experimental (32300).
I.C. 95%
Modelos Parâmetros Média LI LS Erro M.C.
Usual λ 175,6000 0,0194 571,6584 3,1719
Interferência λ 345,000 50,1233 719,1767 2,7532
λ
i
8,9480 0,1818 22,2554 0,0974
Tóxico λ 391,1000 108,5183 726.1933 2,3997
T
1
0,7435 0,4591 0.9696 0,0021
T
2
0,5579 0,2420 0.8673 0,0025
T
3
0,2907 0,0422 0.5557 0,0021
T
4
0,1905 0,0100 0.4191 0,0016
T
5
0,0972 0,0000 0.2802 0,0013
Completo λ 386,8000 75,9087 741,4752 5,0052
λ
i
16,0100 0,0021 42,7326 0,2208
T
1
0,5037 0,0406 0,9892 0,0031
T
2
0,2535 0,0000 0,6962 0,0028
T
3
0,1266 0,0000 0,4419 0,0021
T
4
0,6413 0,0000 0,2600 0,0012
T
5
0,3305 0,0000 0,1416 0,0081
79
TABELA 9A Fator de Bayes aproximado e grau de evidência em favor do modelo da
margem esquerda, em comparação ao modelo da margem superior, consi-
derando o resultado experimental (32300).
Denominador
Numerador Toxidez Interferência Usual
Completo 5, 8283 × 10
12
N
1, 5502 × 10
8
N
1, 0645 × 10
11
MF
Toxidez 2, 6598 × 10
3
F
1, 8264 × 10
22
MF
Interferência 6, 8668 × 10
18
MF
N: Negativa; Forte; MF: Muito forte
TABELA 10A Convergência dos parâmetros dos modelos para o resultado experimental
(032).
Raftery & Lewis Heidelberg & Welch
Modelo Parâmetros FD estacionária Half-width
Usual λ 1,0152 sim sim
Interferência λ 0,9989 sim sim
λ
i
0,9989 não sim
Tóxico λ 0,9986 sim sim
T
1
1,0491 sim sim
T
2
6,2135 sim sim
T
3
3,9540 não sim
Completo λ 0,9986 sim sim
λ
i
1,0152 sim sim
T
1
1,3151 sim sim
T
2
0,9986 sim sim
T
3
1,0320 sim sim
80
TABELA 11A. Estimativas das médias a posteriori e intervalos de credibilidade
(I.C.95%), dos parâmetros dos modelos e erro de Monte Carlo (Erro
M.C.) para o resultado experimental (032).
I.C. 95%
Modelos Parâmetros Média LI LS Erro M.C.
Usual λ 165,7000 0,5299 413,0243 1,7210
Interferência λ 325,3666 74,7787 668,4529 2,0836
λ
i
68,2019 4,2126 156,3212 0,7208
Tóxico λ 315,8596 57,2441 607,4859 2,0355
T
1
0,9891 0,9645 0,9999 0,0002
T
2
0,3280 0,0176 0,6926 0,0041
T
3
0,1346 0,0000 0,3865 0,0024
Completo λ 353,8036 70,6718 675,5638 1,9989
λ
i
86,0443 8,9600 188,7087 0,9794
T
1
0,5009 0,0154 0,9643 0,0036
T
2
0,2514 0,0000 0,6982 0,0031
T
3
0,1261 0,0000 0,4421 0,0018
TABELA 12A Fator de Bayes aproximado e grau de evidência em favor do modelo da
margem esquerda, em comparação ao modelo da margem superior, consi-
derando o resultado experimental (032).
Denominador
Numerador Toxidez Interferência Usual
Completo 1, 2199 × 10
1
F
4, 4141 × 10
2
N
6, 2904 × 10
119
MF
Toxidez 3, 6184 × 10
1
N
5, 1564 × 10
120
MF
Interferência 1, 4251 × 10
121
MF
N: Negativa; Forte; MF: Muito forte
81
ANEXOS B
PROGRAMAS Páginas
PROGRAMA 1B. Rotina para gerar cadeias via MCMC do modelo Usual. . . 83
PROGRAMA 2B. Rotina para gerar cadeias via MCMC do modelo de Interferên-
cia. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 85
PROGRAMA 3B. Rotina para gerar cadeias via MCMC do modelo de Toxidez.
. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 88
PROGRAMA 4B. Rotina para gerar cadeias via MCMC do modelo Completo. 91
82
PROGRAMA 1 B Rotina para gerar cadeias via MCMC do modelo Usual.
# Definição do delineamento #
y <- c(0,0,3,3,2)
r <- 3
d <- length(y)
n <- c(rep(r,length(y)))
p <- length(y)+1
z <- 10^(1-c(2:p))
q <- 2
T <- 1000
Nef <- 10000
# queima (B), salto (J) e tamanho da amostra (N)
B <- 50000
J <- 200
N <- B + Nef
*
J
# hiperparâmetros para iniciar as candidatas #
alphaL <- 2
betaL <- 10
acL <- 2
bcL <- 10
a <- rnorm(T,alphaL
*
betaL,sqrt(alphaL
*
betaL^ 2))
l <- 1
lc <- 1
l <- rgamma(1,alphaL,1/betaL)
pa <- 1
npa <- 0
for(iter in 1:N)
{
lc <- rgamma(1,acL,1/bcL)
#Definindo partes da verossimilhança#
f <- gamma(n+1)/(gamma(y+1)
*
gamma(n-y+1))
pno <- exp(-l
*
z)
pnoc <- exp(-lc
*
z)
#Verossimilhança#
v <-(1-pno)^y
*
pno^(n-y)
vc <- (1-pnoc)^y
*
pnoc^(n-y)
# Razão das probabilidades nas candidatas
pc <- dgamma(lc,alphaL,1/betaL)
p <- dgamma(l ,alphaL,1/betaL)
hc <- dgamma(lc,acL,1/bcL)
h <- dgamma(l ,acL,1/bcL)
critL <- exp(sum(log(vc)-log(v)))+h/hc + pc/p
83
ifelse(critL>0, pa > 1, pa <- exp(critL))
# Troca ponto atual por candidato aceito #
if(pa > runif(1))
{
l <- lc
npa <- npa +1
}
f <- gamma(n+1)/(gamma(y+1)
*
gamma(n-y+1))
v <- f
*
(1-pno)^y
*
pno^(n-y)
vero <- exp(sum(log(v)))
#Inclui T valores da posteriori#
if(iter<=T)
{
a[iter] <- l
}
if(iter>B && (iter-B)%%J==0)
{
a[2:T] <- a[1:(T-1)]
a[1] <- l
mL <- mean(a)
vL <- var(a)+.1
bcL <- vL/mL
acL <- mL/bcL
write(t(c(l,vero)),file= "arquivo", q, append=TRUE)
}
}
84
PROGRAMA 2 B Rotina para gerar cadeias via MCMC do modelo de Interferência.
# Definição do delineamento #
y <- c(0,0,3,3,2)
r <- 3
d <- length(y)
n <-c(rep(r,length(y)))
p <- length(y)+1
z <- 10^(1-c(2:p))
q <- 2
T <- 1000
Nef <- 10000
# queima (B), salto (J) e tamanho da amostra (N)
B <- 50000
J <-200
N <- B + Nef
*
J
# hiperparâmetros para INICIAR as candidatas #
alphaLI <- 5
betaLI <- 5
alphaL <- 5
betaL <- 50
bcLI <- 2
acLI <- 10
bcL <- 2
acL <- 100
a <- rnorm(T,alphaLI
*
betaLI,sqrt(alphaLI
*
betaLI^2))
b <- rnorm(T,alphaL
*
betaL,sqrt(alphaL
*
betaL^2))
################ Valores iniciais para candidatas
l <- 1:2
lc <- 1:2
lc[1] <- rgamma(1,acLI,1/bcLI) # chute para lambda_I
lc[2] <- rgamma(1,acL,1/bcL) # chute para lambda
#Inicia a saida aceitando o chute
pa <- 1:2
pa[1] <- 1
pa[2] <- 1
l <- lc
npa <- 1:2
npa[1] <- 0
npa[2] <- 0
for(iter in 1:N)
{
# Atualizando lambda_I #
lc[1] <- rgamma(1,acLI,1/bcLI)
85
#Definindo partes da verossimilhança #
f <- gamma(n+1)/(gamma(y+1)
*
gamma(n-y+1))
po <- exp(-l[1]
*
z)
*
(1-exp(-l[2]
*
z))
poc <- exp(-lc[1]
*
z)
*
(1-exp(-l[2]
*
z))
# Verossimilhança
v <- po^y
*
(1 - po)^(n-y)
vc <- poc^y
*
(1 - poc)^(n-y)
# Razão das probabilidades nas candidatas
pc <- dgamma(lc[1],alphaLI,1/betaLI)
p <- dgamma(l[1] ,alphaLI,1/betaLI)
hc <- dgamma(lc[1],acLI,1/bcLI)
h <- dgamma(l[1] ,acLI,1/bcLI)
critLI <- sum(log(vc)-log(v))+h/hc + pc/p
ifelse(critLI>0, pa[1] > 1, pa[1] <- exp(critLI))
# Troca ponto atual por candidato aceito #
if(pa[1] > runif(1))
{
l[1] <- lc[1]
npa[1] <- npa[1] +1
}
# Atualizando lambda #
lc[2] <- rgamma(1,acL,1/bcL)
#Definindo partes da verossimilhança #
f <- gamma(n+1)/(gamma(y+1)
*
gamma(n-y+1))
po <- exp(-l[1]
*
z)
*
(1-exp(-l[2]
*
z))
poc <- exp(-l[1]
*
z)
*
(1-exp(-lc[2]
*
z))
# Verossimilhança
v <- po^y
*
(1 - po)^(n-y)
vc <- poc^y
*
(1 - poc)^(n-y)
# Razão das probabilidades nas candidatas
pc <- dgamma(lc[2],alphaL,1/betaL)
p <- dgamma(l[2] ,alphaL,1/betaL)
hc <- dgamma(lc[2],acL,1/bcL)
h <- dgamma(l[2] ,acL,1/bcL)
critL <- sum(log(vc)-log(v))+h/hc + pc/p
ifelse(critL>0, pa[2] > 1, pa[2] <- exp(critL))
86
# Troca ponto atual por candidato aceito #
if(pa[2] > runif(1))
{
l[2] <- lc[2]
npa[2] <- npa[2] +1
v <- vc
}
f <- gamma(n+1)/(gamma(y+1)
*
gamma(n-y+1))
po <- exp(-l[1]
*
z)
*
(1-exp(-l[2]
*
z))
v <- f
*
po^y
*
(1 - po)^(n-y)
vero <- exp(sum(log(v)))
# Armazena na cadeia atual #
if(iter<=T)
{
a[iter] <- l[1]
b[iter] <- l[2]
}
if((iter>B)&&((iter-B)%%J==0))
{
a[2:T] <- a[1:(T-1)]
b[2:T] <- b[1:(T-1)]
a[1] <- l[1]
b[1] <- l[2]
mLI <- mean(a)
vLI <- var(a)+.1
mL <- mean(b)
vL <- var(b)+.1
bcLI <- vLI/mLI
acLI <- mLI/bcLI
bcL <- vL/mL
acL <- mL/bcL
write(t(c(l[1],l[2],vero)),file= "arquivo",q, append=TRUE)
}
}
87
PROGRAMA 3 B Rotina para gerar cadeias via MCMC do modelo de Toxidez.
# Definindo delineamento #
y <- c(0,1,3,0,0)
r <- 3
d <- length(y)
n <- c(rep(r,length(y)))
p <-length(y)+1
z <- 10^(1-c(2:p))
T <- 1-sort(runif(d))
Tc <- T
s <- c(1,T[1:(d-1)])
sc <- s
q <- 2
Nef <- 10000
# queima (B), salto (J), tamanho da amostra (N)
B <- 5000
J <- 100
N <- B + Nef
*
J
M <- 1000
# parâmetros iniciais para as candidatas #
alphaL <- 2
betaL <- 100
acL <- 2
bcL <- 100
a<- rnorm(M,alphaL
*
betaL,sqrt(alphaL
*
betaL^ 2))
# Chutes iniciais para candidatas #
l <- 1
lc <- 1
l <- rgamma(1,alphaL,1/betaL) # chute para lambda
#Inicia a saida aceitando o chute
pa <- 1:2
npa <- 1:2
npa[1] <- 0
npa[2] <- 0
for(iter in 1:N)
{
lc <- rgamma(1,acL,1/bcL)
# Definindo partes da verossimilhança #
f <- gamma(n+1)/(gamma(y+1)
*
gamma(n-y+1))
88
pno <- exp(-l
*
z
*
(1-T))
pnoc <- exp(-lc
*
z
*
(1-T))
# Verossimilhança #
v <- (1-pno)^y
*
pno^(n-y)
vc <- (1-pnoc)^y
*
pnoc^(n-y)
# Razão das probabilidades nas candidatas
pc <- dgamma(lc,alphaL,1/betaL)
p <- dgamma(l ,alphaL,1/betaL)
hc <- dgamma(lc,acL,1/bcL)
h <- dgamma(l ,acL,1/bcL)
critL <- sum(log(vc/v))+ h/hc+ pc/p
ifelse(critL>0, pa[1] > 1, pa[1] <- exp(critL))
# Troca ponto atual por candidato aceito #
if(pa[1] > runif(1))
{
l <- lc
npa[1] <- npa[1] +1
v <- vc
}
# Atualizando T (razão toxicante) #
# Verossimilhança #
Tc[1] <- runif(1)
sc[1] <- 1
for (i in 2:d)
{
Tc[i] <- runif(1,0,Tc[(i-1)])
sc[i] <- Tc[(i-1)]
}
pno <- exp(-l
*
z
*
(1-T))
pnoc <- exp(-l
*
z
*
(1-Tc))
v <- (1-pno)^y
*
pno^(n-y)
vc <- (1-pnoc)^y
*
pnoc^(n-y)
critT <- exp(sum(log(vc)-log(v)+(1/s)/(1/sc)
ifelse(critT>0, pa[2] > 1, pa[2] <- exp(critT))
if(critT > runif(1))
{
T <- Tc
npa[2] <- npa[2] +1
89
}
f <- gamma(n+1)/(gamma(y+1)
*
gamma(n-y+1))
pno <- exp(-l
*
z
*
(1-T))
v <- f
*
(1-pno)^y
*
pno^(n-y)
vero <- exp(sum(log(v)))
# Armazena na cadeia atual M valores da posteriori #
if(iter <= M)
{
a[iter] <- l
}
if(iter>B && (iter-B)%%J==0)
{
a[2:M] <- a[1:(M-1)]
a[1] <- l
mL <- mean(a)
vL <- var(a)+ .1
bcL <- vL/mL
acL <- mL/bcL
write(t(c(l,T,vero)),file=" arquivo.txt" ,q, append=TRUE)
}
}
90
PROGRAMA 4 B Rotina para gerar cadeias via MCMC do modelo de Completo.
#Definindo delineamento#
y <- c(0,0,3,3,2)
r <- 3
d <-length(y)
n <-c(rep(r,length(y)))
p <-length(y)+1
z <- 10^(1-c(2:p))
T <- 1-sort(runif(d))
Tc <- T
s <- c(1,T[1:(d-1)])
sc <- s
q <- 4
M <- 1000
Nef <- 10000
# queima(B), salto(J) e tamanho da amostra(N)
B <- 50000
J <- 200
N <- B + Nef
*
J
# parâmetros iniciais para as candidatas #
alphaLI <- 2
betaLI <- 100
alphaL <- 2
betaL <- 100
bcLI <- 2
acLI <- 100
bcL <- 2
acL <- 100
a <- rnorm(M,alphaLI
*
betaLI,sqrt(alphaLI
*
betaLI^ 2))
b <- rnorm(M,alphaL
*
betaL,sqrt(alphaL
*
betaL^ 2))
# Chutes iniciais para candidatas
l <- 1:2
lc <- 1:2
lc[1] <- rgamma(1,acLI,1/bcLI) # chute para lambda_I
lc[2] <- rgamma(1,acL,1/bcL) # chute para lambda
#Inicia a saida aceitando o chute
pa <- 1:3
npa <- 0
*
c(1:3)
v <- 1
91
vc <- 1
vero <- 1
for(iter in 1:N)
{
# Atualizando lambda_I #
lc[1] <- rgamma(1,acLI,1/bcLI)
#Definindo partes da verossimilhança
f <- gamma(n+1)/(gamma(y+1)
*
gamma(n-y+1))
po <- exp(-l[1]
*
z
*
(1-T))
*
(1-exp(-l[2]
*
z
*
(1-T)))
poc <- exp(-lc[1]
*
z
*
(1-T))
*
(1-exp(-l[2]
*
z
*
(1-T)))
# Verossimilhança
v <- po^y
*
(1 - po)^(n-y)
vc <- poc^y
*
(1 - poc)^(n-y)
# Razão das probabilidades nas candidatas
pc <- dgamma(lc[1],alphaLI,1/betaLI)
p <- dgamma(l[1] ,alphaLI,1/betaLI)
hc <- dgamma(lc[1],acLI,1/bcLI)
h <- dgamma(l[1] ,acLI,1/bcLI)
critLI <- sum(log(vc)-log(v))+ h/hc + pc/p
ifelse(critLI>0, pa[1] > 1, pa[1] <- exp(critLI))
# Troca ponto atual por candidato aceito #
if(pa[1] > runif(1))
{
l[1] <- lc[1]
npa[1] <- npa[1]+1
}
# Atualizando Lambda #
lc[2] <- rgamma(1,acL,1/bcL)
# Definindo partes da verossimilhança #
f <- gamma(n+1)/(gamma(y+1)
*
gamma(n-y+1))
po <- exp(-l[1]
*
z
*
(1-T))
*
(1-exp(-l[2]
*
z
*
(1-T)))
poc <- exp(-l[1]
*
z
*
(1-T))
*
(1-exp(-lc[2]
*
z
*
(1-T)))
# Verossimilhança
v <- po^y
*
(1 - po)^(n-y)
vc <- poc^y
*
(1 - poc)^(n-y)
# Razão das probabilidades nas candidatas
pc <- dgamma(lc[2],alphaL,1/betaL)
p <- dgamma(l[2] ,alphaL,1/betaL)
hc <- dgamma(lc[2],acL,1/bcL)
92
h <- dgamma(l[2] ,acL,1/bcL)
critL <- sum(log(vc)-log(v))+log(h)-log(hc)+log(pc)-log(p)
ifelse(critL>0, pa[2] > 1, pa[2] <- exp(critL))
# Troca ponto atual por candidato aceito #
if(pa[2] > runif(1))
{
l[2] <- lc[2]
npa[2] <- npa[2]+1
}
# Atualizando c (razão toxicante) #
# Verossimilhança #
Tc[1] <- runif(1)
sc[1] <- 1
for (i in 2:d)
{
Tc[i] <- runif(1,0,Tc[(i-1)])
sc[i] <- Tc[(i-1)]
}
#Definindo partes da verossimilhança #
f <- gamma(n+1)/(gamma(y+1)
*
gamma(n-y+1))
po <- exp(-l[1]
*
z
*
(1-T))
*
(1-exp(-l[2]
*
z
*
(1-T)))
poc <- exp(-lc[1]
*
z
*
(1-Tc))
*
(1-exp(-l[2]
*
z
*
(1-Tc)))
# Verossimilhança
v <- po ^y
*
(1 - po)^(n-y)
vc <- poc^y
*
(1 - poc)^(n-y)
# Razão das probabilidades nas candidatas
critT <- exp(sum(log(vc)-log(v)+log(1/s)-log(1/sc)))
ifelse(critT>0, pa[3] > 1, pa[3] <- exp(critT))
# Troca ponto atual por candidato aceito #
if(pa[3] > runif(1))
{
T <- Tc
npa[3] <- npa[3]+1
}
f <- gamma(n+1)/(gamma(y+1)
*
gamma(n-y+1))
po <- exp(-l[1]
*
z
*
(1-T))
*
(1-exp(-l[2]
*
z
*
(1-T)))
v <- f
*
po^y
*
(1 - po)^(n-y)
vero <- exp(sum(log(v)))
# Armazena na cadeia atual #
if(iter<=M)
93
{
a[iter] <- l[1]
b[iter] <- l[2]
}
if(iter>B && (iter-B)%%J==0)
{
a[2:M] <- a[1:(M-1)]
b[2:M] <- b[1:(M-1)]
a[1] <- l[1]
b[1] <- l[2]
mLI <- mean(a)
vLI <- var(a)+.1
mL <- mean(b)
vL <- var(b)+.1
bcLI <- vLI/mLI
acLI <- mLI/bcLI
bcL <- vL/mL
acL <- mL/bcL
write(t(c(l[1],l[2],T,vero)),file= "arquivo",q, append=TRUE)
}
}
94
ANEXOS C
TABELA DE NÚMERO MAIS PROVÁVEL Páginas
TABELA 1C. Tabela de NMP para 3 tubos utilizando-se as diluições 0,1; 0,01 e
0,001 e intervalo de confiança de 95%. . . . . . . . . . . . . . . . . . . . . . 96
95
TABELA 1C.Tabela de NMP para 3 tubos utilizando-se as diluições 0,1; 0,01 e
0,001 e intervalo de confiança de 95%.
Tubos I.C. 95 %
0,1 0,01 0,001 NMP LI LS
0 0 0 <3.0 - 9.5
0 0 1 3.0 0.15 9.6
0 1 0 3.0 0.15 11
0 1 1 6.1 1.2 18
0 2 0 6.2 1.2 18
0 3 0 9.4 3.6 38
1 0 0 3.6 0.17 18
1 0 1 7.2 1.3 18
1 0 2 11 3.6 38
1 1 0 7.4 1.3 20
1 1 1 11 3.6 38
1 2 0 11 3.6 42
1 2 1 15 4.5 42
1 3 0 16 4.5 42
2 0 0 9.2 1.4 38
2 0 1 14 3.6 42
2 0 2 20 4.5 42
2 1 0 15 3.7 42
2 1 1 20 4.5 42
2 1 2 27 8.7 94
Fonte: Bacteriological Analytical Manual, 2001
...continua...
96
Tabela 1C: Cont.
Tubos I.C. 95 %
0,1 0,01 0,001 NMP LI LS
2 2 0 21 4.5 42
2 2 1 28 8.7 94
2 2 2 35 8.7 94
2 3 0 29 8.7 94
2 3 1 36 8.7 94
3 0 0 23 4.6 94
3 0 1 38 8,7 110
3 0 2 64 17 180
3 1 0 43 9 180
3 1 1 75 17 200
3 1 2 120 37 420
3 1 3 160 40 420
3 2 0 93 18 420
3 2 1 150 37 420
3 2 2 210 40 430
3 2 3 290 90 1000
3 3 0 240 42 1000
3 3 1 460 90 2000
3 3 2 1100 180 4100
3 3 3 >1100 420 -
Fonte: Bacteriological Analytical Manual, 2001
97
Livros Grátis
( http://www.livrosgratis.com.br )
Milhares de Livros para Download:
Baixar livros de Administração
Baixar livros de Agronomia
Baixar livros de Arquitetura
Baixar livros de Artes
Baixar livros de Astronomia
Baixar livros de Biologia Geral
Baixar livros de Ciência da Computação
Baixar livros de Ciência da Informação
Baixar livros de Ciência Política
Baixar livros de Ciências da Saúde
Baixar livros de Comunicação
Baixar livros do Conselho Nacional de Educação - CNE
Baixar livros de Defesa civil
Baixar livros de Direito
Baixar livros de Direitos humanos
Baixar livros de Economia
Baixar livros de Economia Doméstica
Baixar livros de Educação
Baixar livros de Educação - Trânsito
Baixar livros de Educação Física
Baixar livros de Engenharia Aeroespacial
Baixar livros de Farmácia
Baixar livros de Filosofia
Baixar livros de Física
Baixar livros de Geociências
Baixar livros de Geografia
Baixar livros de História
Baixar livros de Línguas
Baixar livros de Literatura
Baixar livros de Literatura de Cordel
Baixar livros de Literatura Infantil
Baixar livros de Matemática
Baixar livros de Medicina
Baixar livros de Medicina Veterinária
Baixar livros de Meio Ambiente
Baixar livros de Meteorologia
Baixar Monografias e TCC
Baixar livros Multidisciplinar
Baixar livros de Música
Baixar livros de Psicologia
Baixar livros de Química
Baixar livros de Saúde Coletiva
Baixar livros de Serviço Social
Baixar livros de Sociologia
Baixar livros de Teologia
Baixar livros de Trabalho
Baixar livros de Turismo