jueves, 27 de octubre de 2016

25-Cadenas de Markov 6

Cadenas de Markov 6 de 10

Estados recurrentes y transitorios.

Supongamos una cadena de Markov con espacio temporal \(\left\{0,1,2,\ldots\right\}\) y con un espacio de estados discreto \(\mathcal{S}\). Un estado \(i \in \mathcal{S}\) es recurrente si su probabilidad de un eventual regreso es uno. Es decir cuando

\[ P(X_n=i \text{ para alguna } n = 1,2,\ldots\mid X_0 = i) = 1.\]

Además cuando un estado no sea recurrente diremos que es un estado transitorio.

Una manera equivalente pero más útil de definir un estado recurrente es la siguiente.

Definición 1.

Se dice que un estado es recurrente si \[f_{ii} = 1,\] donde \[f_{ii} = \sum_{n=1}^\infty f_{ii}(n),\] y \[f_{ii}(n) = P(X_n = i, X_{n-1}\neq i, \ldots, X_{1}\neq i\mid X_0=i).\]

Tiempo medio de recurrencia.

Cuando un estado \(i\) es recurrente tenemos dos casos, que en promedio, la cadena regrese a \(i\) en un tiempo finito o en un tiempo infinito.

Definición 2.

Definimos el tiempo de regreso a un estado \(i\) como la variable

\[T_i = \min\left\{ n\ge 1: X_n = i \text{ suponiendo que } X_0 = i\right\}.\] Si el conjunto \(\left\{n\ge 1: X_n = i \text{ suponiendo que } X_0 = i\right\} = \emptyset\) entonces definimos \(T_i = \infty\). Llamamos tiempo medio de recurrencia a la siguiente esperanza: \[\mu_i := E(T_i) = \sum_{n = 1}^\infty n P(T_i = n) = \sum_{n = 1}^\infty nf_{ii}(n).\] Se dice que un estado \(i\) es recurrente positivo si cumple ser recurrente y además \[ \mu_i < \infty,\] en otro caso decimos que el estado \(i\) es recurrente nulo.

Se puede demostrar que basta tener una cadena cuya matriz de probabilidades de transición sea finita e irreducible para concluir que todos los estados son recurrentes positivos, el lector interesado puede ver la demostración, por ejemplo, en las páginas 66 a 69 del libro de Luis Rincón Introducción a los procesos estocásticos http://lya.fciencias.unam.mx/lars/Publicaciones/procesos2012.pdf.

Proposición 1 (F + I = RP).

Si una cadena de Markov es finita e irreducible, entonces es recurrente positiva, es decir, todos sus estados son recurrentes positivos.

Cálculo de \(\mu_i\) por simulación.

Supongamos que una cadena de Markov tiene matriz de probabilidades de transición \(P\), entonces con el siguiente código se calculan aproximaciones a los tiempos medios de recurrencia, que por la Proposición 1, sabemos que son cantidades finitas.

## La función tm_recurrencia aproxima el valor
## de mu_i (tiempo medio de recurrencia) simulando
## muchas cadenas y contando el número de pasos que
## tardan en regresar al estado inicial.
## 
## X0 es el estado inicial
## P la matriz de probabilidades de transición
## n el número de cadenas a simular
## estadocero = TRUE si el espacio de estados inicia
## en cero y estadocero = FALSE si comienza en 1

tm_recurrencia <- function(X0,P,n,estadocero = T){
        dim <- length(P[2,]) 
        tau <- numeric(n)
        if(estadocero == T){
                X0 <-  X0 + 1
        }else{ X0 <- X0}
        for(x in 1:n){
                aux <- 0
                cadena <- X0
                k <- 1
                while(aux != X0){
                        cadena[k+1] <- sample(1:dim,1,r = T,P[cadena[k],])
                        if(cadena[k+1] == X0){
                                tau[x] <- k
                                aux <- X0
                        }else{
                                k <- k+1 
                                aux <- 0
                        }
                }
        }
        mean(tau)
}

## Ejemplo Cadena de Ehrenfest.
p0 <- c(0,1,0,0,0,0,0,0,0,0)
p1 <- c(1/9,0,8/9,0,0,0,0,0,0,0)
p2 <- c(0,2/9,0,7/9,0,0,0,0,0,0)
p3 <- c(0,0,3/9,0,6/9,0,0,0,0,0)
p4 <- c(0,0,0,4/9,0,5/9,0,0,0,0)
p5 <- c(0,0,0,0,5/9,0,4/9,0,0,0)
p6 <- c(0,0,0,0,0,6/9,0,3/9,0,0)
p7 <- c(0,0,0,0,0,0,7/9,0,2/9,0)
p8 <- c(0,0,0,0,0,0,0,8/9,0,1/9)
p9 <- c(0,0,0,0,0,0,0,0,1,0)
P <- rbind(p0,p1,p2,p3,p4,p5,p6,p7,p8,p9)
colnames(P) <- (0:9)
rownames(P) <- (0:9)
library(MASS)
fractions(P)
##   0   1   2   3   4   5   6   7   8   9  
## 0   0   1   0   0   0   0   0   0   0   0
## 1 1/9   0 8/9   0   0   0   0   0   0   0
## 2   0 2/9   0 7/9   0   0   0   0   0   0
## 3   0   0 1/3   0 2/3   0   0   0   0   0
## 4   0   0   0 4/9   0 5/9   0   0   0   0
## 5   0   0   0   0 5/9   0 4/9   0   0   0
## 6   0   0   0   0   0 2/3   0 1/3   0   0
## 7   0   0   0   0   0   0 7/9   0 2/9   0
## 8   0   0   0   0   0   0   0 8/9   0 1/9
## 9   0   0   0   0   0   0   0   0   1   0
###########
## Cálculo de valores exactos (véase página 79 de libro de Luis Rincón)
###########
mu <- 0
mu_0 <- 2^9
for(i in 1:9){
        mu[i] <- mu_0/choose(9,i)
}
mu <- c(mu_0,mu)

###########
## Cálculo de valores por simulación usando la función
##         tm_recurrencia
###########
set.seed(123) #para reproducibilidad
mu_sim <- 0
for(i in 0:9){
        mu_sim[i+1] <- tm_recurrencia(i,P,10000)
}

Exactos <- round(mu,2)
Por_simulacion <- round(mu_sim,2)
Diferencia <- round(abs(mu-mu_sim),2)

Resultados <- cbind(Exactos,Por_simulacion, Diferencia)
aux <- paste('mu_',0)
for(i in 1:9){
        aux[i+1] <- paste('mu_',i) 
}
row.names(Resultados) <- aux
Resultados
##       Exactos Por_simulacion Diferencia
## mu_ 0  512.00         518.28       6.28
## mu_ 1   56.89          57.35       0.46
## mu_ 2   14.22          14.21       0.01
## mu_ 3    6.10           6.12       0.02
## mu_ 4    4.06           4.05       0.02
## mu_ 5    4.06           4.13       0.06
## mu_ 6    6.10           6.10       0.00
## mu_ 7   14.22          13.88       0.34
## mu_ 8   56.89          57.54       0.65
## mu_ 9  512.00         511.31       0.69

No hay comentarios:

Publicar un comentario