Temos por objetivo aqui apresentar uma breve introdução a dois conceitos fundamentais que são necessários para compreender a estrutura básica dos Modelos Ocultos de Markov (HMM). Em primeiro lugar, é apresentada uma breve panorâmica das distribuições de misturas, dado que a distribuição marginal de um Modelo Oculto de Markov (HMM) é um modelo de mistura discreta. Posteriormente, são introduzidas Cadeias de Markov, que são de grande importância, uma vez que o processo de parâmetros de um HMM é modelado por uma Cadeia de Markov.



1 Mistura independente de distribuições


Reconsiderar a série de vendas semanais de sabonetes introduzida na Introdução (Figura 1) e vamos assumir que não há correlação serial na série, ou seja, as vendas semanais são contagens independentes. O modelo básico para contagens independentes não vinculadas é a distribuição de Poisson com função de probabilidade \[ P(S)=\lambda^s e^{-\lambda}/s! \] e a propriedade restritiva que a variância iguala a esperança \(\mbox{E}(S)=\mbox{Var}(S)\).

Para a série de vendas de sabão, tem-se

mean(soap); var(soap)
## [1] 5.442149
## [1] 15.4012

indicando uma forte dispersão em relação à distribuição de Poisson e, portanto, inapropriação desta última. Isto é apresentado na Figura 1 que apresenta um histograma das observações e a distribuição de Poisson teórica. Além disso, a distribuição das vendas semanais de sabonetes é bimodal, enquanto a distribuição de Poisson tem apenas uma moda.

library(ggplot2)
Resposta = c(rep("Teórica", 23), rep("Estimada", 23))
Poisson.est = hist(soap, breaks = 0:23, plot = FALSE)
Histograma = c(dpois(0:22, lambda = mean(soap)), Poisson.est$density)
Pontos = 0:22
dados = data.frame(Resposta, Histograma, Pontos)
ggplot(dados, aes(Pontos, Histograma, fill = Resposta))+
    geom_bar(stat = "identity",position = "dodge") + 
    theme(legend.position=c(.8,0.8))

Figura 1: Histograma das vendas semanais de um produto de sabão específico e a distribuição Poisson estimada. 

ggplot(data = dados1, aes(x = semanas, y = soap)) + 
  geom_line(color = "#00AFBB", size = 1) + geom_hline(yintercept = mean(soap)) +
  labs(x = "Semanas", y = "Unidades vendidas") + 
  annotate("text", x = semanas[55], y = mean(soap) + 1.5, 
           label = "paste(hat(lambda) == 5.44)", parse = TRUE, size = 7)

Figura 2: Série de vendas semanais de um produto de sabão específico.

Um método de superar as deficiências da distribuição Poisson em tal caso é usar um modelo de mistura independente. As distribuições de misturas podem ser aplicadas em muitos campos diferentes. Eles são muito úteis no tratamento de observações ou dados superdispersos com distribuições bimodais ou multimodais mais gerais. Os modelos de mistura assumem que a superdispersão e/ou multimodalidade das observações de uma variável podem ser devidas à heterogeneidade não observada na população, ou seja, a população consiste em grupos diferentes com distribuições diferentes para essa variável. Imagine, por exemplo, a distribuiçõo do número de pacotes de cigarros comprados por um número de clientes em um supermercado. Os clientes podem ser separados em vários grupos: não-fumantes, pessoas que fumam ocasionalmente, etc. Isso leva à superdispersão em relação à Poisson e talvez até a uma distribuição multimodal de suas compras.

Suponha que cada contagem da série de vendas de sabão seja gerada por uma das duas distribuições de Poisson com esperanças \(\lambda_1\) e \(\lambda_2\), em que a escolha da esperança é feita por algum outro mecanismo aleatório que chamamos de processo de parâmetros. Aqui, não temos heterogeneidade em um grupo de clientes, mas heterogeneidade ao longo do tempo. Existem fases com baixo índice de vendas e fases com alta taxa de vendas. Isso pode ser devido a covariáveis, como o preço que não consideramos aqui. Por exemplo, se \(\lambda_1\) é selecionado com probabilidade \(\delta_1\) e \(\lambda_2\) com probabilidade \(\delta_2=1-\delta_1\), então a variância da série resultante é maior que sua esperança, como será mostrado mais adiante. Se o processo de parâmetros é uma série de variáveis aleatórias independentes, então as contagens também são independentes. Este é um exemplo de uma mistura independente.

Em geral, uma distribuição de mistura independente consiste em um certo número de distribuições componentes. Este é o caso de uma mistura discreta. Em muitas aplicações, não é razoável supor que a heterogeneidade da população seja representada por um número finito de componentes com distribuiçãção discreta (Böhning, 2000). Em vez disso, pode-se pensar em um modelo de mistura contínua que podem ser interpretados como um caso especial de uma mistura discreta com inúmeros componentes. Um exemplo proeminente para uma mistura contínua é a distribuição binomial negativa que pode ser derivada como uma mistura de distribuições Poisson, onde o parâmetro \(\lambda\) têm distribuição Gama. Como as misturas contínuas não são relevantes para os HMMs, a seguir, consideramos apenas distribuições discretas de misturas. Mais detalhes sobre misturas contínuas podem pode ser encontrados em Zucchini, Böker, and Stadie (2001) ou Böhning (2000).

No caso de duas distribuições componentes, a distribuição de mistura é caracterizada funções de probabilidade ou funções de densidade associadas às variáveis \(X_1\) e \(X_2\), respectivamente: \[ \begin{array}{lcc} \mbox{Variável aleatória} & X_1 & X_2 \\ \hline \mbox{função de probabillidade} & p_1 & p_2 \\ \mbox{função de densidade} & f_1 & f_2 \\\hline \end{array} \]

