Um Esquema Eficiente de Amostragem em Modelos Dinâmicos Generalizados com Aplicações em Funções de T por Romy Elena Rodríguez Ravines - Versão HTML

/ Inicio / Estatística / Um Esquema Eficiente de Amostragem em Modelos Dinâmicos Generalizados com Aplicações em Funções de T

ATENÇÃO: Esta é apenas uma visualização em HTML e alguns elementos como links e números de página podem estar incorretos.
Faça o download do livro em PDF, ePub para obter uma versão completa.

0.98

1.00

0.4

0.5

0.6

0.7

0.8

0.9

0.2

0.4

0.6

0.8

1.0

(a) k - KT

(b) λ - KT

3.5

3.0

2.5

2.0

1.5

1.0

0.5

0.0

0.90

0.92

0.94

0.96

0.98

1.00

0.3

0.4

0.5

0.6

0.7

0.8

0.9

0.4

0.5

0.6

0.7

0.8

0.9

1.0

(d) k - TF

(e) λ - TF

(f) ρ - TF

Figura 2.3: Amostras a posteriori dos parâmetros k e λ sob o modelo (2.15) com Erro II. A linhas pontilhadas correspondem à média a posteriori. (KT = usando a transformac¸˜ao de Koyck, TF usando

a parametrizac¸˜ao como func¸˜ao de transferência).

à escolha da distribuic¸˜ao a priori e testar a hip ótese ρ = λ.

2.5.4

Comparaç˜ao com uma Abordagem Clássica

Com a finalidade de comparar os resultados das sec¸ ˜oes anteriores com os de uma abordagem

clássica, seguimos uma estratégia bastante utilizada para trabalhar com séries de tempo n˜ao esta-

cionárias. Especificamente seguimos Pesaran & Shin (1999) e Bentzen & Engsted (2001) que mostraram que a aproximac¸˜ao com modelos ARDL é perfeitamente válida quando (algumas das)

as variáveis s˜ao n˜ao-estacionárias integradas de ordem 1, I(1), com a única restric¸˜ao de que as

Tabela 2.3: Critérios de Comparaç˜ao de Modelos: DIC (Deviance Information Criteria) e EPD (Expected

Predictive Deviance)

Transformac¸˜ao de Koyck (KT)

Func¸˜ao de Transferência (FT)

Modelo (2.15)

DIC

312,49

311,40

com Erro I

EPD

1455,57

1453,29

Modelo (2.15)

DIC

280,37

279,20

com Erro II

EPD

760,16

767,91

variáveis I(1) sejam cointegradas e exista apenas uma equac¸˜ao de co-integrac¸˜ao.

Usando a transformac¸˜ao de Koyck, o modelo em (2.15) pode ser visto como um modelo

ARDL(1, 0), assim k e k(1 − λ) podem ser interpretados como elasticidades de longo e curto pra-

zos, respectivamente. No contexto bayesiano, a inferência sobre func¸ ˜oes dos parâmetros pode ser

realizada, de forma direta, a partir das amostras obtidas das distribuic¸ ˜oes a posteriori dos mes-

mos. Neste caso, basta usar as amostras do MCMC de k e λ para ter uma amostra a posteriori de

k(1−λ). A primeira linha da tabela 2.4 mostra as elasticidades de longo prazo estimadas com o mo-

delo selecionado na sec¸˜ao 2.5.3 (com Erro II). Ela também apresenta a elasticidade de longo prazo

obtida com a representac¸˜ao de func¸˜ao de transferência (com Erro II) incluindo um termo constante

(intercepto) no modelo (segunda linha), com o procedimento de m´ınimos quadrados ordinários

(OLS) 1 e com o método de co-integrac¸˜ao de Johansen 2. Observamos que as estimativas pontuais para o parâmetro de longo prazo obtidas com todos os modelos s˜ao bastante similares. A figura

2.4(a) apresenta o histograma da amostra da posteriori de k obtida com o MCMC e a distribuic¸˜ao

normal obtida para ˆk com o procedimento OLS. A figura 2.4(b) mostra as distribuic¸ ˜oes emp´ırica e

assint ótica obtidas para k(1 − λ). Embora as estimativas pontuais para k sob ambos procedimentos

sejam similares, na figura observamos que a distribuic¸˜ao assint ótica é muito diferente da distribui-

c¸˜ao a posteriori. Notamos que na distribuic¸˜ao assint ótica P(ˆk(H) ≥ 0.905) ≥ 0.50 enquanto que na

distribuic¸˜ao a posteriori, P(k ≥ 0.905|H) = 0.025, onde H representa a base de dados completa. Este

fato pode ser explicado porque o último resultado é baseado unicamente na amostra observada.

1sugerido por Engle & Granger (1987)

2sugerido por Johansen (1991)

Tabela 2.4: Elasticidades de longo e curto prazo

Model

Longo prazo

Curto prazo

média

d.p.

média

d.p.

Func¸˜ao de Transferência a

0,959

0,018

0,251

0,116

Func¸˜ao de Transferência b

0,910

0,042

0,262

0,108

OLS

0,878

0,014

0,266

0,114

Johansen (VECM) c

0,884

0,022

–

a Modelo em (2.16).

b Modelo em (2.16) incluindo um termo constante ou intercepto.

c Considerando uma tendência linear nos dados e 0 defasagens no VAR, este modelo foi escolhido seguindo o critério de Schwartz

3.5

3.0

2.5

2.0

1.5

1.0

0.5

0.0

0.70

0.80

0.90

1.00

−0.2

0.0

0.2

0.4

0.6

0.8

k(1 − λ)

(a) Longo Prazo

(b) Curto Prazo

Figura 2.4: Elasticidades de longo e curto prazo estimadas com o modelo em (2.16) com Erro II e um intercepto na equaç˜ao das observaç ˜oes. Em ambos painéis, os histogramas correspondem à

amostra da distribuic¸˜ao a posteriori, a linha vertical corresponde à média a posteriori e as linhas

pontilhadas a distribuic¸˜ao normal assint ótica estimada via OLS.

2.5.5

Bondade de Ajuste e Previs˜ao

Com o objetivo de examinar a bondade de ajuste e a capacidade preditiva de nossa abordagem, re-

estimamos os parâmetros do modelo em (2.16) com Erro II, incluindo um intercepto na primeira

equac¸˜ao, utilizando somente as primeiras 49 observac¸ ˜oes. A figura 2.5(a) mostra o intervalo de

95% de credibilidade estimado para a amostra entre 1947.I and 1959.II. Observamos um bom