Além disso, para o processo de parâmetros, é necessária uma variável aleatória discreta \(C\) para realizar a mistura: \[\begin{equation} C \, = \, \left\{ \begin{array}{cl} 1, & \mbox{com probabilidade} \; \delta_1 \\ 2, & \mbox{com probabilidade} \; \delta_2 \, = \, 1 \,- \,\delta_1 \end{array}\right.\cdot \end{equation}\]

Pode-se imaginar \(C\) como o resultado de jogar uma moeda: Se \(C\) assume o valor 1, obtenha uma observação de \(X_1\). De acordo com isso, se \(C\) assumir o valor 2, obtenha uma observação de \(X_2\). A estrutura desse processo para o caso de duas distribuições de componentes contínuos é mostrada na Figura I.3.

Figura 3: Estrutura do processo de uma distribuição de mistura de dois componentes.
Figura retirada do livro Zucchini & MacDonald (2009).

Em situações práticas, não sabemos para que lado a moeda pousou. Apenas as observações geradas por \(X_1\) ou \(X_2\) podem ser observadas e, na maioria dos casos, elas não podem ser atribuídas a variáveis aleatória distintas. Dada a probabilidade, cada componente e as respectivas distribuições de probabilidade, a função de probabilidade ou de densidade da mistura pode ser facilmente calculada.

Seja \(X\) a realização da mistura. Então, a função de probabilidade ou de densidade é dada por \[\begin{equation} \begin{array}{ccccl} p(x) & = & \delta_1 p_1(x) \, + \, \delta_2 p_2(x) & & \mbox{(caso discreto)}, \\ f(x) & = & \delta_1 f_1(x) \, + \, \delta_2 f_2(x) & & \mbox{(caso contínuo)}\cdot \end{array} \end{equation}\]

A demonstração para o caso discreto é direta: \[\begin{equation} \begin{array}{ccl} P(X=x) & = & P\Big( \big\{\{X_1=x\}\cap \{C=1\}\big\}\cup \big\{\{X_2=x\}\cap \{C=2\} \big\} \Big) \\ & = & P(X_1=x,C=1) \, + P(X_2=x,C=2) \\ & = & P(X_1=x|C=1)P(C=1) \, + \, P(X_2=x|C=2)P(C=2) \\ & = & p_1(x)\delta_1 \, + \, p_2(x)\delta_2\cdot \end{array} \end{equation}\]

Deste resultado para o caso contínuo pode ser comprovado de forma análoga. Dois exemplos gráficos de misturas de dois componentes são apresentados na Figura I.4.

Figura 4: Distribuição da mistura de dois componentes com componentes discretos ou contínuos.
Figura retirada do livro Zucchini & MacDonald (2009).

A extensão para o caso de \(m\) componentes é simples. Sejam \(\delta_1,\delta_2,\cdots,\delta_m\) os pesos atribuídos aos diferentes componentes, \(p_1,p_2,\cdots,p_m\) ou \(f_1, f_2, \cdots, f_m\) denotando suas funções de probabilidade ou de densidade, respectivamente. Então, a distribuição da mistura dada pela variável aleatória \(X\) pode ser facilmente calculada por uma combinação linear das distribuições componentes: \[\begin{equation} \begin{array}{ccccl} p(x) & = & \displaystyle \sum_{i=1}^m \delta_i p_i(x) & & \mbox{(caso discreto)}, \\ f(x) & = & \displaystyle \sum_{i=1}^m \delta_i f_i(x) & & \mbox{(caso contínuo)}\cdot \end{array} \end{equation}\]

Além disso, para uma mistura de \(m\) distribuições discretas, obtém-se \[\begin{equation} \mbox{E}(X) \, = \, \displaystyle \sum_x xp(x) \, = \, \sum_x x \sum_{i=1}^m \delta_i p_i(x) \, = \, \sum_{i=1}^m \delta_i \sum_x x p_i(x) \,= \, \sum_{i=1}^m \delta_i \mbox{E}(X_i)\cdot \end{equation}\]

Um resultado análogo pode ser derivado para uma mistura de distribuições contínuas. Assim, no caso de dois componentes discretos ou contínuos, pode-se afirmar que: \[\begin{equation} \mbox{E}(X) \, = \, \delta_1\mbox{E}(X_1) \, + \, \delta_2\mbox{E}(X_2) \, = \, \delta_1\mbox{E}(X_1) \, + \, (1-\delta_1)\mbox{E}(X_2)\cdot \end{equation}\]

Mais geral, o \(k\)-ésimo momento \(\mbox{E}(X^k)\) de uma mistura é simplesmente uma combinação linear dos respectivos momentos de seus componentes: \[\begin{equation} \mbox{E}(X^k) \, = \, \sum_{i=1}^m \delta_i\mbox{E}(X_i^k), \qquad k=1,2,\cdots\cdot \end{equation}\] Agora, observe que isto não é válido para os momentos centrais, por exemplo, a variância de uma mistura: \[\begin{equation} \mbox{Var}(X) \, \neq \, \sum_{i=1}^m \delta_i\mbox{Var}(X_i)\cdot \end{equation}\] Por exemplo, no caso de dois componentes, a variância da mistura pode ser calculada a partir da seguinte expressão: \[\begin{equation} \mbox{Var}(X) \, = \, \delta_1\mbox{Var}(X_1) \, + \, \delta_2\mbox{Var}(X_2) \, + \, \delta_1\delta_2\big(\mbox{E}(X_1)-\mbox{E}(X_2)\big)^2\cdot \end{equation}\]

Em geral, a variância de um modelo de mistura pode ser calculada usando a expressão: \[\begin{equation} \mbox{Var}(X) \, = \, \mbox{E}(X^2) - \big(\mbox{E}(X)\big)^2, \end{equation}\] em combinação com o resultado acima para \(\mbox{E}(X^k)\).


Exemplo 1.


Seja a distribuição da mistura de \(X\) Poisson. Então \[\begin{equation} \mbox{E}(X) \, = \, \displaystyle \sum_{i=1}^m \delta_i\lambda_i \qquad \mbox{e} \qquad \mbox{E}(X^2)=\sum_{i=1}^m \delta_i (\lambda_i+\lambda_ i^2), \end{equation}\] e segue-se que \[\begin{equation} \mbox{Var}(X) \, = \, \sum_{i=1}^m \delta_i(\lambda_i+\lambda_ i^2) \, - \, \Big( \sum_{i=1}^m \delta_i\lambda_i\Big)^2\cdot \end{equation}\]


A estimação dos parâmetros de uma distribuição de mistura é geralmente realizada utilizando o método
de máxima verossimilhança (ML). Em geral, a verossimilhança de um modelo de mistura com \(m\) componentes é dada por \[\begin{equation} L(\theta,\delta;x) \, = \, \prod_{j=1}^n \sum_{i=1}^m \delta_ip_i(x_j;\theta_i), \qquad \mbox{caso discreto} \end{equation}\] ou \[\begin{equation} L(\theta;,\delta;x) \, = \, \prod_{j=1}^n \sum_{i=1}^m \delta_if_i(x_j;\theta_i), \qquad \mbox{caso contínuo} \end{equation}\] onde \(\theta=(\theta_1,\cdots,\theta_m)\) são os vetores de parâmetros das distribuições componentes, \(\delta=(\delta_1,\cdots,\delta_m)\) são os parâmetros de mistura e \(x=(x_1,\cdots,x_n)\) são as observações. Assim, no caso de distribuições componentes com um parâmetro, \(2m-1\) parâmetros devem ser estimados, ou seja, \(m\) parâmetros para as distribuições componentes e \(m-1\) parâmetros de mistura, isso porque os parâmetros de mistura devem somar para um. A maximização da verossimilhança nâo é trivial, uma vez que não é possível resolver o problema de maximização analiticamente. Isso vai ser demonstrado no seguinte exemplo.


Exemplo 2.


Suponha que \(X_1\sim Poisson(\lambda_1)\) e \(X_2\sim Poisson(\lambda_2)\) e que \(\delta_1\) e \(\delta_2\) sejam os parâmetros de mistura com \(\delta_1+\delta_2 = 1\). Isto significa que a distribuição da mistura \(p(x)\) é dada por \[\begin{equation} p(x)=\delta_1 p_1(x)+\delta_2 p_2(x), \end{equation}\] onde \(p_1(x)=\lambda_1e^{-\lambda_1}/x!\) e \(p_2(x)=\lambda_2e^{-\lambda_2}/x!\).

Como \(\delta_2 = 1 -\delta_1\), três parâmetros devem ser estimados: \(\lambda_1\)¸ \(\lambda_2\) e \(\delta_1\). Vamos denotar por \(x=(x_1,x_2,\cdots,x_n)\) os valores observados da distribuição da mistura com a função de distribuição \(p(x)\).

Então, a função de verossimilhança é dada pela seguinte expressão: \[\begin{equation} \begin{array}{ccl} L(\lambda_1,\lambda_2,\delta_1;x) & = & p(x_1)\times p(x_2)\times \cdots \times p(x_n) \, = \, \displaystyle \prod_{i=1}^n p(x_i) \\ & = & \displaystyle \prod_{i=1}^n \left( \delta_1\dfrac{\lambda_1^{x_i}e^{-\lambda_1}}{x!}+(1-\delta_1)\dfrac{\lambda_2^{x_i}e^{-\lambda_2}}{x!}\right)\cdot \end{array} \end{equation}\]

A maximização em relação a \(\lambda_1\), \(\lambda_2\) e \(\delta_1\) não é trivial, pois esta função é o \(n\)-ésimo produto de uma soma. Tomando o logaritmo e diferenciando não simplifica a expressão.

Para a estimação dos parâmetros podemos utilizar o pacote de funções FlexMix. Nele estão implementadas funções com uma estrutura geral para misturas finitas de modelos de regressão. A estimação dos parâmetros é realizada usando o algoritmo EM: o passo-E é implementado pela função flexmix, enquanto o usuário pode especificar o passo-M.

library(flexmix)
# Mistura de duas distribuições Poisson
ajuste1 = flexmix(soap ~ 1, k = 2, model = FLXglm(family = "poisson"))
summary(ajuste1)
## 
## Call:
## flexmix(formula = soap ~ 1, k = 2, model = FLXglm(family = "poisson"))
## 
##        prior size post>0 ratio
## Comp.1 0.148   26    233 0.112
## Comp.2 0.852  216    238 0.908
## 
## 'log Lik.' -632.971 (df=3)
## AIC: 1271.942   BIC: 1282.409
parameters(ajuste1, which = "concomitant")
##               1         2
## prior 0.1477053 0.8522947
parameters(ajuste1, which = "model")
## Comp.1.coef.(Intercept) Comp.2.coef.(Intercept) 
##                2.525426                1.439177
exp(parameters(ajuste1, which = "model"))
## Comp.1.coef.(Intercept) Comp.2.coef.(Intercept) 
##               12.496216                4.217223
EIC(ajuste1)
## [1] 0.8213469
AIC(ajuste1)
## [1] 1271.942
BIC(ajuste1)
## [1] 1282.409
Teorico.mix = parameters(ajuste1, which = "concomitant")[1]*
   dpois(0:22, lambda = exp(parameters(ajuste1, which = "model")[1])) +
   parameters(ajuste1, which = "concomitant")[2]*
   dpois(0:22, lambda = exp(parameters(ajuste1, which = "model")[2]))
Histograma1 = c(Teorico.mix, Poisson.est$density)
dados.hist1 = data.frame(Resposta, Histograma1, Pontos)
ggplot(dados.hist1, aes(Pontos, Histograma1, fill = Resposta)) + 
  labs(y="Frequência", x = " ") + 
  ggtitle("contagens e mistura de dois componentes")+
  geom_bar(stat = "identity", position = "dodge") + 
  theme(legend.position=c(.8,0.8))

# Mistura de três distribuições Poisson
ajuste2 = flexmix(soap ~ 1, k = 3, model = FLXglm(family = "poisson"))
summary(ajuste2)
## 
## Call:
## flexmix(formula = soap ~ 1, k = 3, model = FLXglm(family = "poisson"))
## 
##        prior size post>0 ratio
## Comp.1 0.112   23    219 0.105
## Comp.2 0.182   23    219 0.105
## Comp.3 0.706  196    239 0.820
## 
## 'log Lik.' -627.4394 (df=5)
## AIC: 1264.879   BIC: 1282.323
parameters(ajuste2, which = "concomitant")
##               1         2         3
## prior 0.1119638 0.1824616 0.7055746
parameters(ajuste2, which = "model")
## Comp.1.coef.(Intercept) Comp.2.coef.(Intercept) Comp.3.coef.(Intercept) 
##               2.6160491               0.7046014               1.6138970
exp(parameters(ajuste2, which = "model"))
## Comp.1.coef.(Intercept) Comp.2.coef.(Intercept) Comp.3.coef.(Intercept) 
##               13.681563                2.023040                5.022345
EIC(ajuste2)
## [1] 0.6456009
AIC(ajuste2)
## [1] 1264.879
BIC(ajuste2)
## [1] 1282.323
Teorico.mix = parameters(ajuste2, which = "concomitant")[1]*
   dpois(0:22, lambda = exp(parameters(ajuste2, which = "model")[1])) +
   parameters(ajuste2, which = "concomitant")[2]*
   dpois(0:22, lambda = exp(parameters(ajuste2, which = "model")[2])) + 
   parameters(ajuste2, which = "concomitant")[3]*
   dpois(0:22, lambda = exp(parameters(ajuste2, which = "model")[3]))
Histograma2 = c(Teorico.mix, Poisson.est$density)
dados.hist2 = data.frame(Resposta, Histograma2, Pontos)
ggplot(dados.hist2, aes(Pontos, Histograma2, fill = Resposta)) + 
   labs(y="Frequência", x = " ") + 
   ggtitle("contagens e mistura de três componentes")+
   geom_bar(stat = "identity", position = "dodge") + 
   theme(legend.position=c(.8,0.8))

Desta forma ajustamos uma distribuição de mistura com duas e três funções de probabilidade Poisson independentes, sendo escolhido como modelo mais adequado àquele com três misturas Poisson independentes.

As estimativas dos parâmetros são \(\widehat{\delta}_1=0.1867675\), \(\widehat{\delta}_2=0.1117240\), \(\widehat{\delta}_3=0.701585\), \(\widehat{\lambda}_1=2.050877\), \(\widehat{\lambda}_2=13.689673\) e \(\widehat{\lambda}_3=5.035153\).

Ainda temos que o Critério de Informação de Entropía (EIC) é 0.6408171, um valor equivalente à medida \(R^2\) no modelo de regressão linear, um valor de \(AIC\) de 1264.896 e \(BIC\) de 1282.341.


Tanto os valores do AIC quanto do BIC podem ser ligeiramente melhorados usando um modelo de mistura de quatro componentes mas a melhoria é demasiado pequena para justificar os parâmetros adicionais.


2 Cadeias de Markov


A teoria das Cadeias de Markov é bem elaborada, no entanto, aqui podemos apenas dar uma breve introdução às Cadeias de Markov e algumas de suas propriedades básicas necessárias para a construção dos HMMs. Para uma descrição detalhada ver Cadeias de Markov.

Considere um processo estocástico, isto é, uma sequência de variáveis aleatórias \(\{C_t \, : \, t\in T\}\) indexada por algum conjunto \(T\) e que aceita valores em algum conjunto \(S\), chamado de espaço de estados. O conjunto do tempo \(T\) pode ser discreto ou contínuo, mas, a seguir, nos concentramos no caso discreto, isto é, nos processos a tempo discreto \(\{C_1,C_2,C_3,\cdots\}=\{C_t \, : \, t\in T\}\) onde \(T = \{1,2,3,\cdots\}\).

Um processo estocástico \(\{C_t \, : \, t = 1,2,\cdots\}\) é um processo de Markov se para cada tempo \(t\), o próximo estado \(C_{t+1}\) depender apenas do estado atual da cadeia, \(C_t\). Ou seja, dados \(C_t, C_{t+1}\) não dependem mais da história da cadeia \(\{C_1,C_2,\cdots,C_{t-1}\}\). Mais matematicamente, pode-se definir um processo de Markov da seguinte maneira.

Seja \(c_1,c_2,\cdots,c_t,c_{t+1}\), \(t\in\{1,2,\cdots\}\) uma sequência de observações de um processo estocástico \(\{C_s, \, s=1,2,\cdots\}\). \(\{C_s\}\) é um processo de Markov se tiver a propriedade Markov \[\begin{equation} P(C_{t+1}=c_{t+1} \, | \, \underbrace{C_t=c_t, C_{t-1}=c_{t-1},\cdots,C_1=c_1}_{\text{"história inteira"}}) \, = \, P(C_{t+1}=c_{t+1} \, | \, C_t=c_t), \end{equation}\] para todo \(t\in\{1,2,\cdots\}\). Ou seja, dado o presente do processo, \(C_t\), seu futuro, \(C_{t+1}\), é independente de seu passado, \(C_{t-1}, C_{t-2},\cdots, C_1\). Essa estrutura de um processo de Markov é demonstrada na Figura I.5.

Figura 5: Um processo de Markov.

Até agora, não especificamos mais o espaço de estados \(S\). Análogo ao conjunto do tempo \(T\), o espaço de estados \(S\) de um processo de Markov pode ser contínuo ou discreto. Um processo de Markov cujo espaço de estados é discreto chama-se Cadeia de Markov. A seguir, trataremos apenas de Cadeias de Markov cujo espaço de estados é um conjunto limitado de números inteiros \(S = \{1,2,\cdots,m\}\), isto é, \(c_t \in \{1,2,\cdots,m\}\) para todos os \(t \in \{1,2,\cdots\}\). Então, se \(C_t = i, \, i\in\{1,2,\cdots,m\}\) dizemos que a Cadeia de Markov está no \(i\)-ésimo estado no momento \(t\).

Para descrever uma Cadeia de Markov, as probabilidades de mudar de um estado \(i\) para outro estado \(j\), \(P(C_{t + 1} = j \, | \, C_t = 1)\), devem ser consideradas mais próximas. Em geral, existem duas possibilidades: essas chamadas “probabilidades de transição” podem mudar no tempo ou podem ser estáveis no tempo. O último caso pode ser caracterizado da seguinte maneira.

Uma Cadeia de Markov é chamada de cadeia homogênea ou de Markov com probabilidades de transição estacionárias se as probabilidades de transição forem independentes de \(t\), isto é, \[\begin{equation} P(C_{t+1}=j \, | \, C_t=i) \, = \, \gamma_{_{i,j}}, \end{equation}\] para todo \(t\in \{1,2,\cdots\}\) e \(i,j\in \{1,\cdots,m\}\). Então, \(\gamma_{i,j}\) denota a probabilidade de transição do estado \(i\) para o estado \(j\). A seguir, trataremos apenas de cadeias homogêneas de Markov.

A função de probabilidades de transição de uma Cadeia de Markov homogênea de \(m\) estados pode ser representada por uma matriz \(m\times m\), a chamada matriz de probabilidades de transição \(\Gamma\): \[\begin{equation} \Gamma \, = \, \begin{pmatrix} \gamma_{_{1,1}} & \cdots & \gamma_{_{1,m}} \\ \vdots & \ddots & \vdots \\ \gamma_{_{m,1}} & \cdots & \gamma_{_{m,m}} \end{pmatrix}, \end{equation}\] com \(\gamma_{_{i,j}}=P(C_{t+1}=j \, | \, C_t=i)\) e \(\displaystyle \sum_{j=1}^m \gamma_{_{i,j}}=1\), \(\forall i\).

Como cada linha \(i\) de \(\Gamma\) descreve a função de probabilidade, ou seja, as probabilidades de transição para mudar do estado \(i\) para o estado \(j\), a linha é igual a um. A matriz de probabilidade de transição \(\Gamma\) contêm as probabilidades de transição em uma etapa e, portanto, descreve o comportamento de curto prazo da Cadeia de Markov. Para descrever o comportamento de longo prazo de uma Cadeia de Markov, podem-se definir as probabilidades de transição em \(k\) etapas, ou seja, as probabilidades de transição do estado \(i\) em \(t\) para o estado \(j\) em \(t + k\): \[\begin{equation} \gamma_{_{i,j}}^{(k)} \, = \, P(C_{t+k}=j \, | \, C_t=i)\cdot \end{equation}\]

As probabilidades de transição em \(k\)-passos ser calculadas por meio de alguma álgebra matricial, pois sustenta que a matriz \(\Gamma^{(k)}\), que contêm as probabilidades de transição \(k\)-passos é simplesmente a \(k\)-ésima potência da matriz de probabilidade de transição \(\Gamma\): \[\begin{equation} \Gamma^{k} \, = \, \begin{pmatrix} \gamma_{_{1,1}}^{(k)} & \cdots & \gamma_{_{1,m}}^{(k)} \\ \vdots & \ddots & \vdots \\ \gamma_{_{m,1}}^{(k)} & \cdots & \gamma_{_{m,m}}^{(k)} \end{pmatrix}\cdot \end{equation}\] Este resultado segue das equações de Chapman-Kolmogorov (para uma prova, ver Grimmett and Stirzaker (2001)): \[\begin{equation} \gamma_{_{i,j}}^{(m+n)} \, = \, \displaystyle \sum_{k=1}^m \gamma_{_{i,k}}^{(m)}\gamma_{_{k,j}}^{(n)}\cdot \end{equation}\]

Nesse contexto, diz-se que o estado \(i\) se comunica com o estado \(j\), escrito \(i\to j\), se a cadeia puder visitar o estado \(j\) com probabilidade positiva, começando pelo estado \(i\). Ou seja, \(i\to j\) se existir algum \(k\in \{1,2,\cdots\}\) com \(\gamma_{_{i,j}}^{(k)}> 0\). Além disso, diz-se que \(i\) e \(j\) intercomunicam e escrevem-se \(i\longleftrightarrow j\) se \(j \to i\). Então, uma Cadeia de Markov é definida como irredutível se \(i\longleftrightarrow j\) para todos os \(i,j \in \{1,2,\cdots,m\}\). A seguir, como na maioria das aplicações, assumimos que a Cadeia de Markov é irredutível.

As probabilidades de transição em \(k\)-passos fornecem as probabilidades condicionais de estar no estado \(j\) no tempo \(t + k\), dado que esteve no estado \(i\) no tempo \(t\). No entanto, em geral, também interessa a probabilidade marginal de estar no estado \(i\) em um determinado momento \(t\), \(\delta_i(t)\). Dada a distribuição de probabilidade inicial para o primeiro estado \[\begin{equation} \delta(1) \, = \, \big(\delta_1(1),\delta_2(1),\cdots,\delta_m(1) \big) \, = \, \big( P(C_1=1), P(C_1=2), \cdots, P(C_1=m)\big), \end{equation}\] a função de probabilidade para o estado \(C_t\), \(\delta(t)\), pode ser calculada como \[\begin{equation} \delta(t) \, = \, \big(\delta_1(t),\delta_2(t),\cdots,\delta_m(t) \big) \, = \, \big( P(C_t=1), P(C_t=2), \cdots, P(C_t=m)\big)\, = \, \delta(1)\Gamma^{k-1}\cdot \end{equation}\] Uma questão então levantada é se qualquer afirmação pode ser feita com relação ao comportamento de \(\delta(t)\) para um \(t\) muito grande.

Pode-se mostrar que, caso a Cadeia de Markov seja homogênea e irredutível, \(\delta(t)\) converge para um vetor fixo, digamos \(\delta = (\delta_1,\delta_2,\cdots,\delta_m)\), para \(t\) grande, e que \(\delta\) é o vetor único de comprimento \(m\) obtido como resolução de \[\begin{equation} \sum_{i=1}^m \delta_i \gamma_{_{i,j}} \, = \, \delta_j, \qquad \forall j\in \{1,2,\cdots,m\} \qquad \mbox{sujeito a} \qquad \sum_{j=1}^m \delta_j \, = \, 1, \end{equation}\] ou equivalente \[\begin{equation} \delta = \delta \Gamma \qquad \mbox{sujeito a} \qquad \delta\pmb{1}^\top \, = \, 1\cdot \end{equation}\]

Para uma prova desse resultado, ver Seneta (1981). O sistema de equações pode ser resolvido facilmente. Por exemplo, para uma Cadeia de Markov de dois estados com matriz de probabilidade de transição dada por \[\begin{equation} \Gamma \, = \, \begin{pmatrix} \gamma_{_{1,1}} & \gamma_{_{1,2}} \\ \gamma_{_{2,1}} & \gamma_{_{2,2}} \end{pmatrix} \end{equation}\] obtém-se \[\begin{equation} \delta \, = \, (\delta_1,\delta_2) \, = \, \dfrac{1}{\gamma_{_{1,2}}+\gamma_{_{2,1}}}(\gamma_{_{2,1}}+\gamma_{_{1,2}})\cdot \end{equation}\]

O vetor \(\delta\) é chamado de distribuição estacionária de \(\Gamma\) e uma Cadeia de Markov é considerada estacionária, se a distribuição estacionária \(\delta\) existe e se vale para cada \(t\), \(t\in \{1,2,\cdots\}\), em particular para a distribuição inicial do primeiro estado, \(\delta(1) = \delta\). Na prática, depende da aplicação se assumir ou não a Cadeia de Markov subjacente de um HMM estacionária; no entanto, além do exemplo a seguir, consideramos apenas Cadeias de Markov estacionárias.


Exemplo 3.


A seguir, consideramos um exemplo simples para esclarecer as propriedades dadas das Cadeias de Markov. Imagine uma sequência de dias chuvosos e ensolarados, começando com um dia ensolarado, modelada de tal maneira que o clima de amanhã dependa apenas do clima de hoje e as probabilidades de transição da respectiva Cadeia de Markov sejam dadas pela tabela a seguir: \[ \begin{array}{l|cc} & \mbox{dia t+1} & & & \\ \mbox{dia t} & \mbox{chuvoso} & \mbox{ensolarado} \\\hline \mbox{chuvoso} & 0.9 & 0.1 \\ \mbox{ensolarado} & 0.6 & 0.4 \\\hline \end{array} \] Essas probabilidades de transição podem ser interpretadas da seguinte maneira. Supondo que hoje seja chuvoso, a probabilidade de obter um dia chuvoso amanhã é de 90%, a probabilidade de um dia ensolarado é de apenas 10%. Por analogia, a probabilidade de um dia ensolarado ser seguido por outro dia ensolarado é de 40, enquanto a probabilidade de um dia chuvoso seguir é de 60%. Isso é monstrado novamente na Figura I.6.

Figura 6: Probabilidades de transição para uma sequência de dias chuvosos e ensolarados.

Essa situção pode ser modelada por uma Cadeia de Markov simples de dois estados, a saber, sua matriz de probabilidade de transição \(\Gamma\), dada por \[\begin{array}{rcl} \Gamma & = & \begin{pmatrix} \gamma_{_{1,1}} & \gamma_{_{1,2}} \\ \gamma_{_{2,1}} & \gamma_{_{2,2}} \end{pmatrix} \\ & = & \begin{pmatrix} P(C_{t+1}=1 \, | \, C_t=1) & P(C_{t+1}=2 \,| \, C_t=1) \\ P(C_{t+1}=1 \, | \, C_t=2) & P(C_{t+1}=2 \, | \, C_t=2) \end{pmatrix} \\ & = & \begin{pmatrix} P(C_{t+1}=\text{chuvoso} \, | \, C_t=\text{chuvoso}) & P(C_{t+1}=\text{ensolarado} \,| \, C_t=\text{chuvoso}) \\ P(C_{t+1}=\text{chuvoso} \, | \, C_t=\text{ensolarado}) & P(C_{t+1}=\text{ensolarado} \, | \, C_t=\text{ensolarado}) \end{pmatrix} \\ & = & \begin{pmatrix} 0.9 & 0.1 \\ 0.6 & 0.4 \end{pmatrix}\cdot \end{array}\] Além disso, como assumimos que a Cadeia de Markov começa com um dia ensolarado, temos a seguinte distribuição de probabilidade canônica do clima de hoje \(\delta(1)\): \[\begin{array}{rcl} \delta(1) & = & \big( \delta_1(1) \quad \delta_2(1)\big) \, = \, \big( P(C_1=1) \quad P(C_1=2)\big) \\ & = & \big(P(C_1=\text{chuvoso}) \quad P(C_1=\text{ensolarado}) \big) \, = \, (0 \quad 1)\cdot \end{array}\] Dado \(\Gamma\) e \(\delta(1)\), as distribuições de probabilidade do clima dos dias seguintes podem ser calculadas da seguinte forma: \[\begin{array}{rcl} \delta(2) & = & \big( P(C_2=1) \quad P(C_2=2)\big) \, = \, \delta(1) \, \Gamma \\ & = & (0 \quad 1)\begin{pmatrix} 0.9 & 0.1 \\ 0.6 & 0.4 \end{pmatrix} \, = \, (0.6 \quad 0.4), \\ \delta(3) & = & \delta(2) \, \Gamma \\ & = & (0.6 \quad 0.4)\begin{pmatrix} 0.9 & 0.1 \\ 0.6 & 0.4 \end{pmatrix} \, = \, (0.78 \quad 0.22), \\ \vdots & \vdots & \vdots \end{array}\]

Assim, por exemplo, a probabilidade de que o dia depois de amanhã seja um dia chuvoso se hoje for um dia ensolarado é igual a 78%.

Além disso, resolvendo \[\begin{equation} \big( \delta_1 \quad \delta_2\big) \, = \, \big( \delta_1 \quad \delta_2\big) \begin{pmatrix} 0.9 & 0.1 \\ 0.6 & 0.4 \end{pmatrix} \qquad \text{restrito a} \qquad \delta_1+\delta_2=1 \end{equation}\] ou usando a fórmula dada acima leva à distribuição estacionária \(\delta=(0.86 \quad 0.14)\).

Isso significa que a probabilidade de um dia no futuro, o que é longe o suficiente para hoje, ser chuvosa &é de cerca de 86%, enquanto a probabilidade de um dia ensolarado é igual a apenas 14%, independentemente do fato de hoje ser um dia ensolarado.

Além disso, se não soubéssemos o estado inicial e, portanto, \(\delta(1)\), e se tivéssemos assumido que a Cadeia de Markov era estacionária, \(\delta\) denotaria a distribuição de probabilidade marginal do tempo para qualquer dia, até hoje, chamada de distribuição estacionária.


2.1 Distribuição estacionária


Diz-se que uma Cadeia de Markov com a matriz de probabilidades de transição \(\Gamma\) tem por distribuição estacionária \(\pi\), um vetor linha com elementos não negativos, se \[\begin{equation} \pi\Gamma \, = \, \pi \qquad \text{ e } \qquad \pi\pmb{1}^\top \, = \, 1\cdot \end{equation}\] O primeiro desses requisitos expressa a estacionariedade e o segundo requisito verifica se \(\pi\) é de fato uma função de probabilidade.

Por exemplo, a Cadeia de Markov com matriz de probabilidades de transição dada por \[\begin{equation} \Gamma \, = \, \begin{pmatrix} 1/3 & 1/3 & 1/3 \\ 2/3 & 0 & 1/3 \\ 1/2 & 1/2 & 0 \end{pmatrix} \end{equation}\] tem como distribuição estacionária \[\begin{equation} \pi \, = \, \dfrac{1}{32}\begin{pmatrix} 15 & 9 & 8 \end{pmatrix}\cdot \end{equation}\]

Um resultado geral que pode ser usado convenientemente para calcular uma distribuição estacionária é o seguinte. O vetor \(\pi\) com elementos não negativos é uma distribuição estacionária da Cadeia de Markov com matriz de probabilidades de transição \(\Gamma\) se, e somente se, \[\begin{equation} \pi (\mbox{I}_m-\Gamma+\text{U}) \, = \, 1, \end{equation}\] onde \(1\) é um vetor de uns, \(\mbox{I}_m\) é a matriz \(m\times m\) identidade e \(\text{U}\) é uma matriz \(m\times m\) de uns. A prova deste resultado é um exercício.


Exemplo 4.


Considerando a matriz de probabilidades de transição dada acima, podemos encontrar a distribuição estacionária como a seguir.

library(markovchain)
G = matrix(c(1/3,1/3,1/3,2/3,0,1/3,1/2,1/2,0), byrow = TRUE, ncol = 3)
ProbT = new("markovchain", states = c("1","2","3"), transitionMatrix = G, name="Exemplo")
ProbT
## Exemplo 
##  A  3 - dimensional discrete Markov Chain defined by the following states: 
##  1, 2, 3 
##  The transition matrix  (by rows)  is defined as follows: 
##           1         2         3
## 1 0.3333333 0.3333333 0.3333333
## 2 0.6666667 0.0000000 0.3333333
## 3 0.5000000 0.5000000 0.0000000
steadyStates(ProbT)
##            1       2    3
## [1,] 0.46875 0.28125 0.25
matriz = diag(3)-G+matrix(1, ncol=3, nrow=3)
solve(t(matriz), matrix(1,ncol=1,nrow=3))
##         [,1]
## [1,] 0.46875
## [2,] 0.28125
## [3,] 0.25000

Este resultado é numéricamente equivalente à \[\begin{equation} \pi \, = \, \dfrac{1}{32}\begin{pmatrix} 15 & 9 & 8 \end{pmatrix}\cdot \end{equation}\]


2.2 Reversibilidade


Uma propriedade das cadeias de Markov e outros processos aleatórios que às vezes interessa é a reversibilidade. Diz-se que um processo aleatório é reversível se suas distribuições de dimensões finitas são invariantes sob a reversão do tempo.

No caso de uma Cadeia de Markov irredutível estacionária com função de probabilidades de transição \(\Gamma\) e distribuição estacionária \(\delta\), é necessário e suficiente para a reversibilidade que as condições detalhadas de equilíbrio \[\begin{equation} \pi(i) \gamma_{i,j} \, = \, \pi(j) \gamma_{j,i} \end{equation}\] sejam satisfeitas por todos os estados \(i\) e \(j\) (Kelly, 1979).

Essas condições são satisfeitas por todas as Cadeia de Markov irredutíveis estacionárias de dois estados, que são assim reversíveis. A Cadeia Markov no exemplo na Seção I.2.1 não é reversível, no entanto, porque \[\begin{equation} \pi(1)\gamma_{1,2} \, = \, \dfrac{15}{32}\times \dfrac{1}{3} \, = \, \dfrac{5}{32} \qquad \text{mas} \qquad \pi(2)\gamma_{2,1} \, = \, \dfrac{9}{32}\times \dfrac{2}{3} \, = \, \dfrac{6}{32}\cdot \end{equation}\]


2.3 Função de autocorrelação


Teremos ocasião de comparar a função de autocorrelação (ACF) de um Modelo Oculto de Markov com a de uma Cadeia Markov. Portanto, discutiremos esta última agora. Assumimos, naturalmente, que os estados da Cadeia de Markov são quantitativos e não meramente categóricos.

A função de autocorrelação de uma Cadeia de Markov \(C_t\) em \(\{1,2,3,\cdots,m\}\) assumida estacionária e irredutível, pode ser obtida da seguinte forma.

Em primeiro lugar, definindo \(v=(1,2,\cdots,m)\) e \(V=\text{diag}(v)\), nós temos, para todos os inteiros não negativos \(k\) \[\begin{equation} \mbox{Cov}\big( C_t,C_{t+k}\big) \, = \, \pi V\Gamma^kv^\top-\big( \pi v^\top\big)^2, \end{equation}\] a prova é um Exercício. Em segundo lugar, se \(\Gamma\) é diagonalizável e os seus valores próprios, diferentes de 1, são indicados por \(\omega_2,\omega_3,\cdots,\omega_m\), então \(\Gamma\) pode ser escrito como \[\begin{equation} \Gamma \, = \, U \Omega U^{-1} \end{equation}\] onde \(\Omega=\mbox{diag}(1,\omega_2,\omega_3,\cdots,\omega_m)\) e as colunas de \(U\) são os vetores próprios correspondentes a \(\Gamma\).

Temos então, para inteiros não negativos \(k\), \[\begin{array} \mbox{Cov}\big( C_t,C_{t+k}\big) & = & \pi VU \Omega^k U^{-1}v^\top-\big( \pi v^\top\big)^2 \\ & = & a\Omega^k b^\top - a_1b_1 \, = \, \displaystyle \sum_{i=2}^m a_ib_i\omega_i^k, \end{array}\]

onde \(a=\pi VU\) e \(b^\top=U^{-1}v^\top\). Portanto \[\begin{equation} \mbox{Var}(C_t) \, = \, \sum_{i=2}^m a_ib_i \end{equation}\] e para os inteiros não negativos \(k\) \[\begin{equation} \rho(k) \, = \, \mbox{Corr}\big(C_t,C_{t+k} \big) \, = \, \dfrac{\displaystyle \sum_{i=2}^m a_ib_i\omega_i^k}{\displaystyle \sum_{i=2}^m a_ib_i}\cdot \end{equation}\] Esta é uma média ponderada da \(k\)-ésima potência dos valores próprios \(\omega_2,\omega_3,\cdots,\omega_m\) e algo semelhante ao ACF de um processo autoregressivo Gaussiano de ordem \(m-1\). Observe que a equação acima implica, no caso \(m=2\) que \(\rho(k)=\rho(1)^k\) para os inteiros não negativos \(k\) e \(\rho(1)\) é o valor próprio diferente de 1 de \(\Gamma\).


2.4 Cadeias de Markov de ordem superior


Nos casos em que as observações sobre um processo com espaço de estados finito parecem não satisfazer a propriedade de Markov, uma possibilidade que se sugere é usar uma Cadeia de Markov de ordem superior, ou seja, um modelo \(\{C_t\}\) que satisfaça a seguinte generalização da propriedade de Markov para alguns \(l\geq 2\): \[\begin{equation} P(C_t \, | \, C_{t-1}, C_{t-2}, \cdots) \, = \, P(C_t\, | \, C_{t-1}, \cdots, C_{t-l})\cdot \end{equation}\]

Um relato dessa Cadeia de Markov de ordem superior pode ser encontrado, por exemplo, em Lloyd (1980). Embora tal modelo não seja no sentido usual uma Cadeia de Markov, ou seja, não é uma Cadeia de Markov de primeira ordem, podemos redefinir o modelo de modo a produzir um processo equivalente, que sim seja.

Se consideramos \(Y_t=(C_{t-l+1},C_{t-l+2},\cdots,C_t)\), então \(Y_t\) é uma Cadeia de Markov de primeira ordem em \(S^l\), onde \(S\) é o espaço de estados de \(C_t\). Embora algumas propriedades possam ser mais difíceis de estabelecer, nenhuma teoria essencialmente nova está antes envolvida na análise de uma Cadeia de Markov de ordem superior em vez de uma cadeia de primeira ordem.

Uma Cadeia de Markov de segunda ordem, se estacionária, é caracterizada pelas probabilidades de transição \[\begin{equation} \gamma(i,j,k) \, = \, P\big( C_t=k \, | \, C_{t-1}=j, C_{t-2}=i\big) \end{equation}\] e a distribuição bivariada estacionária \[\begin{equation} u(j,k) \, = \, P\big( C_{t-1}=j, C_t=k\big) \end{equation}\] satisfazendo \[\begin{equation} u(j,k) \, = \, \sum_{i=1}^m u(i,j)\gamma(i,j,k) \qquad \mbox{e} \qquad \displaystyle \sum_{j=1}^m\sum_{k=1}^m u(j,k) \, = \, 1\cdot \end{equation}\]

Por exemplo, a Cadeia de Markov \(\{C_t\}\) de segunda ordem estacionária mais geral nos dois estados 1 e 2 é caracterizada pelas seguintes quatro probabilidades de transição: \[\begin{array}{rcl} a & = & P\big( C_t=2 \, | \, C_{t-1}=1, C_{t-2}=1\big), \\ b & = & P\big( C_t=1 \, | \, C_{t-1}=2, C_{t-2}=2\big), \\ c & = & P\big( C_t=1 \, | \, C_{t-1}=2, C_{t-2}=1\big), \\ d & = & P\big( C_t=2 \, | \, C_{t-1}=1, C_{t-2}=2\big)\cdot \end{array}\]

O processo \[\begin{equation} \{Y_t\} \, = \, \{(C_{t-1}, C_t)\} \end{equation}\] é então uma Cadeia de Markov de primeira ordem nos quatro estados \((1,1), (1,2), (2,1), (2,2)\), com uma matriz de probabilidade de transição \[\begin{equation} \begin{pmatrix} 1-a & a & 0 & 0 \\ 0 & 0 & c & 1-c \\ 1-d & d & 0 & 0 \\ 0 & 0 & b & 1-b \end{pmatrix}\cdot \end{equation}\]

Observe os zeros estruturais que aparecem nesta matriz. Não é possível, por exemplo, fazer uma transição diretamente de \((2,1)\) para \((2,2)\); daí o zero na linha 3 e na coluna 4 na matriz de probabilidades de transição acima. Os parâmetros \(a, b, c\) e \(d\) são delimitados por 0 e 1, mas de outra forma não são limitados.

A distribuição estacionária de \(\{Y_t\}\) é proporcional ao vetor \[\begin{equation} \big( b(1-d), ab, ab, a(1-c)\big), \end{equation}\] da qual se segue que a matriz \(\big(u(j,k)\big)\) de probabilidades bivariadas estacionárias para \(\{C_t\}\) é \[\begin{equation} \dfrac{1}{b(1-d)+2ab+a(1-c)}\begin{pmatrix} b(1-d) & ab \\ ab & a(1-c) \end{pmatrix}\cdot \end{equation}\]

É claro que o uso de uma Cadeia de Markov de ordem superior, em vez de uma de primeira ordem, aumenta o número de parâmetros do modelo; uma Cadeia de Markov de ordem \(l\) com \(m\) estados tem \(m^l(m-1)\) probabilidades de transição independentes. Pegram (1980) e Raftery (1985) propuseram, portanto, certas classes de modelos parcimoniosos para cadeias de ordem superior. Os modelos de Pegram têm \(m+l-1\) parâmetros e os de Raftery \(m(m-1)+l-1\). Para \(m=2\) os modelos são equivalentes, mas para \(m>2\) os de Raftery são mais gerais e podem representar uma gama mais ampla de padrões de dependência e estruturas de autocorrelação. Em ambos os casos, um aumento de um na ordem da Cadeia de Markov requer apenas um parâmetro adicional. Os modelos de Raftery, que ele denomina modelos mistura de distribuições de transiçã, definidos da seguinte forma.


Definição 1 (Cadeias de Markov de ordem superior)


O processo \(\{C_t\}\) é uma Cadeia de Markov de ordem superior ou um mistura de distribuições de transição se assume valores em \(S=\{1,2,\cdots,m\}\) e satisfaz \[\begin{equation} P\big( C_t=j_0 \, | \, C_{t-1}=j_1, \cdots, C_{t-l}=j_l\big) \, = \, \sum_{i=1}^l \lambda_i q(j_i,j_0), \end{equation}\] onde \(\sum_{i=1}^l \lambda_i=1\) e \(Q=\big( q(j,k)\big)\) é uma matriz \(m\times m\) com estradas não negativas e soma por linha igual a 1, de tal forma que o lado direito da equação acima é limitado por zero e um para todos os \(j_0,j_1,\cdots,j_l\in S\).


Este último requisito, que gera \(m^{l+1}\) pares de restrições não lineares sobre os parâmetros, assegura que as probabilidades condicionais na equação na Definição I.1 são mesmo probabilidades e a condição de que as linhas somem 1 de \(Q\) assegura que a soma sobre \(j_0\) destas probabilidades condicionais é um. Note que Raftery não assume que os parâmetros \(\lambda\) sejam não negativos.

Uma variedade de aplicações são apresentadas por Raftery (1985) e por Raftery and Tavaré (1994). Em vários dos modelos ajustados existem estimativas negativas de alguns dos coeficientes \(\lambda_i\). Para mais informações sobre esta classe de modelos, ver Haney (1993), Berechtold (2001) e Berchtold and Raftery (2002).


3 Exercícios


  1. Seja \(X\) uma variável aleatória com distribuição uma mistura de duas distribuições com esperanças \(\mu_1, \mu_2\) e variâncias \(\sigma_1^2, \sigma_2^2\), respectivamente, onde os parâmetros de mistura são \(\delta_1\) e \(\delta_2\) com \(\delta_1+\delta_2=1\).

    1. Prove que \(\mbox{Var}(X)=\delta_1\sigma_1^2+\delta_2\sigma_2^2+\delta_1\delta_2(\mu_1-\mu_2)^2\).

    2. Mostre que a mistura de duas distribuições Poisson, \(P(\lambda_1), P(\lambda_2)\), com \(\lambda_1\neq \lambda_2\), é superdispersa, ou seja, \(\mbox{Var}(X)>\mbox{E}(X)\).


  2. Uma distrinbuição de Poisson com inflação de zeros é às vezes usada como modelo para contagens sem limites exibindo sobredispersão em relação à Poisson. Tal modelo é uma mistura de duas distribuições: uma é uma Poisson e a outra é identicamente zero.

    1. É possível que um modelo deste tipo exiba uma subdispersão?

    2. Agora considere a binomial inflacionada de zeros. É possível em tal modelo que a variância seja menor que a média?


  3. Considere uma Cadeia de Markov estacionária de dois estados e matriz de probabilidades de transição dada por \[\begin{equation} \Gamma \, = \, \begin{pmatrix} \gamma_{_{1,1}} & \gamma_{_{1,2}} \\ \gamma_{_{2,1}} & \gamma_{_{2,2}} \end{pmatrix}\cdot \end{equation}\]

    1. Mostre que a distribuição estacionária é \[\begin{equation} \big(\pi(1), \pi(2)\big) \, = \, \dfrac{1}{\gamma_{_{1,2}}+\gamma_{_{2,1}}}(\gamma_{_{2,1}}, \gamma_{_{1,2}})\cdot \end{equation}\]

    2. Considere o caso \[\begin{equation} \Gamma \, = \, \begin{pmatrix} 0.9 & 0.1 \\ 0.2 & 0.8 \end{pmatrix} \end{equation}\] e as duas sequências de observações seguintes, que se supõe serem geradas pela Cadeia de Markov acima

      \[\begin{array}{ccccccc} \text{Sequência 1:} & 1 & 1 & 1 & 2 & 2 & 1 \\ \text{Sequência 2:} & 2 & 1 & 1 & 2 & 1 & 1 \end{array}\]

      Calcule a distribuição estacionária de cada uma das sequências. Note que cada sequência contém o mesmo número de uns e dois. Porquê estas sequências não são igualmente prováveis?


  4. Considere uma Cadeia de Markov estacionária de dois estados e matriz de probabilidades de transição dada por \[\begin{equation} \Gamma \, = \, \begin{pmatrix} \gamma_{_{1,1}} & \gamma_{_{1,2}} \\ \gamma_{_{2,1}} & \gamma_{_{2,2}} \end{pmatrix}\cdot \end{equation}\] Mostre que a \(k\)-ésima matriz de transição, ou seja, \(\Gamma^k\) é dada por \[\begin{equation} \Gamma^k \, = \, \begin{pmatrix} \pi(1) & \pi(2) \\ \pi(1) & \pi(2) \end{pmatrix}+\omega^k \begin{pmatrix} \pi(2) & -\pi(2) \\ -\pi(1) & \pi(1) \end{pmatrix}, \end{equation}\] onde \(\omega=1-\gamma_{_{1,2}}-\gamma_{_{2,1}}\) e \(\pi(1)\) e \(\pi(2)\) foram definidos no exercício anterior. Uma maneira de mostrar isto é diagonalizar a matriz de probabilidade de transição.

    1. Esta é uma das várias abordagens possíveis para encontrar a distribuição estacionária de uma Cadeia de Markov; proposta por Grimmett and Stirzaker (2001). Suponha que \(\Gamma\) seja a matriz de probabilidades de transição de uma Cadeia de Markov homogênea com \(m\) estados e seja \(\pi\) um vetor com \(m\) componentes todos não negetivos. Mostre que \(\pi\) é a distribuição estacionária da Cadeia de Markov se, e somente se, \[\begin{equation} \pi (\mbox{I}_m-\Gamma+\text{U}) \, = \, 1, \end{equation}\] onde \(1\) é um vetor de uns, \(\mbox{I}_m\) é a matriz \(m\times m\) identidade e \(\text{U}\) é uma matriz \(m\times m\) de uns.

    2. Encontre a distribuição estacionária correspondente às seguintes matrizes de probabilidades de transição. Alguma delas causa problema?

    1. \[\begin{pmatrix} 0.7 & 0.2 & 0.1 \\ 0 & 0.6 & 0.4 \\ 0.5 & 0 & 0.5 \end{pmatrix}\]
    2. \[\begin{pmatrix} 0 & 1 & 0 \\ 1/3 & 0 & 2/3 \\ 0 & 1 & 0 \end{pmatrix}\]
    3. \[\begin{pmatrix} 0 & 0.5 & 0 & 0.5 \\ 0.75 & 0 & 0.25 & 0 \\ 0 & 0.75 & 0 & 0.25 \\ 0.5 & 0 & 0.5 & 0 \end{pmatrix}\]
    4. \[\begin{pmatrix} 0.25 & 0.25 & 0.25 & 0.25 \\ 0.25 & 0.25 & 0.5 & 0 \\ 0 & 0 & 0.25 & 0.75 \\ 0 & 0 & 0.5 & 0.5 \end{pmatrix}\]
    5. \[\begin{pmatrix} 1 & 0 & 0 & 0 \\ 0.5 & 0 & 0.5 & 0 \\ 0 & 0.75 & 0 & 0.25 \\ 0 & 0 & 0 & 1 \end{pmatrix}\]

  6. Sejam as quantidades \(a_i\) não-negativas, com \(\sum_i a_i>0\). Utilizando os multiplicadoes de Lagrange, maximize \[\begin{equation} S \, = \, \sum_{i=1}^m a_i \log\big( \delta_i\big), \end{equation}\] sobre \(\delta_i\geq 0\), restrito a \(\sum_i \delta_i=1\). Verifique as condições que tanto a primeira como a segunda derivadas devem satisfazer.

  7. Exemplo em Bisgaard and Travis (1991). Considere a seguinte sequência de 21 observações, assumidas como resultantes de uma Cadeia de Markov homogênea de dois estados \[\begin{array}{ccccccccccccccccccccccccc} 1 & 1 & 1 & 0 & 1 & & 1 & 0 & 1 & 1 & 1 & & 1 & 0 & 1 & 1 & 0 & & 1 & 1 & 1 & 1 & 1 & & 1\cdot \end{array}\]

    1. Estimar a matriz de probabilidades de transição por máxima verossimilhança condicionada à primeira observação.

    2. Estimar a matriz de probabilidades de transição por máxima verossimilhança incondicional assumindo estacionaridade da Cadeia de Markov.

    3. Use as funções R, contour e persp para produzir gráficos de contorno e perspectiva do logaritmo da verossimilhança incondicional, como uma função das duas probabilidades de transição fora da diagonal principal.


  8. Considere as duas seguintes matrizes de probabilidades de transição, nenhuma das quais é diagonalizável:

    1. \[\begin{pmatrix} 1/3 & 1/3 & 1/3 \\ 2/3 & 0 & 1/3 \\ 1/2 & 1/2 & 0 \end{pmatrix}\]
    2. \[\begin{pmatrix} 0.9 & 0.08 & 0 & 0.02 \\ 0 & 0.7 & 0.2 & 0.1 \\ 0 & 0 & 0.7 & 0.3 \\ 0 & 0 & 0 & 1 \end{pmatrix}\]

    Em cada caso, escreva \(\Gamma\) na forma canônica de Jordan e assim encontre uma expressão explícita para as probabilidades de transição em \(t\)-passos, \(t=1,2,\cdots\).


  9. Exemplo em Singh (2003). Considere a seguinte, muito curta, sequência de ADN: \[\begin{array}{ccccccccc} AACGT & & CTCTA & & TCATG & & CCAGG & & ATCTG \end{array}\]

    Ajuste uma Cadeia de Markov homogênea a estes dados por:

    1. Máxima verossimilhança condicionada à primeira observação;

    2. Assumindo estacionariedade e maximizando a verossimilhança não condicional de todas as 25 observações.