ajuste, pois somente duas das observac¸ ˜oes ficaram fora do intervalo obtido. A figura 2.5(b) apre-

senta o intervalo de 95% de credibilidade estimado para as observac¸ ˜oes deixadas fora da amostra.

Observamos que todos os valores reais est˜ao contidos no intervalo estimado.

305

280

260

300

Y 240

Y 295

220

290

200

285

1947 1949 1951 1953 1955 1957 1959

1959.2

1959.5

1959.8

1960.1

1960.4

(a) Valores replicados (1947.I-1959.II)

(b) Previs˜ao (1959.III-1960.IV)

Figura 2.5: Valores replicados e previs˜ao do consumo obtido com o modelo em (2.16) com Erro II e incluindo um intercepto. Em ambos painéis a linha s ólida corresponde à média, as linhas pon-

tilhadas ao primeiro e terceiro quartis da distribuic¸˜ao a posteriori. As cruzes correspondem aos

valores observados.

2.5.6

Tratamento da Endogeneidade

Até agora temos considerado a renda como variável predeterminada somente com o objetivo de

ilustrar nossa proposta e compará-la com Zellner & Geisel (1970); contudo, é bastante conhecida

a interdependência presente num modelo macro-econométrico de renda e consumo. Para lidar

com esta situac¸˜ao, nesta sec¸˜ao adotamos a abordagem de Congdon (2001, exemplo 7.38). Primeiro

ajustamos um modelo AR(1) para X com a finalidade de formar uma variável instrumental e depois

estimamos o modelo em (2.16) com Erro II.

A tabela 2.5 apresenta os resultados obtidos, incluindo as estimativas para o coeficiente auto-

regressivo do modelo AR(1), α. Basicamente, comparamos as estat´ısticas a posteriori para k e

k(1 − λ) com as apresentadas na segunda linha da tabela 2.4. Observamos que, como era esperado,

as médias obtidas s˜ao menores ao usar variáveis instrumentais enquanto que os desvios padr ˜oes

s˜ao maiores.

Tabela 2.5: Estat´ısticas a posteriori associadas com os parâmetros do modelo em (2.16) com Erro II e modelo AR(1) para a renda

média

d.p.

2, 5%

50%

97, 5%

0,902

0,102

0,606

0,922

0,995

1,001

0,768

0,145

0,407

0,808

0,952

1,003

k(1 − λ)

0,211

0,135

0,029

0,172

0,548

1,003

0,753

0,108

0,549

0,748

0,973

1,002

σ2

9,819

2,011

6,587

9,602

14,33

1,001

0,998

0,013

0,972

0,998

1,024

1,003

b ˆ

R é o fator de reduc¸˜ao de escala (na convergência, ˆR = 1).

2.6

Consideraç ˜oes Finais

Embora os primeiros modelos de defasagens distribu´ıdas tenham sido propostos durante a década

de 50, eles s˜ao bastante úteis e muito utilizados na atualidade. Certamente esta classe de modelos

envolve vários t ópicos que merecem atenc¸˜ao: a forma funcional da relac¸˜ao entre os seus coefi-

cientes, a estrutura considerada para o termo de erro, a n˜ao estacionariedade das séries temporais,

variáveis end ógenas, dados de painel, etc.

Neste cap´ıtulo discutimos a inferência dos modelos de defasagens distribu´ıdas do ponto de

vista bayesiano. Especificamente, mostramos que estes modelos podem ser representados dentro

de uma classe muito ampla: os modelos lineares dinâmicos. Este fato permite-nos usar todas as

vantagens oferecidas pela abordagem bayesiana, como a incorporac¸˜ao de informac¸˜ao subjetiva e a

descric¸˜ao da incerteza, através da distribuic¸˜ao a posteriori dos parâmetros de interesse.

A inferência sobre modelos de defasagens distribu´ıdas pode ser realizada com métodos MCMC.

Neste cap´ıtulo mostramos que esta tarefa é bastante simplificada com o uso do WinBUGS, um pa-

cote de livre disponibilidade. Mostramos também que a comparac¸˜ao de modelos pode ser rea-

lizada utilizando o DIC, critério que está implementado no WinBUGS, e o EPD. Além disso, as

cadeias obtidas através do WinBUGS podem ser exploradas com a biblioteca CODA (Best et al.,

1997), que também é de livre disponibilidade, e é útil para obter indicadores (diagn ósticos) de

convergência das cadeias.

Na verdade, outros algoritmos particulares podem ser utilizados, mas o ganho na eficiência

obtido com eles, é compensado pelo WinBUGS devido à facilidade de implementac¸˜ao dos c ódigos.

Este fato permite-nos explorar diferentes alternativas de modelagem, avaliar o efeito da escolha

de prioris, o uso de variáveis transformadas ou a inclus˜ao de termos estocásticos na equac¸˜ao do

sistema, como foi apresentado aqui.

Modelos mais complexos podem ser analisados do ponto de vista bayesiano. Por exemplo,

Migon & Harrison (1985) mostram um caso n˜ao normal e n˜ao linear usando estimac¸˜ao Linear Bayes, e Alves (2005) prop ˜oe um esquema MCMC neste mesmo contexto.

Um t ópico de importância te órica, que é parte do nosso trabalho atual, é a escolha da forma

funcional da func¸˜ao de transferência e como isto afeta as previs ˜oes e incertezas associadas. Uma

possibilidade é o uso de uma func¸˜ao de transferência de segunda ordem.

2.7

Apêndice: Exemplo de c ódigo usado em WinBUGS

Este é o c ódigo utilizado para ajustar o modelo de Koyck em (2.12) com erros n˜ao correlacionados

e homoscedásticos, isto é ( ind

N(0, σ2)).

∼

model #(Ct=consumo, Yt=renda){

e[1]

<- 0

for(t in 2:T){

mean.Y[t] <- lambda*Y[t-1] + k*(1-lambda)*X[t] - lambda*e[t-1]

Y[t]

dnorm(mean.Y[t], tau.y)

Y.hat[t]

dnorm(mean.Y[t], tau.y)

e[t]

<- Y[t] - mean.Y[t]

}

# Prior

lambda

˜ dbeta(1,1)

tau.y

˜ dgamma(0.0001,0.0001)

var.y

<- 1/tau.y

phi

<- k*(1-lambda)

}

Cap´ıtulo 3

UM ESQUEMA DE AMOSTRAGEM EFICIENTE PARA MODELOS DIN ˆ

AMICOS

GENERALIZADOS

Neste cap´ıtulo apresentamos um esquema de amostragem m últipla dos parâmetros de estado

de modelos dinâmicos n˜ao lineares, dentro de um amostrador de Gibbs. O esquema proposto

combina a abordagem do Conjugate Updating para modelos lineares generalizados dinâmicos com

o Backward Sampling utilizado nos modelos dinâmicos normais para amostrar do vetor de estados.

Denominamos este esquema de Conjugate Updating Backward Sampling (CUBS). Observamos que

o CUBS reduz significativamente o tempo computacional necessário para atingir a convergência

de todas as cadeias dos modelos, sendo simples de implementar. As amostras das distribuic¸ ˜oes a

posteriori dos outros parâmetros dos modelos podem ser obtidas com outros algoritmos. Apresen-

tamos os resultados de um estudo de Monte Carlo extensivo realizado para comparar a eficiência

do CUBS versus outros algoritmos previamente propostos na literatura bayesiana e, também, duas

aplicac¸ ˜oes com dados reais.

3.1

Introduç˜ao

Os modelos dinâmicos n˜ao lineares e n˜ao gaussianos s˜ao amplamente utilizados na estat´ıstica

aplicada, em particular, sob a abordagem bayesiana. Neste cap´ıtulo consideramos que um modelo

dinâmico de tempo discreto definido na fam´ılia exponencial é tal que

yt|ηt, φ ∼ exp[φ{ytηt − a(ηt)}]b(yt, φ),

t = 1, 2, . . . .

(3.1a)

g(µt) = Ft(ψ ) θ

(3.1b)

θt = Gt(ψ )θ

t−1 + wt,

wt ∼ N(0, W)

(3.1c)

θ0|Y0 ∼ N(m0, C0),

onde os vetores θt s˜ao denominados parâmetros de estado e est˜ao relacionados através do tempo

via (3.1c), a equac¸˜ao do sistema; θ0|Y0 é a informac¸˜ao dispon´ıvel anterior ao per´ıodo em que se tem informac¸˜ao através das observac¸ ˜oes; ηt é o parâmetro natural e E[yt|ηt] = µt = ˙a(ηt) está relacionado

a θt via uma func¸˜ao de ligac¸˜ao conhecida, g(·); ψ e ψ denotam todas as quantidades desconheci-

das envolvidas na definic¸˜ao de Ft(ψ ) e G

), respectivamente. Também, w

t(ψ2

t é o termo de erro

do sistema com variância (de evoluc¸˜ao) W. Freq üentemente W é uma matriz diagonal e pode

variar ao longo do tempo. Do ponto de vista bayesiano, a especificac¸˜ao do modelo fica completa

depois da atribuic¸˜ao da distribuic¸˜ao a priori dos parâmetros. A distribuic¸˜ao a posteriori resultante

para os parâmetros em (3.1) n˜ao tem uma forma fechada conhecida, por esta raz˜ao a inferência

n˜ao pode ser realizada de forma anal´ıtica. Logo, métodos de simulac¸˜ao estocástica s˜ao necessários

para gerar as amostras de interesse, mas a alta correlac¸˜ao existente entre os θts pode dificultar a

gerac¸˜ao de amostras independentes. No contexto dos métodos MCMC, nosso objetivo é propor

um algoritmo com passo de Metropolis-Hastings independente e eficiente (pouca demanda de

tempo computacional e fácil de implementar) para amostrar da distribuic¸˜ao condicional completa

de θt, t = 1, . . . , T.

Antes do uso intensivo das técnicas MCMC, West et al. (1985), Kitagawa (1987) e Fahrmeir

(1992), entre outros, propuseram diferentes métodos para obter aproximac¸ ˜oes da distribuic¸˜ao a

posteriori em quest˜ao. Desde a década passada, como pode ser visto em Migon et al. (2005), vários

algoritmos MCMC têm sido propostos para amostrar dos estados dos modelos dinâmicos. Em

particular, Gamerman (1998) e Geweke & Tanizaki (2001) consideraram abordagens bayesianas com passos de Metropolis-Hastings dentro do amostrador de Gibbs. Dado que o ponto cr´ıtico de

um algoritmo de Metropolis-Hastings é a escolha da distribuic¸˜ao proposta, eles est˜ao interessados

na construc¸˜ao de propostas eficientes. Por exemplo, Gamerman (1998) sugeriu o uso de um mo-

delo dinâmico normal ajustado com a finalidade de construir uma boa densidade proposta. Esta

idéia pode ser implementada de três maneiras diferentes: na amostragem individual ou em blocos

dos estados, ou na amostragem individual dos termos de erro da equac¸˜ao de evoluc¸˜ao. Gamer-

man concluiu que os movimentos individuais s˜ao prefer´ıveis aos movimentos em blocos e que

amostrar dos erros é mais eficiente porque as cadeias resultantes s˜ao menos autocorrelacionadas,

o que ajuda a acelerar a convergência do algoritmo. Entretanto, esta proposta requer um esforc¸o

computacional bastante significativo pois o c ódigo, geralmente, é complicado de escrever e pode

levar “muito” tempo para completar uma única iterac¸˜ao. Por outro lado, Geweke & Tanizaki (2001)

propuseram várias formas de construir densidades propostas para amostrar os estados em passos

individuais e realizaram um estudo de Monte Carlo para avaliar a sensibilidade dos resultados

a esta escolha. Eles conclu´ıram que os resultados s˜ao bastante robustos à escolha da densidade

proposta para o passo de Metropolis-Hastings.

Nossa alternativa é amostrar do vetor de estados {θ1, . . . , θT} de (3.1), em bloco, isto é, usando

um movimento m últiplo num esquema análogo ao FFBS (do inglês Forward Filtering and Backward

Sampling), de Fr ühwirth-Schnater (1994) e de Carter & Kohn (1994) proposto para os modelos lineares dinâmicos normais. Ambos trabalhos conclu´ıram que amostrar todos os estados em forma

simultânea acelera a convergência das cadeias em comparac¸˜ao ao caso de movimentos indivi-

duais. A idéia principal do FFBS é o uso da decomposic¸˜ao da distribuic¸˜ao a posteriori conjunta

dos estados, no produt ório das distribuic¸ ˜oes retrospectivas (ou suavizadas). No caso gaussiano,

as quantidades necessárias para a especificac¸˜ao dos momentos das distribuic¸ ˜oes suavizadas s˜ao

obtidas a partir do Filtro de Kalman. Logo, para obter amostras da distribuic¸˜ao condicional com-

pleta conjunta, primeiro calculam-se os momentos “para frente” (on-line) ou passo de filtragem, e

depois se amostra das distribuic¸ ˜oes normais suavizadas (que corresponde ao passo de Backward

Sampling).

Neste cap´ıtulo, nossa proposta para modelos n˜ao normais é substituir o primeiro passo do FFBS

(o passo de filtragem) pelo Conjugate Updating apresentado em West et al. (1985), para modelos di-

nâmicos na fam´ılia exponencial. Em outras palavras, nossa proposta consiste em aproximar as

distribuic¸ ˜oes on-line introduzindo uma priori conjugada para o parâmetro can ônico da distribui-

c¸˜ao da observac¸˜ao, e depois usar um passo de amostragem retrospectiva como no caso gaussiano

(o passo Backward Sampling). Denominamos este algoritmo CUBS, abreviac¸˜ao do inglês Conjugate

Updating Backward Sampling. A grande forc¸a deste método reside em considerar a estrutura de

correlac¸˜ao temporal dos θt. As principais vantagens s˜ao: facilidade de implementac¸˜ao e reduc¸˜ao

do tempo e esforc¸o computacional quando comparado com outros algoritmos. O CUBS pode ser

usado dentro de outros algoritmos MCMC e, portanto, as distribuic¸ ˜oes a posteriori dos demais

parâmetros do modelo em (3.1) podem ser amostradas de maneira usual.

Este cap´ıtulo está organizado da seguinte forma. Na sec¸˜ao 3.2 descrevemos o algoritmo MCMC

proposto, destacando as vantagens de sua utilizac¸˜ao. Na sec¸˜ao 3.3 comparamos os resultados

obtidos com o CUBS com outros algoritmos previamente propostos na literatura bayesiana. A

comparac¸˜ao é feita com base num estudo de simulac¸˜ao extensivo realizado com dados artificiais.

Também apresentamos os resultados do uso do CUBS em duas aplicac¸ ˜oes a dados reais. Final-

mente, apresentamos algumas conclus ˜oes e projetos futuros na sec¸˜ao 3.4.

3.2

Esquema de Amostragem Proposto

Nesta sec¸˜ao propomos um esquema de amostragem MCMC baseado na combinac¸˜ao da filtragem

conjugada (Conjugate Updating) com a amostragem retrospectiva (Backward Smoothing). Iniciamos

esta sec¸˜ao revisando os principais aspectos da filtragem conjugada e a amostragem retrospec-

tiva introduzidos, respectivamente, no contexto dos modelos lineares generalizados dinâmicos

(DGLM) e modelos lineares dinâmicos (DLM) normais. Depois apresentamos nossa proposta, o al-

goritmo CUBS, que combina as idéias descritas anteriormente. Terminamos esta sec¸˜ao ressaltando

as vantagens do CUBS. Para definir a nossa notac¸˜ao, daqui em diante Θ = {θ1, . . . , θT} representa o

vetor de estados, Yt = {y1, . . . , yt} representa a informac¸˜ao até o tempo t e Φ = (φ, ψ , ψ , W) denota

o vetor de todas as quantidades desconhecidas em (3.1), exceto Θ. De acordo com a especificac¸˜ao

em (3.1) e assumindo distribuic¸ ˜oes a priori independentes para φ, Ψ1, Ψ2, W, e para a informac¸˜ao inicial θ0, a distribuic¸˜ao a posteriori é

p(θ1, . . . , θT, θ0, ψ , ψ , W, φ|y

p(y

1, . . . , yT) ∝

t|θt, ηt, φ)

t=1

p(θt|θt−1, Ft(ψ ), G

), W)p(θ

)p(ψ )p(W)p(φ). (3.2)

t(ψ2

0)p(ψ1

Como a distribuic¸˜ao em (3.2) n˜ao tem forma fechada conhecida, propomos usar métodos MCMC

para amostrá-la. Nosso interesse é amostrar da distribuic¸˜ao condicional completa de {θ1, . . . , θT}.

Especificamente, usamos um algoritmo h´ıbrido, um amostrador de Gibbs com passos de Metropolis-

Hastings. Nosso foco é a construc¸˜ao de uma distribuic¸˜ao proposta eficiente.

Nos anos 80, antes do uso dos métodos MCMC, West et al. (1985) propuseram um algoritmo

seq üencial para aproximar a distribuic¸˜ao a posteriori de um modelo linear generalizado dinâmico

(DGLM). Esse algoritmo é denominado Conjugate Updating e baseia-se nos procedimentos de Linear

Bayes, uma aproximac¸˜ao realizada ao n´ıvel das distribuic¸ ˜oes a priori. Nos DGLM, a priori dos esta-

dos é especificada parcialmente através dos seus primeiro e segundo momentos, digamos at e Rt.

Estes s˜ao utilizados para determinar rt e st, os parâmetros de uma distribuic¸˜ao a priori conjugada

apropriada para ηt, isto é, ηt|Yt−1 ∼ CP(rt, st). Estes momentos s˜ao determinados resolvendo um

conjunto relativamente simples de duas equac¸ ˜oes tal como é mostrado no apêndice (ver tabela 3.6).

Os parâmetros da posteriori conjugada para o parâmetro natural s˜ao obtidos facilmente depois que

yt é observado. A informac¸˜ao nesta posteriori é repassada para a distribuic¸˜ao a posteriori dos es-

tados. Como esta distribuic¸˜ao n˜ao é completamente conhecida, West et al. (1985) propuseram o

uso de estimac¸˜ao via Linear Bayes, considerando uma func¸˜ao perda quadrática para derivar uma

aproximac¸˜ao para a distribuic¸˜ao a posteriori do vetor de estados. A aproximac¸˜ao consiste em

p(θt|Yt, Φ) ∝ p(θt|Yt−1, Φ)p(yt|θt, Φ)

p(θt|ηt, Yt−1, Φ)p(ηt|Yt−1, Φ)p(yt|ηt, Φ)dηt,

onde ηt é o parâmetro canônico de p(yt|ηt, Φ) no tempo t. Como o primeiro termo na integral acima

é parcialmente conhecido, usamos o principio de Linear Bayes para obter

p(θt|Yt, Φ) ∝

p(θt | ηt, Yt−1, Φ)p(ηt|Yt, Φ)dηt

[mt, Ct],

onde [mt, Ct] indica que p(θt|Yt, Φ) é parcialmente conhecida através dos seus dois primeiros mo-

mentos, dados por

mt = at + RtFt( f ∗ −

ft) qt

q∗ 1

= Rt − RtFtF

tRt 1 − qt qt

onde at = Gtmt−1, Rt = GtCt−1G +

Wt e ( ft, qt) e ( f ∗ q∗) s˜ao, respectivamente, a média e variância,

a priori e a posteriori de ηt. Vale a pena destacar que o conhecido Filtro de Kalman é recuperado

sob a hip ótese de normalidade para yt e a func¸˜ao identidade como func¸˜ao de ligac¸˜ao.

O quadro 1 apresenta um esquema da análise seq üencial, via estimac¸˜ao Linear Bayes, para

modelos lineares dinâmicos generalizados, em paralelo com a análise seq üencial, via Filtro de

Kalman, para modelos lineares dinâmicos normais. Observamos que ηt faz o papel de uma variável

latente que facilita a aproximac¸˜ao de p(θt|Dt). Por conveniência assumimos que ηt|Dt−1 segue uma

distribuic¸˜ao conjugada com p(yt|ηt), mas esta hipótese pode ser relaxada.

Quadro 1: Análise Seq üencial

• Nos modelos lineares dinâmicos normais:

Evoluc¸˜ao

. . .

Atualizac˜ao

(θ

t−1|Yt−1)

(θt|Yt−1)

(θt|Yt) . . .

• Nos modelos lineares generalizados dinâmicos:

Evoluc¸˜ao

. . . (θ

t−1|Yt−1)

(θt|Yt−1)

(θt|Yt) . . .

Atualizac˜ao

(η

t|Yt−1)

(ηt|Yt)

No contexto dos modelos lineares dinâmicos normais é conhecido que as implementac¸ ˜oes de

movimentos individuais nos métodos MCMC s˜ao bastante ineficientes devido à alta correlac¸˜ao

temporal entre os parâmetros de estado. Para contornar esta situac¸˜ao, Fr ühwirth-Schnater (1994)

e Carter & Kohn (1994) propuseram amostrar todos os elementos de Θ em (3.1) num passo de amostragem m últipla. A abordagem proposta por eles está baseada na conhecida decomposic¸˜ao

da distribuic¸˜ao condicional completa conjunta dos estados:

T−1

p(Θ|YT, Φ) = p(θT|YT, Φ)

p(θt|θt+1, Yt, Φ).

(3.4)

t=1

Como p(θt|θt+1, Yt, Φ) ∝ p(θt+1|θt, Yt, Φ)p(θt|Yt, Φ), é fácil mostrar que p(θt|θt+1, Yt, Φ) é N(ms, Cs t

t ),

onde

ms = m

+ W

t + CtGt(GtCtGt

t)−1(θt+1 − Gtmt)

Cs = C

+ W

t − CtGt(GtCtGt

t)−1GCt,

onde mt e Ct s˜ao os primeiros e segundos momentos obtidos fazendo uso do Filtro de Kalman.

Daqui, o algoritmo FFBS consiste principalmente em amostrar os elementos de Θ de forma seq üen-

cial.

Com base nas idéias descritas acima, descrevemos agora como obter uma boa distribuic¸˜ao pro-

posta para o passo de Metropolis-Hastings na amostragem de Θ em (3.1). Um valor candidato

para Θ, Θ∗, será amostrado de uma distribuic¸˜ao multivariada obtida combinando as idéias de

West et al. (1985) e Fr ühwirth-Schnater (1994). Especificamente, a distribuic¸˜ao normal multivariada N(Θ∗|ms, Cs), com média e variância on-line - m, C - aproximadas pelo Conjugate Updating, ao

invés do Filtro de Kalman. Cada θ∗t é amostrado seq üencialmente de t = T até t = 1, de suas

distribuic¸ ˜oes retrospectivas, dadas pela fatorizac¸˜ao de N(Θ∗|ms, Cs) em T densidades condicionais

univariadas, de maneira análoga ao caso gaussiano.

A principal idéia de nosso esquema é descrita da seguinte maneira. Sejam m(i) e C(i), t =

1, . . . , T, os dois primeiros momentos das distribuic¸˜oes on-line obtidas via o Conjugate Updating, em

(3.3), e sejam ms(i) e Cs(i) os dois primeiros momentos das distribuic¸ ˜oes retrospectivas na i−ésima t

iterac¸˜ao; ent˜ao, uma amostra de Θ pode ser obtida nas seguintes três etapas:

(1). Amostre Θ∗ de N(m(i), C(i));

(2). Amostre Θ∗

t de N(ms(i) Cs(i)),

t = T − 1, . . . , 1; onde

ms(i) = m(i) + C(i)G

G + W)−1(θ∗

− G

)

t(GtC(i)

t+1

tm(i)

Cs(i) = C(i) − C(i)G

G + W)−1G

t(GtC(i)

tC(i)

(3). Fac¸a Θ(i) = Θ∗ com probabilidade π e Θ(i) = Θ(i−1) com probabilidade 1 − π, onde π =

min(1, A) e

p(Θ∗|YT, Φ) q(Θ|ms(i), Cs(i))

A =

(3.5)

p(Θ|YT, Φ) q(Θ∗|ms(i), Cs(i))

onde q(·) é a densidade proposta com momentos m e C.

Em outras palavras, com o CUBS exploramos o uso de um algoritmo seq üencial, o Conjugate

Updating, em forma conjunta com um método de amostragem, o FFBS. Basicamente, com o objetivo

de construir nossa densidade proposta, introduzimos, a cada instante de tempo t, em (3.1), uma

priori conjugada (CP) para o parâmetro natural da fam´ılia exponencial. Da´ı que o CUBS é parecido

com um algoritmo de dados aumentados (data augmentation) onde os parâmetros can ônicos atuam

como variáveis latentes que ajudam na simulac¸˜ao eficiente dos estados. Contudo, dado que várias

aproximac¸ ˜oes s˜ao realizadas para tentar o uso de uma amostragem em bloco, este resultado é

uma alternativa para a construc¸˜ao de uma densidade proposta independente para um passo do

Metropolis-Hastings. Em particular, neste cap´ıtulo enfatizamos o uso de uma distribuic¸˜ao normal

multivariada como densidade proposta, porém, outras distribuic¸ ˜oes, como a t−Student, podem

ser utilizadas.

Vale a pena ressaltar que a distribuic¸˜ao proposta obtida via CUBS fornece uma boa aproximac¸˜ao

para a distribuic¸˜ao objetivo p(Θ|·), portanto, as taxas de aceitac¸˜ao s˜ao sempre bastante razoáveis.

Além disso, o CUBS é fácil de implementar, caracter´ıstica que ajuda a diminuir o n úmero de erros

de programac¸˜ao involuntários. Finalmente, embora em alguns casos sejam necessários métodos

numéricos para calcular os valores exatos de rt and st, na prática podemos usar aproximac¸ ˜oes de

func¸ ˜oes que contribuam para diminuir o tempo computacional.

Dado o vetor de parâmetros, Θ, amostras da distribuic¸˜ao a posteriori dos outros parâmetros do

modelo podem ser obtidas seguindo as técnicas padr˜ao como Metropolis-Hastings, Slice Sampling

(Neal, 2003), ARMS, abreviac¸˜ao do inglês Adaptive Rejection Metropolis Sampling (Gilks & Wild,

1992; Gilks et al., 1995), etc. Estes métodos n˜ao afetam a convergência das cadeias geradas com

nosso esquema de amostragem, mas podem afetar o tempo computacional necessário para a con-

vergência de todas as cadeias (do algoritmo completo). Entretanto, pela nossa experiência, esse

tempo será menor que o utilizado pelo algoritmo de Gamerman (1998). Para fins de exposic¸˜ao

resumimos nosso esquema MCMC no quadro 2.

Quadro 2: MCMC + CUBS

(1). Inicializaç˜ao: dar valores iniciais θ(0), ψ(0) e iniciar as iterac¸˜oes, i = 1;

(2). Amostrar Θ(i) usando CUBS:

(a) Calcular os momentos de p(θ |tDt, ψ(i−1)), m(i) e C(i), com o Conjugate Updating;

(b) Amostrar θ∗ com o Backward Sampling.

i. Amostrar θ∗ de Normal(m(i), C(i))

ii. Amostrar θ∗,

|θ∗ , ψ(i−1)

t t = T − 1, . . . , 1, de p(θt

)

t+1

min(1, A) e A é a raz˜ao de aceitac¸˜ao do Metropolis-Hastings:

ω(θ∗)

p(θ∗)

A = min 1,

, ω

(θ∗) =

;

(θ)

q(θ∗)

(3). Amostrar ψ(i) usando, em geral, um passo de Metropolis-Hastings;

(4). Amostrar φ(i) usando, em geral, um passo de Metropolis-Hastings;

(5). Atualizaç˜ao: fac¸a i = i + 1 e volte para (2)c até a convergência.

3.3

Comparaç˜ao de Esquemas de Amostragem

3.3.1

Estudo de Monte Carlo

Nesta sec¸˜ao descrevemos um estudo de Monte Carlo que teve como objetivo comparar o uso do

CUBS dentro de um algoritmo MCMC, para fazer inferência sobre todos os parâmetros de estado

de um particular modelo linear dinâmico, versus alguns dos esquemas de amostragem propostos

por Gamerman (1998) e Geweke & Tanizaki (2001).

− Geraç˜ao dos Dados Artificiais

Geramos dados artificiais segundo um modelo relativamente simples porque nosso interesse prin-

cipal é avaliar a eficiência na gerac¸˜ao de amostras somente dos parâmetros de estado. O modelo

dinâmico Poisson de primeira ordem mostrado em (3.6), com θ0 = 0, 50 e W = 0, 01, foi utilizado

para gerar três conjuntos de dados, cada um com 100 séries temporais. A diferenc¸a entres os três

conjuntos foi o tamanho das séries temporais. Desta forma, assumimos que

yt|ηt ∼ Poisson(ηt),

t = 1, . . . , T

(3.6a)

log(ηt) = θt

(3.6b)

θt = θt−1 + wt

wt ∼ N(0, W)

(3.6c)

θ0|Y0 ∼ N(m0, C0).

(3.6d)

Em resumo, geramos valores aleat órios de θt com base em (3.6c), a seguir geramos YT com base

em (3.6a) para T = 50, 100, 300. Depois disso, rodamos nossas rotinas MCMC para obter amostras

da distribuic¸˜ao a posteriori de Θ e W, e comparamos os resultados com os valores verdadeiros. As

distribuic¸ ˜oes a priori de W e θ0 foram, respectivamente, uma gama invertida com ambos parâme-

tros iguais a 0, 001, isto é, IG(0, 001; 0, 001), e uma normal com média zero e variância grande, isto

é, N(0; 103). Em outras palavras, utilizamos prioris n˜ao informativas.

− Esquemas de Amostragem Estudados

As quantidades desconhecidas em (3.6) s˜ao (θ0, θ1 . . . , θT, W). Como a distribuic¸˜ao condicional

completa de W é gama invertida, amostramos diretamente dela, entretanto para o vetor de estados

usamos passos de Metropolis-Hastings.

Como mencionamos na sec¸˜ao 3.1, Gamerman (1998) sugere a construc¸˜ao de um modelo linear dinâmico gaussiano ajustado (ou auxiliar) tal que

yt = F θ

t t + vt,

vt ∼ N(0, ˜

Vt)

(3.7a)

θt = Gtθt−1 + wt, wt ∼ N(0, Wt),

(3.7b)

onde as observac¸ ˜oes ajustadas ou pseudo-observac¸ ˜oes ˜

yt e suas variâncias associadas ˜

Vt s˜ao

yt = ηt + (yt − µt)g (µt)

Vt = ä(ηt)[g (µt)]2.

Com base nesse “novo” modelo, os seguintes esquemas de amostragem foram propostos:

• Proposta I. Amostrar w∗, um valor candidato para os erros do sistema, w

t, da distribuic

¸ ˜ao

condicional completa obtida a partir da seguinte re-parametrizac¸˜ao de (3.7):

˜yt = Ft

Gt−jwj + vt,

Gt = G ∀t,

t = 2, . . . , T,

j=1

vt ∼ N(0, ˜

Vt),

wt ∼ N(0, W),

w1 ∼ N(m1, C1).

• Proposta II. Amostrar θ∗t da condicional completa de θt em (3.7) (densidade normal).

Na realidade há uma terceira proposta que consiste em usar o algoritmo FFBS em (3.7). Con-

tudo, esta proposta n˜ao é recomendada por Gamerman pois conduz a taxas de aceitac¸˜ao muito

baixas, por isso n˜ao é implementada neste cap´ıtulo.

Geweke & Tanizaki (2001) prop ˜oem diferentes distribuic¸ ˜oes propostas e recomendam a amos-

tragem individual dos estados, isto é, amostrar os estados um a um e n˜ao em bloco. Três dos

esquemas propostos s˜ao:

• Proposta I. Amostre θ∗t da densidade obtida a partir da equac¸˜ao do sistema. Em outras

palavras, amostrar θ∗t de N(Gt(ψ )θ

t−1, Wt).

• Proposta II. Amostrar θ∗t da densidade normal com média e variância iguais às estimativas

retrospectivas, no instante t, obtidas da aplicac¸˜ao do Filtro de Kalman Estendido (EKS). Esse

filtro consiste na aplicac¸˜ao das recurs ˜oes do Filtro de Kalman e recurs ˜oes de suavizac¸˜ao a

um modelo dinâmico normal ajustado, como em (3.7). Em outras palavras, também utiliza

pseudo-observac¸ ˜oes como em Gamerman (1998).

• Proposta III. Amostrar θ∗t de uma densidade normal com base num passeio aleatório, isto é,

a média é o valor corrente da cadeia, e variância obtida como na Proposta II.

Neste cap´ıtulo utilizamos as duas alternativas de amostragem individual proposta por Gamer-

man (1998), e as três alternativas descritas acima de Geweke & Tanizaki (2001). Adicionalmente, consideramos os momentos retrospectivos obtidos depois de aplicar o Conjugate Updating, para

construir densidades propostas individuais para θt, em forma análoga à proposta II de Geweke &

Tanizaki (2001), com este esquema avaliamos a amostragem de movimentos simples que decorre

naturalmente de nossa proposta. Em resumo, os sete esquemas que comparamos s˜ao:

I. CUBS: amostragem em bloco como apresentada na sec¸˜ao 3.2.

II. Conjugate Updating: amostragem individual. A proposta é uma densidade normal com média

e variância iguais aos momentos obtidos com a análise retrospectiva no instante t, ap ós o uso

do Conjugate Updating em toda a série.

III. De Gamerman (1998) - Proposta I: amostragem individual dos erros do sistema;

IV. De Gamerman (1998) - Proposta II: amostragem individual direta dos estados;

V. De Geweke & Tanizaki (2001) - Proposta I: a proposta é obtida a partir da equac¸˜ao do sistema; VI. De Geweke & Tanizaki (2001) - Proposta II: a proposta é uma densidade normal com média e variância iguais aos momentos retrospectivos obtidos com o Filtro de Kalman Estendido.

VII. De Geweke & Tanizaki (2001) - Proposta III: a proposta é uma densidade normal com média igual ao valor corrente da cadeia e variância igual ao segundo momento retrospectivo obtido

com o Filtro de Kalman Estendido.

− Critérios de Comparaç˜ao

Com a finalidade de comparar os esquemas de amostragem mencionados acima, calculamos o

fator de ineficiência e a raiz quadrada do erro quadrático médio para Θ e YT. Seja INEF o fator de

ineficiência, dado por

n−1 n − j

INEF = 1 + 2

j=1

onde ρj é a autocorrelac¸˜ao com defasagem j para os valores de uma cadeia de n valores (n → ∞).

INEF é uma medida de qu˜ao longe está a cadeia para representar uma amostra independente (veja

Gamerman & Lopes (2006) para maiores detalhes). Quanto maior é INEF, menor é a eficiência do

esquema.

Seja RMSEθ a raiz quadrada do erro quadrático médio para θ

rep,t, definida por

2 1/2

RMSEθ = 1

θ(k) − θ(k)

rep,t

k=1

onde ˆ

θ(k) é a média da amostra da distribuic¸˜ao a posteriori e θ(k) é o valor verdadeiro do t−ésimo

rep,t

estado do k−ésimo conjunto de dados artificiais (K = 100). Ent˜ao a média e desvio padr˜ao de

RMSEθ s˜ao dados, respectivamente, por

média(RMSEθ) = 1

RMSEθ

t=1

2 1/2

d.p.(RMSEθ) = 1

RMSEθ − média(RMSEθ)

t=1

Em outras palavras, média(RMSEθ) e d.p.(RMSEθ) s˜ao as medidas resumo do erro quadrático

médio obtidos considerando as 100 replicac¸ ˜oes de cada conjunto de dados. Os valores de Yrep,t e

RMSEY s˜ao obtidos de forma análoga.

Os outros critérios de comparac¸˜ao apresentados s˜ao: tempo computacional e taxas de aceitac¸˜ao.

Os dados foram gerados no pacote R vers˜ao 2.4 (R Development Core Team, 2005) e as rotinas

MCMC foram escritas em Ox vers˜ao 3.2 (Doornik, 2002).

− Resultados

Os resultados apresentados nesta sec¸˜ao provêm da análise de amostras de tamanho 10000 das

cadeias que foram iteradas 50000 vezes, tendo as primeiras 40000 sido descartadas como per´ıodo

de aquecimento. A única excec¸˜ao s˜ao os resultados do esquema III cujas cadeias foram iteradas somente por 15000 vezes, descartando-se as primeiras 10000 iterac¸ ˜oes. Adotamos esta estratégia

com base nos resultados de Gamerman (1998) sobre a velocidade de convergência e a alta demanda

de tempo computacional desse algoritmo. Além disso, n˜ao fazemos uso do espac¸amento entre

observac¸ ˜oes (thinning interval) para obter as estat´ısticas a posteriori, pois um dos nossos objetivos

é medir a autocorrelac¸˜ao que cada método produz em iterac¸ ˜oes consecutivas.

Em forma global os resultados indicam que sob o esquema I, as ineficiências exibidas pelas

cadeias do vetor de parâmetros s˜ao independentes da posic¸˜ao no tempo, isto é, as ineficiências das

cadeias de θt para t no in´ıcio, meio e final da série de tempo s˜ao praticamente as mesmas. Outro

resultado interessante é que n˜ao dependem do n úmero de observac¸ ˜oes da série. Embora seja uma

amostragem em bloco, as taxas de aceitac¸˜ao s˜ao sempre razoáveis, estando, em média, entre 26%

e 46%. Os esquemas de amostragem individual (II, V, VII e VII) foram os menos eficientes. Adicionalmente, as func¸ ˜oes de autocorrelac¸˜ao emp´ıricas (que n˜ao s˜ao apresentadas aqui) do esquema

I foram t˜ao boas quanto as do esquema III.

A tabela 3.1 mostra o RMSEYT, RMSEΘ, a média e desvio padr˜ao das taxas de aceitac¸˜ao, consi-

derando todos os θt, assim como a média e desvio padr˜ao dos tempos computacionais utilizados

com cada algoritmo para completar todas as iterac¸ ˜oes. Nessa tabela observamos que os valores

médios dos RMSEs, (média(RMSE)), de todos os esquemas s˜ao bastante similares, como esperado,

tanto para a resposta (YT) como para os estados (Θ). As taxas de aceitac¸˜ao dos esquemas III e

IV s˜ao as maiores, mais de 90% em todos os casos, entretanto as taxas de aceitac¸˜ao de V, VI e

VII est˜ao entre 30% e 50%. As taxas de aceitac¸˜ao do esquema I s˜ao, em média, parecidas com as taxas dos esquemas individuais. Vale a pena mencionar que para alguns conjuntos de dados, a

taxa de aceitac¸˜ao foi a mais baixa entre todos os esquemas, o que é esperado por se tratar de um

Tabela 3.1: Erro Quadrático Médio (RMSE), Taxa de Aceitaç˜ao (média de T × K taxas) e Tempo Com-

putacional (média em segundo das K = 100 amostras) para as séries temporais de tamanhos

T = 50, 100, 300.

RMSEYT

RMSEΘ

Taxa de aceitac¸˜ao

Tempo(seg)

Esquema

média

d.p.

média

d.p.

média

d.p.

T = 50

1,252

0,109

0,236

0,036

42,629

–

280,30

1,276

0,116

0,234

0,026

33,525

7,157

214,98

III

1,259

0,110

0,249

0,051

97,252

1,228

629,08

1,244

0,114

0,258

0,029

98,144

0,681

120,79

1,315

0,132

0,237

0,034

51,423

4,531

89,04

1,240

0,113

0,261

0,028

37,796

6,301

248,25

VII

1,236

0,110

0,259

0,030

44,613

5,816

219,65

T = 100

1,309

0,133

0,224

0,027

38,454

–

334,25

1,327

0,143

0,226

0,024

32,350

6,986

257,16

III

1,313

0,134

0,233

0,040

97,942

0,851

1533,09

1,434

0,382

0,254

0,049

98,712

3,228

147,62

1,383

0,152

0,286

0,038

51,461

4,041

104,54

1,312

0,142

0,235

0,019

34,271

6,864

301,94

VII

1,307

0,140

0,230

0,021

42,016

6,042

274,08

T = 300

1,608

0,303

0,228

0,035

31,352

–

247,76

1,615

0,306

0,379

0,221

32,177

9,470

314,09

III

1,610

0,298

0,227

0,037

98,428

0,687

2158,93

4,522

3,753

0,545

0,252

95,515

14,121

163,87

1,711

0,329

0,471

0,133

51,365

4,108

74,44

1,607

0,301

0,236

0,034

32,881

8,985

371,31

VII

1,606

0,302

0,273

0,087

40,361

7,802

338,89

esquema em bloco. Contudo, é diferente de zero, resultado contrário ao obtido pela amostragem

em bloco comentada em Gamerman (1998). Além disso, como esperado, observamos que a taxa

de aceitac¸˜ao diminui conforme o tamanho da série temporal aumenta. Os tempos computacionais

para os esquemas I, II, V, VI eVII foram muito parecidos (entre 200 e 300 segundos em média), porém, o esquema III utilizou o mesmo tempo para completar apenas 10% do n úmero de iterac¸ ˜oes completadas pelos outros métodos.

Considerando as 100 séries temporais, a figura 3.1 mostra um resumo, através de um diagrama

de caixa (box plot), dos logaritmos das ineficiências das cadeias amostradas dos estados θt, nos

instantes t = 5, 25, 45 para T = 50, t = 5, 45, 95 para T = 100, t = 5, 145, 295 para T = 300 e W

para todos os casos. Nesta figura observamos que para as séries de tamanho T = 50 e T = 100,

os esquemas V, VI e VII exibem as maiores ineficiências. Para o esquema III notamos que as ineficiências diminuem quando o valor de t aumenta, ou seja, as cadeias s˜ao geradas de forma

mais eficiente conforme a informac¸˜ao aumenta. Este resultado é esperado uma vez que a equac¸˜ao

de estados é obtida em func¸˜ao dos termos de erro; logo, com o passar do tempo mais observac¸ ˜oes

s˜ao utilizadas para a construc¸˜ao da distribuic¸˜ao proposta. Os diagramas de caixa do esquema I

(CUBS) mostram que com este método obtivemos as menores ineficiências em t = 5 e t = 25, e as

segundas melhores em t = 45, quando T = 50. Para T = 300 o esquema III resultou ser um pouco

melhor (em termos de eficiência) que o esquema I no meio da série temporal, por exemplo em

t = 145, e muito melhor no final da mesma, por exemplo em t = 295.

A figura 3.2 exibe os diagramas de caixa dos logaritmos das medianas das amostras das distri-

buic¸ ˜oes a posteriori de W. As estimativas s˜ao razoáveis para todos os esquemas usados. Contudo,

segundo estes resultados, o esquema V tende sempre a subestimar a variância do sistema. Este fato pode ser explicado pela pouca mistura das cadeias. Adicionalmente, vale a pena mencionar que

alguns resultados, n˜ao apresentados aqui, mostraram que se o valor de W usado para gerar dados

é “grande”, por exemplo W = 0, 1, os esquemas III e IV podem ter problemas para aproximar as distribuic¸ ˜oes a posteriori de interesse causando taxas de aceitac¸˜ao pr óximas de zero.

3.3.2

Dados Reais: Precipitaç˜ao em Tokyo

Utilizamos todos os esquemas de amostragem descritos na sec¸˜ao anterior para ajustar um modelo

dinâmico ao n úmero de ocorrências de chuva acima de 1 mm em Tokyo em cada dia do ano entre

1983-1984. Este exemplo foi analisado previamente por Kitagawa (1987) e Gamerman (1998), entre outros. O objetivo é estimar a probabilidade de sucesso (ocorrência de chuva) para cada dia calen-dário. O modelo considerado em (3.8) tem uma resposta binomial com probabilidade πt, tal que

●

4.4

●

4.2

●

4.0

●

3.8

●

3.6

●

log(inefficiency)

3.4

●

3.2

●

3.0

III

VII

III

VII

III

VII

III

VII

(a) θ5

(b) θ25

(d) W

4.5

●

4.0

●

3.5

●

3.0

●

2.5

●

log(inefficiency)

2.0

●

1.5

●