CNSF

Procesos Estocásticos I Capacitación técnica especializada en el nuevo marco de Solvencia Gerónimo Uribe Bravo Instituto de Matemáticas Universidad Nacional Autónoma de México CAPÍTULO 1 Introducción 1. Definiciones y clasificación básica de procesos estocásticos Un proceso estocástico es una colección de variables aleatorias (Xt )t∈T indexadas por un conjunto T y definidas en algún espacio de probabilidad (Ω, F , P). Interpretamos al conjunto de ı́ndices T como un parámetro temporal; para nosotros T será {0, . . . , n}, N, algún intervalo [0, t] ó [0, ∞). Interpretamos a un proceso estocástico como la evolución en el tiempo de algún fenómeno cuya dinámica se rige por el azar. Un ejemplo sencillo de esto es la cantidad de soles que vamos acumulando al participar en un juego de volados. Otro ejemplo es la evolución en el tiempo de la reserva de una compañı́a de seguros. En el primer ejemplo, se puede indexar al proceso por algún intervalo de naturales, en cuyo caso hablaremos de un proceso estocástico a tiempo discreto. Además, dicho proceso toma valores en los naturales, por lo que también se trata de un proceso con espacio de estados discreto. En el segundo caso, se puede pensar en un modelo indexado por un subintervalo de [0, ∞) y hablaremos de un proceso estocástico a tiempo continuo. Además, en principio el valor de la reserva podrı́a ser cualquier real no-negativo y por lo tanto hablamos de un proceso con espacio de estados continuo Uno de los primeros resultados generales dentro de la teorı́a de los procesos estocásticos es el teorema de consistencia de Kolmogorov que nos permite construir procesos estocásticos a partir de colecciones vectores aleatorios (que satisfacen la condición técnica de ser consistentes). La prueba de este teorema se puede hacer basándose en la existencia de una sucesión de variables aleatorias uniformes. Antes de analizar por qué existe una sucesión de variables uniformes independientes, ejemplificaremos cómo se pueden construir algunos de los procesos estocásticos que analizaremos en este curso. Ejemplo 1.1 (Caminatas aleatorias simples y el problema de la ruina). Imaginemos la siguiente situación: tengo un capital de 20 pesos al tiempo cero y cada instante de tiempo apuesto un peso en un volado, ganándo si cae águila. ¿cómo puedo estudiar matemáticamente a la evolución de mi capital en el tiempo? De particular interés es la variable aleatoria que nos indica el instante en que me arruino por primera vez, misma que a priori podrı́a ser infinita si jamás me arruino. 1 1. Definiciones básicas 2 El modelo matemático es el siguiente: consideremos variables aleatorias U 1 , U2 , . . . uniformes en (0, 1) e independientes. A la variable aleatoria 1Ui ≤1/2 , que toma los valores cero y uno, la interpretaremos como indicándonos si el resultado del i-ésimo volado es águila (cuando toma el valor uno) y por lo tanto, la variable 21Ui ≤1/2 − 1 toma los valores 1 si cae águila y −1 si cae sol. Ejercicio 1.1. Con el objeto de verificar que comprendemos la noción de independencia, probar que las variables aleatorias 1U1 ≤1/2 , 1U2 ≤1/2 , . . . son independientes y con distribución Bernoulli de parámetro 1/2. Finalmente, podemos definir X0 = 20 y Xn+1 = Xn + 21Un+1 ≤1/2 − 1. El siguiente código en R simula la evolución de mi fortuna. C <- 20 # C es un vector cuya entrada i será mi capital al tiempo i aux <-C # Esta variable me dice cuál es el último valor de mi capital while ( aux >0) { # Mientras no me haya arruinado aux <- aux +2 * ( runif (1) <1 / 2) -1 # actualizo mi capital al sumarle una variable que toma valores -1 y 1 con probabilidad 1 / 2 C <-c (C , aux ) # Agrego el último valor de mi fortuna al vector C } plot ( C ) Listing 1.1. Ruina.R En la Figura 1 podemos apreciar un ejemplo del resultado de correr el código anterior. Ejemplo 1.2 (Apostando con prisa). Modificaremos el ejemplo anterior como sigue: tengo un capital de 20 pesos al tiempo cero y cada instante de tiempo apuesto en un volado ya sea la mitad de mi fortuna si tengo más de 6 pesos o 2 pesos si mi fortuna es menor o igual a 6, ganando si cae águila. Un modelo matemático es el siguiente: consideremos variables aleatorias U1 , U2 , . . . uniformes en (0, 1) e independientes y definamos (
⌊Xn /2⌋ Xn > 6 X0 = 20 y Xn+1 = Xn + 2 ∗ 1Un+1 ≤1/2 − 1 . 2 Xn ≤ 6 El modelo anterior se puede implementar fácilmente en R con el códio siguiente. C <- 20 # C es un vector cuya entrada i será mi capital al tiempo fortuna <-C # Esta variable me dice cuál es el último valor de mi while ( fortuna >0) { # Mientras no me haya arruinado monto <-2 * ( fortuna <=6) + floor ( fortuna / 2) * ( fortuna >6) # Calculo apostaré , que es la mitad de mi fortuna cuando tengo más si no , dos pesos . i capital el monto que de 6 pesos y 3 0 5 10 C 15 20 1. Definiciones básicas 0 20 40 60 80 100 Index Figura 1. Trayectoria de una caminata aleatoria simple que comienza en 20 y es detenida al llegar a cero fortuna <- fortuna + monto * (2 * ( runif (1) >1 / 2) -1) # actualizo mi capital al sumarle una variable que toma valores -1 y 1 con probabilidad 1 / 2 multiplicada por el monto de la apuesta C <-c (C , fortuna ) # Agrego el último valor de mi fortuna al vector C } plot ( C ) Listing 1.2. Prisa.R Por supuesto, esperamos que esta estrategia nos llege más rápido a la ruina. En la Figura 2 podemos apreciar dos ejemplos de trayectorias simuladas de la evolución de la fortuna bajo este esquema de apuestas. Los dos ejemplos anteriores corresponden a procesos con tiempo y espacio discreto. Ahora analizaremos un modelo a tiempo continuo y espacio discreto. Ejemplo 1.3 (Contéos aleatorios ). Imaginemos que queremos modelar los tiempos sucesivos en que cambiamos un foco en nuestro lugar de trabajo. Supondremos que en cuanto se funde un foco lo cambiamos (instantáneamente) por uno nuevo. Es natural asumir que podemos modelar los tiempos de vida de los sucesivos focos mediante una sucesión de variables aleatorias independientes. El supuesto 4 C 0 0 10 50 20 100 30 C 40 150 50 200 60 250 1. Definiciones básicas 2 4 6 8 10 12 0 20 40 Index 60 80 Index Figura 2. Dos trayectorias que muestran la evlución de un capital al someterlo a un esquema de apuestas arriesgadas adicional que impondremos es que éstas tienen distribución exponencial de tasa λ > 0, donde la tasa es el recı́proco de la media. Sean U1 , U2 , . . . variables independientes de distribución uniforme en (0, 1) y definamos a Ti = − log(Ui ) /λ. Entonces T1 , T2 , . . . son variables aleatorias exponenciales independientes de tasa λ. La variable Ti la interpretamos como el tiempo de vida del i-ésimo foco. Puesto que la distribución exponencial está caracterizada por la propiedad de pérdida de memoria P(Si > t + s |Si > t) = P(Si > s) , el suponer que el tiempo de vida de un foco tiene distribución exponencial puede ser cuestionable puesto que debe haber un efecto de desgaste en su tiempo de vida. Sin embargo, lo que se espera con el modelo es que capture la escencia del fenómeno que queremos modelar. El proceso estocástico de interés es el que va midiendo la cantidad de focos que hemos cambiado en el intervalo de tiempo [0, t] que mediremos en años. Sean T0 = 0, Tn+1 = Tn + Sn+1 y Nt = ∞ X 1Ti ≤t . i=1 Entonces se interpreta a Tn como el instante de tiempo en el que cambiamos el n-ésimo foco y a Nt como la cantidad de focos que hemos cambiado en [0, t]. Se puede simular a la función aleatoria t 7→ Nt en el intervalo [0, 1] mediante el siguiente código. lambda =24 * 360 / 1000 xi = rexp (1 , lambda ) foco # Media , en a \ ~ nos , del tiempo de vida de un foco # xi representa el tiempo en el que se cambió el último 5 0 2 20 40 c(1:(N + 2)) 8 6 4 c(1:(N + 2)) 60 10 12 80 1. Definiciones básicas 0.0 0.2 0.4 0.6 0.8 1.0 0.0 0.2 T 0.4 0.6 0.8 1.0 T Figura 3. Los tiempos sucesivos en los que se cambia un foco a lo largo de 1 y 8 años T = c (0 , xi ) # El vector T irá acumulando los tiempos en que vamos cambiando los focos N =0 # N nos dirá cuantos focos hemos cambiado al final de un a~ no while ( xi <1) { # Mientras no haya pasado un a~ no N <-N +1 # Aumentamos el número de focos cambiados en uno xi <- xi + rexp (1 , lambda ) # Vemos el tiempo en el que debemos cambiar el siguiente foco T = c (T , xi ) # Aumentamos un evento temporal } plot (T , c (1:( N +2) ) ) Listing 1.3. Poisson.R En la Figura 3 podemos observar dos trayectorias con los tiempos en los que se cambian los focos, primero en un año y luego en ocho. Podemos observar una cierta regularidad, como si hubiera cierta tendencia determinista (una linea recta) y unas fluctuaciones aleatorias que capturan los sucesivos tiempos de cambio de foco. Ejemplo 1.4 (Tiempos de espera). Imaginemos la fila de un banco. Supongamos que los clientes van llegando a tiempos aleatorios y que cada uno requiere un servicio que es también una variable aleatoria. Lo que se quiere medir es: al correr el sistema, si un cliente llega al momento t, ‘?Cuánto debe esperar para salir del banco? Un modelo posible se enfoca en los tiempos entre los arribos de los clientes y supone que éstos son variables aleatorias exponenciales de algún parámetro λi . Además, podemos suponer que los tiempos de servicio son variables aleatorias independientes con distribución común, que fijaremos como la exponencial de tasa 1. Definiciones básicas 6 λs para fijar ideas. Se supone que todas las variables en cuestón son independientes. Este modelo, aunque sea criticable en su supuesto de pérdida de memoria heredado de las variables exponenciales, tiene la particularidad de que se pueden hacer cálculos explı́citos que no son posibles en modelos más generales. Además, ejemplifica algunas caracterı́sticas de los modelos más generales. Sean S1 , S2 , . . . variables exponenciales independientes de parámetro λi y ξ1 , ξ2 , . . . variables exponenciales independientes (entre si y de las Si ) de parámetro λs . Si T0 = 0, Tn+1 = Tn + Sn+1 , Nt = ∞ X 1Tn ≤t , R0 = 0 y Rn+1 = Rn + ξn+1 , n=1 definimos entonces los procesos X t = RN t − t y Qt = Xt − min Xs . s≤t Entonces Qt representa el tiempo de servicio necesario para atender a los clientes que se encuentran presentes en el banco al tiempo t. Por otra parte, podemos simular al modelo matemático de la cola mediante el siguiente código en R. li <-1 # Tasa interarribo ls <-2 # Recı́proco de la media de servicio T =8 * 60 # Tiempo de la simulación t <-c (0 , rexp (1 , li ) ) # Inicia lización del vector de eventos temporales q <-c (0 , rexp (1 , ls ) ) # Inicia lización del vector de estado de la cola while ( tail (t ,1) <T ) { # Mientras no haya sobrepasado el umbral temporal taux = rexp (1 , li ) # Me fijo en cuánto falta para la llegada del próximo cliente if ( taux < tail (q ,1) ) { # En particular si el próximo cliente llega antes de que la cola se vacı́e t <-c (t , tail (t ,1) + taux ) # En cuyo caso agrego el evento de llegada q <-c (q , tail (q ,1) - taux + rexp (1 , ls ) ) # Junto con el tiempo de servicio que requiere menos el que ya he realizado } else { # Si el próximo cliente llega después de que la cola se vacı́e t <-c (t , tail (t ,1) + tail (q ,1) , tail (t ,1) + taux ) # Agrego dos eventos temporales : cuando se vacı́a la cola y cuando llega el próximo cliente q <-c (q ,0 , rexp (1 , ls ) ) # Agrego además un estado de cola =0 más el servicio del próximo cliente que llega } } plot (t , q ) Listing 1.4. Cola.R Al ejecutar el código se obtienen gráficos como los de la Figura 4. 7 0 0 5 1 10 2 15 q q 3 20 25 4 30 5 2. Variables uniformes independientes 0 100 200 300 t 400 0 100 200 300 400 t Figura 4. Estado de la cola cuando λi = 1 y λs = 1, 2 2. La construcción fundamental de una sucesión de variables aleatorias independientes Como vimos en los ejemplos anteriores, y es cierto en gran generalidad, podemos construir procesos estocásticos muy generales a partir de sucesiones de variables aleatorias inependientes. En cierto sentido, dichas sucesiones son los ejemplos más sencillos de procesos estocásticos, en los que no hay realmente una evolución. Al ser, sin embargo, los bloques fundamentales con los que se construyen todos los demás, nos detendremos en su construcción matemática. 2.1. El modelo matemático de una sucesión de volados. Primero se ejemplificará la construcción de una sucesión de variables aleatorias independientes a partir de una sola variable. Consideremos al espacio de probabilidad (Ω, F , P) en el que Ω = (0, 1], F = B(0,1] y P es la medida de Lebesgue restringida a Ω. Definamos a las variables aleatorias dn : Ω → R como sigue: a cada ω ∈ Ω, se le puede asignar su expansión diádica con colas infinitas de tal forma que ∞ X dn (ω) , ω= 2n n=1 donde cada dn es cero o uno. Aunque la expansión diádica todavı́a no esté bien definida, puesto que por ejemplo a 1/2 se le podrı́a asociar ya sea (1, 0, 0, . . .) o (0, 1, 1, . . .), la expansión diádica con colas infinitas, que es la segunda en nuestro ejemplo, sı́ lo está. Más formalmente, definamos ( 0 si ω ∈ (0, 1/2] d1 (ω) = . 1 si ω ∈ (1/2, 1] 2. Variables uniformes independientes 8 Notemos que si ω1 = 2ω − d1 (ω), entonces ω1 ∈ (0, 1]; recursivamente, definimos dn+1 (ω) = d1 (ωn ) y ωn+1 = 2ωn − d1 (ωn ) ∈ (0, 1]. Es fácil ver que de hecho, d2 (ω) = 1(1/4,2/4] + 1(3/4,4/4] , d3 (ω) = 1(1/8,2/8] + 1(3/8,4/8] + 1(5/8,6/8] + 1(7/8,8/8] y en general dn (ω) = n−1 2X 1((2i−1)/2n ,2i/2n ] . i=1 Esto implica inmediatamente que si u1 , . . . , un ∈ {0, 1} entonces el conjunto {d1 = u1 , . . . , dn = un } es un intervalo de longitud 1/2n y que por lo tanto d1 , . . . , dn son variables aleatorias independientes de distribución Bernoulli de parámetro 1/2. 2.2. Una sucesión de variables aleatorias uniformes independientes. Ahora demostraremos que si (Xi )i≥1 son variables aleatorias independientes con distribución Bernoulli de parámetro 1/2 (definidas en otro espacio de probabilidad) P entonces U = i≥1 Xi /2i tiene distribución uniforme. En efecto, puesto que P(X1 = u1 , . . . , Xn = un ) = 1/2n = P(d1 = u1 , . . . , dn = un ) , vemos que P n X i=1 Por otro lado, i Xi /2 < x ! U = lim n =P n X i di /2 < x . i=1 n X ! Xi /2i , i=1 y de hecho la sucesión de variables aleatorias es creciente. Por lo tanto U es variable aleatoria y ! ! n n X X di /2i < x = P((0, x]) = x. Xi /2i < x = lim P P(U < x) = lim P n→∞ i=1 n→∞ i=1 Ası́, vemos que U es una variable uniforme. Ahora utilizaremos lo anterior para mostrar que existe un espacio de probabilidad en el que están definidas una sucesión de variables aleatorias uniformes independientes. De hecho el espacio de probabilidad que consideraremos es el mismo (Ω, F , P) que en la Subsección 2.1. Como Z+ y Z2+ tienen la misma cardinalidad, 2. Variables uniformes independientes 9 consideremos una biyección de φ : Z2+ → Z+ . Definamos dni = dφ(n,i) y para cada n ∈ Z+ , sea X dn i . Un = 2i i≥1 (dni )i≥1 Como son variables aleatorias independientes de distribución Bernoulli de parámetro 1/2, se sique que Un tiene distribución uniforme para cada n ∈ N. Se afirma ahora que las variables (Un )n≥1 son independientes. En efecto, esto es consecuencia del siguiente lema, un tanto más general. Notemos que Un es medible respecto de la σ-álgebra generada por (dni )i≥1 , a la cual llamaremos Fn . Lema 1. Sean Fn,i , i ≥ 1, n ≥ 1 σ-álgebras independientes y definamos Fn = σ(Fn,1 , Fn,2 , . . .) . Entonces Fn , n ≥ 1 son σ-álgebras independientes. Demostración. Debemos mostrar que para todo A1 ∈ F1 , . . . , An ∈ Fn , se tiene que P(A1 ∩ · · · ∩ An ) = P(A1 ) · · · P(An ) . (1) Sea Cn = {A1 ∩ · · · ∩ Am : m ≥ 1 y Aj ∈ ∪i≥1 Fn,i para j = 1, . . . , m} . Puesto que [ i≥1 vemos que Fn,i ⊂ Cn ⊂ Fn , σ(Cn ) = Fn . Por otra parte, es fácil ver que Cn es un π-sistema. Consideremos ahora la clase M1 = {A ∈ Fk : P(A ∩ B) = P(A) P(B) si B = B1 ∩ · · · ∩ Bn con Bj ∈ Cj } . Es fácil ver que Mk es un λ-sistema que contiene, por hipótesis a C1 . Por lo tanto M1 = F 1 . Ahora consideramos a M2 = {A ∈ F2 : P(A ∩ B) = P(A) P(B) si B = B1 ∩ B3 ∩ · · · ∩ Bn con B1 ∈ F1 y Bj ∈ Cj para j ≥ 3} . Se prueba entonces que M2 es un λ-sistema que por hipótesis y la igualdad M1 = F1 contiene a C2 . Al aplicar este razonamiento sucesivamente, obtenemos la igualdad (1).  2. Variables uniformes independientes 10 2.3. Una sucesión de variables aleatorias independientes con distribuciones arbitrarias. Ahora utilizaremos la construcción de la sucesión de variables aleatorias uniformes independientes para demostrar el siguiente resultado: Teorema 1.1. Sean µn , n ≥ 1 medidas de probabilidad en R. Entonces existe un espacio de probabilidad (Ω, F , P) y una sucesión de variables aleatorias independientes Xn : Ω → R, n ≥ 1 tales que la distribución de Xn es µn . La herramienta principal de la construcción será el siguiente lema: (tomado de Billingsley p. 190). Recordemos que una función F : R → [0, 1] es la función de distribución de una variable aleatoria real si y sólo si es no decreciente, continua por la derecha (y con lı́mites por la izquierda) tal que limx→−∞ F (x) = 0 y limx→∞ F (x) = 1. Definición. La función de cuantiles de una función de distribución F es la función φ : (0, 1) → R dada por φ(u) = inf {x ∈ R : u ≤ F (x)} . La función de cuantiles satisface la igualdad φ(u) ≤ x ⇔ u ≤ F (x) que se demostrará posteriormente. De esta igualdad se deduce la medibilidad de φ. Prueba del Teorema 1.1. Sabemos que existe un espacio de probabilidad (Ω, F , P) en el que existe una sucesión de variables aleatorias independientes (Un )n≥1 uniformes en (0, 1). Sea Fn la función de distribución asociada a la medida de probabilidad µn y φn la función de cuantiles de Fn . Como φn es una función medible, Xn = φn (Un ) es una variable aleatoria. Además, como las variables Un , 6= 1 son independientes, también lo son las variables Xn , n ≥ 1:  −1 {X1 ∈ A1 , . . . , Xn ∈ An } = U1 , ∈ φ−1 1 (A1 ) , . . . , Un , ∈ φn (An ) . Finalmente:



P Xi−1 ((−∞, x]) = P Ui−1 φ−1 = P Ui−1 ((0, Fi (x)]) = Fi (x) , i ((−∞, x]) por lo que Xi tiene distribución µi .  Ahora demostremos las propiedades de la función de cuantiles φ asociada a la función de distribución F . Sea u ∈ (0, 1); entonces el conjunto {x ∈ R : u ≤ F (x)} es no vacı́o y como F es no decreciente, es un intervalo ya sea de la forma [φ(u) , ∞) o (φ(u) , ∞), ya que φ(u) es el ı́nfimo del conjunto considerado. La segunda opción se descarta al notar que F es continua por la derecha. Por lo tanto, u ≤ F (x) si y sólo si φ(u) ≤ x. CAPÍTULO 2 Cadenas de Markov a tiempo discreto La Figura 1 representa un laberinto. Imaginemos que colocamos una rata en la esquina inferior izquierda y un plato de comida en la esquina superior derecha. Para modelar la trayectoria que sigue la rata hasta encontrar la comida, supongamos que cuando se encuentra en un cuarto del laberinto, la rata va a cualquier otro con la misma probabilidad. Un modelo matemático para esta situación es el siguiente. Enumeremos los cuartos del laberinto de izquierda a derecha, de abajo a arriba, por lo que la rata comienza en el cuarto 1 y encuentra la comida en el cuarto 9. Definamos Pi,j como la probabilidad con que la rata pasa del cuarto i al cuarto j; por ejemplo, vemos que P5,j = ( 1/4 j ∈ {2, 4, 6, 8} . 0 j∈ 6 {2, 4, 6, 8} Figura 1. Laberinto para un experimento aleatorio 11 12 Esta información se puede organizar de  0 1/2 0 1/2 1/3 0 1/3 0   0 1/2 0 0  1/3 0 0 0  0 1/4 0 1/4 P =   0 0 1/3 0   0 0 0 1/2   0 0 0 0 0 0 0 0 forma matricial de la siguiente manera:  0 0 0 0 0 1/3 0 0 0 0   0 1/2 0 0 0   1/3 0 1/3 0 0   0 1/4 0 1/4 0   1/3 0 0 0 1/3  0 0 0 1/2 0   1/3 0 1/3 0 1/3 0 1/2 0 1/2 0 Entonces la probabilidad de que la rata siga la trayectoria 1, 4, 7, 8, 9 para encontrar la comida es P1,4 P4,7 P7,8 P8,9 . Notemos que (1) P Pi,j ≥ 0 para todo i y j y (2) j Pi,j = 1 para toda i. A una matriz con estas dos caracterı́sticas se le llama matriz estocástica. La segunda condición nos dice que Pi,1 , . . . , Pi,9 es una distribución de probabilidad sobre el conjunto {1, . . . , 9}. Si definimos a φi como la función de quantiles asociada, tendremos que φi (Uj ) es una variable aleatoria con la misma distribución que el cuarto al que pasa la rata si está en el cuarto j. Es por esto que si definitmos X0 = 1 y Xn+1 = φXn (Un+1 ) , las variables X0 , X1 , . . . nos modelan el movimiento de la rata por el laberinto. Para poder obtener la trayectoria de la rata detenida hasta que encuentre la comida, podrı́amos modificar la matriz P en   0 1/2 0 1/2 0 0 0 0 0 1/3 0 1/3 0 1/3 0 0 0 0     0 1/2 0 0 0 1/2 0 0 0    1/3 0 0 0 1/3 0 1/3 0 0    1/4 0 1/4 0 1/4 0 1/4 0  P̃ =   0   0 0 1/3 0 1/3 0 0 0 1/3    0 0 0 1/2 0 0 0 1/2 0     0 0 0 0 1/3 0 1/3 0 1/3 0 0 0 0 0 0 0 0 1 que difiere de P salvo en el último renglón, en el cual especificamos que una vez que la rata llegue al cuarto 9 se quede ahı́. Si φ̃i son las funciones de cuantiles asociadas a los renglones de P̃ , podemos entonces modelar la trayectoria de la rata, detenida cuando alcanza la comida mediante la sucesión X̃0 = 1 y X̃n+1 = φ̃X̃n (Un+1 ) . 13 Se presenta a continuación un código en R para simular la trayectoria de la rata. P = matrix ( c (0 ,1 / 2 ,0 ,1 / 2 ,0 ,0 ,0 ,0 ,0 ,1 / 3 ,0 ,1 / 3 ,0 ,1 / 3 ,0 ,0 ,0 ,0 ,0 ,1 / 2 ,0 ,0 ,0 ,1 / 2 ,0 ,0 ,0 ,1 / 3 ,0 ,0 ,0 ,1 / 3 ,0 ,1 / 3 ,0 ,0 ,0 ,1 / 4 ,0 ,1 / 4 ,0 ,1 / 4 ,0 ,1 / 4 ,0 ,0 ,0 ,1 / 3 ,0 ,1 / 3 ,0 ,0 ,0 ,1 / 3 ,0 ,0 ,0 ,1 / 2 ,0 ,0 ,0 ,1 / 2 ,0 ,0 ,0 ,0 ,0 ,1 / 3 ,0 ,1 / 3 ,0 ,1 / 3 ,0 ,0 ,0 ,0 ,0 ,1 / 2 ,0 ,1 / 2 ,0) ,9) # Genera la matriz de transición para la rata en un laberinto X <-1 # El vector X acumulará la trayectoria que sigue la rata ; comienza en el cuarto 1. N <-0 # Paso en el que vamos while ( tail (X ,1) ! = 9) { # Mientras la rata no encuentre la comida del cuarto 9 X <-c (X , sample ( c (1:9) ,1 , prob = P [ tail (X ,1) ,]) ) # Escogemos un cuarto al azar a partir del que se encuentra N <-N +1 # Especificamos que se ha dado un paso más } Listing 2.1. Rata.R Como un ejemplo, se obtuvieron las siguientes dos trayectorias simuladas • 123212121414141236565232363236321212369 • 1456369 A continuación presentamos una serie de preguntas para las cuales la teorı́a subsecuente encuentra una respuesta. • ¿Cuánto tarda la rata en promedio en encontra la comida si comienza en el cuarto i? • Si quitamos la comida y nada más seguimos la trayectoria de la rata, ¿Cuál es la probabilidad de que se encuentre en el cuarto j en el paso n si comienza en i? Parte de la teorı́a que veremos nos dice que la probabilidad se estabiliza conforme n → ∞. • Si de nuevo seguimos sólamente la trayectoria sin comida, ¿estamos seguros de regresar al punto inicial? • Si agregamos la comida, ¿Cuántas veces regresará la rata al cuarto inicial antes de encontrar la comida? A continuación daremos un marco teórico que permite generalizar al modelo anterior. Se trata de las cadenas de Markov cuyo estudio abarcará este capı́tulo. Sea E un conjunto a lo más numerable al que llamaremos espacio de estados. Consideremos a una colección numérica P = (Px,y )x,y∈E a la que pensaremos como una matriz indexada por E. Supongamos que (1) P Px,y ≥ 0 para todo x y y y (2) y Px,y = 1 para toda x. A P le llamamos matriz estocástica. Consideremos también una distribución de probabilidad π sobre E, que podemos pensar como un vector (digamos rengón) π = (πx )x∈E . Definición. Una cadena de Markov con matriz de transición P y distribución inicial π es un proceso estocástico (Xn )n∈N con valores en E tal que si 14 x0 , . . . , x1 ∈ E entonces P(X0 = x0 , . . . , Xn = xn ) = πx0 Px0 ,x1 · · · Pxn−1 ,xn . Teorema 2.1. Dada una matriz de transición P y una distribución inicial π existe un espacio de probabilidad en el que están definidas una sucesión de variables aleatorias (Xn )n∈N definidas en él que conforman una cadena de Markov con matriz de transición P y distribución inicial π. La demostración del teorema es importante pues nos provée de un algoritmo de simulación para cadenas de Markov. Representa otra ilustración del hecho de que cualquier proceso estocástico se puede construir mediante variables uniformes independientes. Demostración. Al enumerar a los elementos de E, podemos pensar que E = {0, . . . , n} ó E = N. Sea φi la función de cuantiles asociada al renglón i de P , φ la función de cuantiles de π y sea (Ω, F , P) un espacio de probabilidad en el que están definidas una sucesión (Ui )i∈N de variables uniformes independientes. Definimos a X0 = φ(U0 ) y Xn+1 = φXn (Un+1 ) . Por definición de función de cuantiles: P(φx (Uj ) = y) = Px,y , por lo que se sigue que si x0 , . . . , xn ∈ E entonces P(X0 = x0 , . . . , Xn = xn ) = P φ(U0 ) = x0 , φx0 (U1 ) = x1 , . . . , φxn−1 (Un ) = xn Y
= P(φ(U0 ) = x0 ) i = 1n−1 P φxi−1 (U1 ) = x1
= πx0 Px0 ,x1 · · · Pxn−1 ,xn .  Una de las caracterı́sticas principales de las cadenas de Markov es la propiedad de Markov: Proposición 2.1 (Propiedad de Markov). Sea X una cadena de Markov de distribución inicial π y transición P . Si P(X0 = x0 , . . . , Xn = xn ) > 0 entonces P(Xn+1 = xn+1 |X0 = x0 , . . . , Xn = xn ) = P(Xn+1 = xn+1 |Xn = xn ) . La interpretación es que la evolución futura de la cadena sólo depende del pasado a través del presente. Demostración. Calculemos el lado izquierdo: por lo cual P(X0 = x0 , . . . , Xn+1 = xn+1 ) = πx0 Px0 ,x1 · · · Pxn ,xn+1 P(Xn+1 = xn+1 |X0 = x0 , . . . , Xn = xn ) = Pxn ,xn+1 . 15 Por otra parte, puesto que {Xn = xn } = [ x0 ,...,xn−1 {X0 = x0 , . . . , Xn−1 = xn−1 , Xn = xn } donde la unión es disjunta, se sigue que X P(Xn = xn ) = πx0 Px0 ,x1 · · · Pxn−1 ,xn x0 ,...,xn y que P(Xn = xn , Xn+1 = xn+1 ) = X x0 ,...,xn πx0 Px0 ,x1 · · · Pxn−1 ,xn Pxn ,xn+1 por lo que P(Xn+1 = xn+1 |Xn = xn ) = Pxn ,xn+1 .  Esta misma técnica de descomposición del espacio de estados nos lleva a lo que se conoce como las ecuaciones de Chapman-Kolmogorov. Proposición 2.2. Sea X una cadena de Markov de distribución inicial π y n transición P . Si Px,y = P(Xn = y |X0 = x) entonces X m n n+m Px,y Py,z . Px,z = y∈E La ecuación anterior recuerda mucho a la de multiplicación de matrices. En n efecto, lo que nos dice es que Px,y es la enésima potencia de la matriz de transición P. Demostración. Al generalizar la idea de la prueba anterior, vemos que si definimos x0 = x y xn+m = z entonces n+m Px,z P(Xn+m = z |X0 = x) X Px0 ,x1 · · · Pxn+m−1 ,xn+m = x1 ,...,xn+m−1 ∈E = X X y∈E x1 ,...,xm−1 ∈E = X m n Px,y Py,z . Px0 ,x1 · · · Pxm−1 ,y X xm+1 ,...,xn+m−1 ∈E Py,xm+1 · · · Pxn+m−1 ,z y∈E  0 -100 -200 X[1:1000] -300 -200 X[1:1000] -300 -30 0 20 40 60 80 100 -400 -500 -500 -50 -400 -40 X[1:100] -20 -100 -10 0 0 16 0 200 Index 400 600 800 1000 Index 0 200 400 600 800 1000 Index Figura 2. 100, 1000 y 10000 pasos de una caminata aleatoria simple con p = 1/4 Por supuesto, en general no es posible calcular explı́citamente las potencias de la matriz de transición. Sin embargo, un paquete como R es capáz de realizar este producto de manera numérica y ası́ poder resolver problemas de orden práctico que se modelen mediante cadenas de Markov. A continuación, se presentan algunos ejemplos de cadenas de Markov. Ejemplo 2.1. La caminata aleatoria simple es una cadena de Markov cuyo espacio de estados es Z y es tal que Pi,i+1 = 1 − Pi,i−1 = p para alguna p ∈ (0, 1). Este es uno de los ejemplos introductorios. Basta entonces mencionar que se puede simular una trayectoria de longitud fija n (de hecho 2 trayectorias) mediante el siguiente código. # Código para simular una trayectoria de n pasos de una caminata aleatoria simple de parámetro p p <-1 / 2 n <- 10000 U <- runif ( n ) Y <-2 * (U < p ) -1 X <- cumsum ( Y ) plot ( X [1:100] , type = " l " ) quartz () # usar x11 () en UNIX y windows () en Windows , esto es para mac . plot ( X [1:1000] , type = " l " ) quartz () plot ( X [1:10000] , type = " l " ) Listing 2.2. CAS1.R Se pueden obtener entonces trayectorias como las de las Figuras 2 y 3 en las que se examinan trayectorias de 100, 1000 y 10000 pasos respectivamente para los parámetros 1/4 y 1/2. En la primera se aprecia la ley fuerte de los grandes números. Por otra parte, se pueden calcular numéricamente las probabilidades de transición a n pasos mediante el código: pa <- .5 n <-6 # Probabilidad de ir de i a i +1 # Cantidad de pasos que daremos -50 -150 -40 -200 -25 -50 -20 -100 X[1:10000] -20 X[1:1000] -30 -10 -15 X[1:100] -5 -10 0 0 0 17 0 20 40 60 80 100 0 200 Index 400 600 800 1000 0 2000 Index 4000 6000 8000 10000 Index Figura 3. 100, 1000 y 10000 pasos de una caminata aleatoria simple con p = 1/2 p <- matrix (0 , n +1 ,2 * n +1) # La entrada P [i , j ] nos da la probabilidad de encontrarnos en j - i al paso i -1 p [1 ,1] <-1 # Inicia lización : comenzamos en 0 al paso 0 con probabilidad 1 for ( i in 1: n ) { # Con cada paso actualizamos nuestras probab ilidade s p [ i +1 ,]=(1 - pa ) * p [i ,]+ pa * c (0 ,0 , p [i ,1:(2 *n -1) ]) } Listing 2.3. CASnPasos.R Podemos entonces obtener la matriz P tal que Pi,j nos da la probabilidad de que una caminata aleatoria simple esté en el estado j − i al paso i − 1. Para que cupiera en la página sólo se corrió con n = 6, pero computacionalmente n = 1000 no representa ningún problema. Una vez almacenados estos datos se pueden utilizar para obtener numericamente la media, varianza, o la esperanza de alguna otra función de la variable aleatoria que nos mide la posición después de n pasos. 1 2 3 4 5 6 7 1 1.00 0.50 0.25 0.12 0.06 0.03 0.02 2 0.00 0.00 0.00 0.00 0.00 0.00 0.00 3 0.00 0.50 0.50 0.38 0.25 0.16 0.09 4 0.00 0.00 0.00 0.00 0.00 0.00 0.00 5 0.00 0.00 0.25 0.38 0.38 0.31 0.23 6 0.00 0.00 0.00 0.00 0.00 0.00 0.00 7 0.00 0.00 0.00 0.12 0.25 0.31 0.31 8 0.00 0.00 0.00 0.00 0.00 0.00 0.00 9 0.00 0.00 0.00 0.00 0.06 0.16 0.23 10 0.00 0.00 0.00 0.00 0.00 0.00 0.00 11 0.00 0.00 0.00 0.00 0.00 0.03 0.09 12 0.00 0.00 0.00 0.00 0.00 0.00 0.00 13 0.00 0.00 0.00 0.00 0.00 0.00 0.02 Ejemplo 2.2 (Cadena de nacimiento y muerte). Se trata de una cadena de Markov cuyo espacio de estados es E = {0, . . . , n} ó N = {0, 1, 2, . . .} con probabilidades de transición son Pi,i+1 = p(i) y Pi,i−1 = q(i) donde 1 − q(i) = p(i) ∈ [0, 1]. (Definimos q(0) = 0 y si E = {0, . . . , n} entonces p(n) = 0.) Ejemplo 2.3 (Cadena de Ehrenfest). En este ejemplo hay dos urnas, con bolas numeradas del 1 al n repartidas entre ambas. A cada instante de tiempo se escoge 680 700 720 X 740 760 780 18 0 2000 4000 6000 8000 10000 Index Figura 4. Trayectoria simulada de 1000 pasos de la cadena de Ehrenfest con n = 1000 un número al azar entre 1 y n y la bola con ese número se cambia de urna. Lo que se mide es la cantidad de bolas en la urna 1 (digamos). Esta será una cadena de Markov con espacio de estados E = {0, . . . , n} y matriz de transición P dada por P0,1 = 1 = Pn,n−1 , Pi,i+1 = 1 − i/n si i < n yPi,i−1 = i/n si i > 0. Este es un caso particular de la cadena de nacimiento y muerte con espacio de estados finito. Se puede simular la cadena mediante un código como el siguiente: # Código para simular una trayectoria de espacio de estados {0 ,... , n }. n <- 1000 m <- 1000 U <- runif ( m ) X <-n / 2 for ( i in 2: m ) { aux <- tail (X ,1) if ( aux ==0) { X <-c (X ,1) } else if ( aux == n ) { X <-c (X ,n -1) } else { X <-c (X , aux +1 -2 * ( U [ i ] < aux / n ) ) } } plot (X , type = " l " ) m pasos de una cadena de Ehrenfest con Listing 2.4. Ehrenfest.R Con él, se obtuvo la Figura 4. Como ejemplo final, el lector puede verificar la liga http://www.r-bloggers. com/basics-on-markov-chain-for-parents/ a un blog en el que se interpreta al juego de serpientes y escaleras en términos de cadenas de Markov con código en R para simular el desarrollo del juego. 1. Clases de comunicación 19 1. Clases de comunicación Sean P una matriz de transición sobre E y X una cadena de Markov con matriz de transición P y distribución inicial ν tal que νx > 0 para toda x ∈ E. Denotaremos por Px a P condicionada por X0 = x. Puesto que X Px (X1 = x1 , . . . , Xn−1 = xn−1 , Xn = y) , Px (Xn = y) = i1 ,...,in−1 ∈E vemos que Px (Xn = y) = X i1 ,...,in Px,x1 · · · Pxn−1 ,y . Por lo tanto, si se introducen a las potencias de la matriz de transición P n , n ≥ 1 0 (y se define Px,y = δx,y ) vemos que n Px (Xn = y) = Px,y . Sean x y y dos estados de E. Diremos que x conduce a y si existe n ≥ 0 tal n que Px,y > 0. Claramente esto ocurre si y sólo si existen x0 , . . . , xn con x0 = x y xn = y tales que Pxk−1 ,xk > 0. Cuando x conduce a y y y conduce a x, diremos que x y y se comunican y lo denotaremos mediante x ∼ y. Proposición 2.3. La relación x ∼ y es una relación de equivalencia en E. A las clases de equivalencia inducidas por la relación ∼ les llamaremos clases de comunicación. Demostración. 0 = 1, vemos que x ∼ x. Reflexividad: Puesto que Px,x Simetrı́a: Por definición x ∼ y si y sólo si x ∼ y. m Transitividad: Si x ∼ y y y ∼ z, sean m y n en N tales que Px,y >0y n Py,z > 0. Puesto que n+m n m Px,z ≥ Px,y Py,z > 0, vemos que x conduce a z y un argumento análogo muestra que entonces x ∼ z.  Se dice que una cadena de Markov es irreducible si tiene una sola clase de comunicación. A la clase de comunicación a la que pertenece el estado x ∈ E la denotamos por Cx ; explı́citamente: Cx = {y ∈ E : x ∼ y} . El concepto de clase de comunicación nos permite dar una primera descomposición del espacio de estados. Ésta se puede refinar al introducir el concepto 2. La propiedad de Markov fuerte 20 de clase de comunicación abierta y cerrada. Este concepto es útil pues se puede reducir el espacio de estados de una cadena de Markov a una clase de comunicación cerrada. Definición. Sea C un subconjunto del espacio de estados E. Decimos que C es un conjunto cerrado si para toda y ∈ E \ C, x no conduce a y. Un conjunto abierto es aquel que no es cerrado. 2. La propiedad de Markov fuerte La propiedad de Markov fuerte es una extensión de la propiedad de Markov a ciertos tiempos aleatorios. Es una herramienta de gran utilidad. En particular nos servirá para estudiar los conceptos de transitoriedad y recurrencia. Antes de pasar a la propiedad de Markov fuerte, veamos la siguiente extensión de la propiedad de Markov. Proposición 2.4. Sea A cualquier subconjunto de E n tal que P(A ∩ {Xn = y}) > 0. Entonces, condicionalmente a A ∩ {Xn = y}, el proceso (Xn+m , m ≥ 0) es una cadena de Markov que comienza en y y tiene matriz de transición P . Demostración. Al descomponer al conjunto A como unión de eventos elementales de la forma {(x0 , . . . , xn−1 )}, vemos que P(A, Xn = y, Xn+1 = y1 , . . . , Xn+m = ym ) X = P(X0 = x0 , . . . , Xn−1 = xn−1 , Xn = y, Xn+1 = y1 , . . . , Xn+m = ym ) (x0 ,...,xn−1 )∈A = X (x0 ,...,xn−1 )∈A Px1 ,x1 · · · Pxn−1 ,y Py,y1 · · · Pym−1 ,ym Ası́, se obtiene P(Xn+1 = y1 , . . . , Xn+m = ym |A, Xn = y) = Pxn−1 ,y Py,y1 · · · Pym−1 ,ym .  Ahora verificaremos que la propiedad de Markov se extiende a ciertos tiempos aleatorios. Un tiempo aleatorio es una variable aleatoria T : Ω → N ∪ {∞}. Dicho tiempo aleatorio es finito si T (ω) ∈ N para toda ω ∈ Ω y es acotado si existe K ∈ N tal que T (ω) ≤ K para todo ω ∈ Ω. Definición. Un tiempo aleatorio T es un tiempo de paro si para toda n ∈ N existe An ⊂ E n+1 tal que {T = n} = {(X0 , . . . , Xn ) ∈ An } . Intuitivamente un tiempo de paro es un tiempo que obtenemos de observar la trayectoria hasta que se cumpla una condición. El instante en que se cumple es el tiempo de paro. 2. La propiedad de Markov fuerte 21 Nuestro primer ejemplo de un tiempo de paro es el tiempo T1 en que una cadena X regresa a su estado inicial. En otras palabras: ( ∞ Xn 6= X0 para toda n T1 = . min {n ≥ 1 : Xn = X0 } en otro caso En efecto es un tiempo de paro puesto que y para n ≥ 2 {T1 = 1} = {X1 = X0 } {T1 = n} = {X1 6= X0 , . . . Xn−1 6= X0 , Xn = X0 } . De igual manera, el tiempo Tn en que ocurre la enésima visita al estado inicial es un tiempo de paro. Esto se prueba por inducción al notar que ya hemos verificado la base inductiva n = 1 y por otro lado [ {Tn+1 = m} = {Tn = l} ∩ {Xl+1 6= X0 , . . . , Xm−1 6= X0 , Xm = X0 } . l<m Otro ejemplo de un tiempo de paro es la primera vez HA en que la cadena accede a un subconjunto A del espacio de estados. En otras palabras: ( ∞ Xn ∈ E \ A para toda n HA = . min {n ≥ 0 : Xn ∈ A} en otro caso Teorema 2.2 (Propiedad de Markov fuerte). Sea A cualquier subconjunto de E n+1 tal que P(A, Xn = y, T = n) > 0. Entonces, condicionalmente a A ∩ {T = n, Xn = y}, el proceso (Xn+m , m ≥ 0) es una cadena de Markov que comienza en y y tiene matriz de transición P . Demostración. Sea An ⊂ E n+1 tal que {T = n} = {(X0 , . . . , Xn ) ∈ An } . Al descomponer al conjunto A como unión de eventos elementales de la forma {(x0 , . . . , xn−1 )}, vemos que P(A, T = n, Xn = y, Xn+1 = y1 , . . . , Xn+m = ym ) X = P(X0 = x0 , . . . , Xn−1 = xn−1 , Xn = y, Xn+1 = y1 , . . . , Xn+m = ym ) (x0 ,...,xn )∈A∩An = X (x0 ,...,xn )∈A∩An Px1 ,x1 · · · Pxn−1 ,y Py,y1 · · · Pym−1 ,ym Ası́, se obtiene P(Xn+1 = y1 , . . . , Xn+m = ym |A, T = n, Xn = y) = Pxn−1 ,y Py,y1 · · · Pym−1 ,ym .  3. Transitoriedad y recurrencia 22 3. Transitoriedad y recurrencia Pasaremos ahora al análisis de dos conceptos que permiten hacer una distinción entre los estados de una cadena de Markov, de acuerdo a si siempre serán revisitados o no. Sea x ∈ E. Definamos a la cantidad de visitas al estado x como la variable aleatoria ∞ X 1Xn =x . Vx = n=0 Esta variable aleatoria podrı́a tomar el valor infinito. Sin embargo, un resultado curioso es que si toma el valor infinito con probabilidad positiva, entonces toma el valor infinito con probabilidad 1. En caso de que Vx sea infinita con probabilidad 1 bajo Px hablamos de un estado recurrente y en caso contrario de un estado transitorio. Analicemos ahora por qué el conjunto {Vx = ∞} tiene probabilidad cero o uno. Para esto, definamos a T0 , T1 , . . . como los instantes sucesivos que X visita al estado x. Bajo la medida Px , T0 = 0, T1 = min {n > 0 : Xn = x} y Tn+1 es la primera vez que la cadena de Markov (XTn +m , m ≥ 0) regresa a x. Hemos visto que cada Tn es un tiempo de paro. Notemos que ∞ X 1Tn <∞ . Vx = n=1 Se afirma ahora que bajo Px , Vx es una variable aleatoria geométrica de parámetro Px (T1 < ∞). En efecto, por una parte se tiene que {Vx ≥ n} = {Tn < ∞} y por otra,la propiedad de Markov fuerte nos permite afirmar que para cada n ≥ 1: Px (Tn+1 < ∞) = Px (Tn+1 < ∞, Tn < ∞) = Ex (1Tn <∞ Px (T1 < ∞)) , por lo cual n Px (Tn < ∞) = Px (T1 < ∞) . El caso en que Px (T1 < ∞) = 1 ocurre si y sólo si Vx es infinita Px casi seguramente. Si no, Vx es geométrica de parámetro Px (T1 < ∞) y por lo tanto su esperanza es finita. Esto nos proporciona una equivalencia, en términos de la matriz de transición, para que un estado sea recurrente. P n Proposición 2.5. El estado x es recurrente si y sólo si x Px,x = ∞. Demostración. La afirmación se sigue de notar que X X n Ex (Vx ) = Ex (1Xn =x ) = Px,x . n n  3. Transitoriedad y recurrencia 23 Ahora veremos que la transitoriedad o recurrencia es de hecho una propiedad de clase. Proposición 2.6. Si x y y se comunican entre si e x es transitorio entonces y es transitorio. m n Demostración. Sean m y n tales que Px,y > 0 y Px,y > 0. Entonces l n m+l+n m Px,x ≥ Px,y Py,y Py,x . Por lo tanto: si X n X m+l+n Px,x < ∞ entonces n l Py,y < ∞.  La conclusión que obtenemos es que en una clase o todos los estados son recurrentes o todos son transitorios y que por lo tanto podemos hablar de clases recurrentes y de clases transitorias. Hay una forma fácil de saber si una clase es transitoria. Proposición 2.7. Sea C ⊂ E una clase abierta. Entonces C es transitoria. Demostración. En efecto, puesto que C es una clase abierta, existe x ∈ C, m n y ∈ E \ C y m ≥ 0 tal que Px,y > 0 mientras que Py,x = 0 para toda n ≥ 0. Por lo tanto ∞ ∞ X X n Py,x =0 Ey (1Xn =x ) = Ey (Vx ) = n=0 n=0 y puesto que Vx es una variable aleatoria no-negativa, entonces Py (Vx = 0) = 1. Ası́, vemos que m Px (Vx < ∞) ≥ Px (Vx (Xm , Xm+1 , . . .) = 0, Xm = y) = Px,y >0 por lo que x es transitorio.  Veremos ahora que la conclusiónes anteriores nos permiten clasificar a las clases de cadenas de Markov con espacio de estados finito. En efecto, Proposición 2.8. Si el espacio de estados es finito, una clase es recurrente si y sólo si es cerrada. Demostración. Sólo hace falta verificar que si C es cerrada entonces es recurrente. Puesto que C es cerrada, vemos que para cualquier x ∈ C, 1 = Px (Xn ∈ C para toda n ≥ 0) . Por otra parte, al ser E finito, lo anterior forza a que exista y ∈ C que se visita infinitas veces bajo Px : 0 < Px (Vy = ∞) . 3. Transitoriedad y recurrencia 24 Si T denota a la primera visita a y, vemos que 0 < Px (Ty < ∞) Py (Vy = ∞) de acuerdo a la propiedad de Markov fuerte. Por lo tanto, vemos que y es recurrente y que ası́ la clase C es recurrente.  La primera conclusión es que una cadena irreducible con espacio de estados finitos tiene a todos los estados recurrentes. Un ejemplo muy concreto serı́a el de la cadena de Ehrenfest. Aún más, en una cadena irreducible y recurrente, de cualquier estado se accede a cualquier otro. Proposición 2.9. Si la cadena es irreducible entonces Px (Vy = ∞) = 1 para toda x, y ∈ E En particular, si recordamos la cadena de Markov que modela el movimiento de una rata por el laberinto ilustrado en la Figura 1, si colocamos comida en algúna celda, entonces la rata la encontrará con probabilidad 1. Demostración. Recordemos que Py (Vy = ∞) = 1 para toda y ∈ E. Por otra parte, al aplicar la propiedad de Markov al instante n, vemos que X n Py,x Px (Vy = ∞) . 1 = Py (Vy = ∞) = x∈E Del lado derecho tenemos un promedio ponderado de las cantidades Px (Vy = ∞) ≤ 1. El promedio es igual a 1 si y sólo sı́ Px (Vy = ∞) = 1 para toda y tal que n n Py,x > 0. Por irreducibilidad, para toda y existe n tal que Py,x > 0 y por lo tanto Px (Vy = ∞) = 1 para toda x, y ∈ E.  Para la cadena de la ruina del jugador, donde el espacio de estados es {0, . . . , N } y la matriz de transición es   i < N, j = i + 1 p Pi,j = 1 − p i > 0, j = i − 1 ,   1 i = 0, N vemos que 0 y N son absorbentes y que de cualquier estado se accede a 0 y a N . Hay por lo tanto 3 clases de comunicación: {0} , {1, . . . , N − 1} , {N }. La primera y la última son cerradas y por lo tanto recurrentes mientras que la segunda es abierta y por lo tanto transitoria. Cabe la pregunta de si en esta cadena alcanzamos alguno de los estados 0 y N con probabilidad 1. La respuesta es por supuesto afirmativa y se puede generalizar como sigue: Proposición 2.10. Sean A = {x ∈ E : Cx es abierta} , C = {x ∈ E : Cx es cerrada} 3. Transitoriedad y recurrencia 25 y HC = ( ∞ min {n ≥ 0 : Xn ∈ C} si {n ∈ N : Xn ∈ C} = ∅ . en caso contrario Si A es finita entonces para todo x ∈ E, Px (HC < ∞) = 1. Demostración. Si x ∈ C entonces 1 = Px (HC = 0) ≤ Px (HC < ∞). Por la propiedad de Markov fuerte, si x ∈ E y y ∈ A, se tiene que Px (Vy = ∞) = Px (Ty < ∞) Py (Vy = ∞) . El segundo factor del lado derecho  que y es transitorio. Puesto que P es cero puesto V = ∞ = 0 para todo x, y ∈ A. E es finito, se concluye que Px y∈A y P Por otra parte, ya que y∈E Vy = ∞, se deduce que para x ∈ A:   X Vy = ∞ = 1. Px  y∈C Como  Px (HC < ∞) ≥ Px  X y∈C se deduce el resultado deseado.  V y = ∞ ,  El análisis de recurrencia y transitoriedad de cualquier cadena con espacio de estados finito es como el de los dos ejemplos anteriores. Sin embargo, cuando el espacio de estados es infinito, el análisis de recurrencia y transitoriedad es mucho más delicado. Se presenta un ejemplo célebre de este tipo de análisis. Para más ejemplos, necesitaremos profundizar en la relación entre la propiedad de Markov y las relaciones de recurrencia, en la siguiente sección. Ejemplo 2.4 (Recurrencia y transitoriedad de la caminata aleatoria simple y el teorema de Pólya). Sea P dada por Pi,i+1 = p = 1 − Pi,i−1 para i ∈ Z, donde p ∈ (0, 1) y sea S = (Sn ) una cadena de Markov con esta matriz de transición. Esta cadena es irreducible y por lo tanto, basta ver la recurrencia o transitoriedad del cero y los demás estados compartirán esta caracterı́stica. Si p 6= 1/2, la ley fuerte de los grandes números nos dice que Sn /n → 2p − 1 conforme n → ∞ casi seguramente. Por lo tanto, para ε < p ∧ 1 − p, existe N tal que para toda n ≥ N , se tiene que (2p − 1 − ε) n < Sn < (2p − 1 + ε) n. Vemos que entonces el número de visitas al estado cero es finito casi seguramente y por lo tanto 0 es transitorio. Note que este argumento no es válido cuando
p = 1/2. −2n 2n . Ahora Si p = 1/2, entonces se puede calcular explı́citamente P0,0 = 2n n 2 utilizamos la fórmula de Stirling, la cual afirma que n! √ = 1, lim n→∞ nn e−n 2πn 4. Ejemplos de utilización de la propiedad de Markov vemos que Por lo tanto 26 n P0,0 √ . 1/ πn X n n P0,0 =∞ y ası́ vemos que S es una cadena recurrente. Consideremos ahora una caminata aleatoria simple y simétrica enP dimensión d, que tomaremos igual a 2 ó 3. Si i, j ∈ Zd , definiremos Pi,j = 1/2d si l |il − jl | = 1. Mediante un argumento combinatorio y la fórmula de Stirling, se puede probar n /nd/2 → Cd donde Cd ∈ (0, ∞). Entonces, vemos que si d = 2 la caminata que P0,0 aleatoria simple y simétrica es recurrente mientras que si d = 3, es transitoria. Éste último resultado se conoce como Teorema de Pólya para la caminata aleatoria simple. 4. Ejemplos de utilización de la propiedad de Markov Ejemplo 2.5 (El problema de la ruina). Consideremos la matriz de transición Pi,i+1 = 1 − Pi,i−1 = p ∈ (0, 1) que da lugar a una caminata aleatoria simple en Z. Sea X una cadena de Markov con transición P y tal que comienza en i bajo Pi . Consideremos N ≥ 2 y m ∈ {1, . . . , N − 1}. Definamos al siguiente tiempo de paro: T = min {n ≥ 0 : Xn ∈ {0, N }} . El objetivo del ejemplo es utilizar la propiedad de Markov, mediante el llamado análisis del primer paso para calcular la probabilidad siguiente: (2) Pm (XT = N ) . La interpretación es la siguiente: un jugador tiene un capital inicial de m pesos y se enrola en un juego de volados en el que gana un peso con probabilidad p (y los juegos son independientes). No deja de jugar hasta que pasan dos cosas: acumula un capital objetivo de N pesos o pierde toda su fortuna. Lo que se quiere determinar es la probabilidad de que termine de jugar con N pesos. Claramente, para este problema, es igual trabajar con la caminata aleatoria simple, y su espacio de estados infinito, que con la caminata aleatoria con absorción en 0 y en N que tiene espacio de estados {0, . . . , n} y matriz de transición Q dada por Q0,0 = QN,N = 1 , Qi,i+1 = p = 1 − Qi,i−1 si 1 ≤ i ≤ N − 1. La cantidad Ql m, 0 nos da la probabilidad haber llegado a cero antes del instante l y de visitar el estado N . El siguiente código R permite hacer un cálculo de esta cantidad cuando N = 10, l = 26, m = 5 y p = 5/11. # Objetivo : encontrar , numéricamente , las p robabil idades de transicion y de ruina para la cadena de la ruina del jugador 4. Ejemplos de utilización de la propiedad de Markov 27 N <- 10 # Capital objetivo p <- 5 / 11 # Probabilidad de ir hacia arriba m <- floor ( N / 2) # Capital inicial P <- matrix (0 , N +1 , N +1) # Matriz de transicion P [1 ,1] <-1 P [ N +1 , N +1] <-1 for ( i in 2: N ) { P [i ,( i -1) :( i +1) ] <-c (1 -p ,0 , p ) } p <- matrix (0 ,1 , N +1) # Distribucion inicial p [ m +1] <-1 aux <-p l <- 26 # Cantidad de pasos que se haran M <- matrix (0 , l +1 , N +1) M [1 ,] <- aux for ( i in 1: l ) { aux <- aux % * % P M [ i +1 ,] <- aux } library ( xtable ) # Cargar el paquete que genera la tabla print ( xtable ( M ) , type = " latex " , file = " Ruina2 . tex " ) # Genera el TeX de la tabla Listing 2.5. Ruina2.R Se escogió n pequeño para poder presentar los resultados en la Tabla 1. Sólo se imprimen dos decimales. Antes que nada, probaremos que T < ∞ Pm casi seguramente. Al notar que si hay N incrementos de X igual a 1 (lo cual sucede con probabilidad pN > 0, entonces T < ∞). Sin embargo, como bajo Pm los incrementos son independientes y toman el valor 1 con probabilidad p, esto sucedera casi seguramente (como se ve por ejemplo al utilizar Borel-Cantelli). Escribamos qm = Pm (XT = N ) y determinemos ahora el valor de qm . Hay dos casos sencillos: q0 = 0 y qn = 1. Por otro lado, observamos que si m ≤ n − 2 entonces 1 < T por lo que qm = pqm+1 + (1 − p) qm−1 . Ası́, la probabilidad de interés queda determinada por una relación de recurrencia con valores de frontera. Afortunadamente, la solución a la relación de recurrencia se conoce. En efecto, escribamos la relación de recurrencia en la forma pqm + qqm = pqm+1 + qqm−1 ó pqm + qqm = pqm+1 + qqm−1 4. Ejemplos de utilización de la propiedad de Markov 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 1 2 3 4 5 6 7 8 9 10 11 0.00 0.00 0.00 0.00 0.00 1.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.55 0.00 0.45 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.30 0.00 0.50 0.00 0.21 0.00 0.00 0.00 0.00 0.00 0.16 0.00 0.41 0.00 0.34 0.00 0.09 0.00 0.00 0.00 0.09 0.00 0.30 0.00 0.37 0.00 0.20 0.00 0.04 0.00 0.05 0.00 0.20 0.00 0.34 0.00 0.28 0.00 0.12 0.00 0.02 0.05 0.11 0.00 0.27 0.00 0.30 0.00 0.19 0.00 0.05 0.02 0.11 0.00 0.20 0.00 0.29 0.00 0.24 0.00 0.12 0.00 0.04 0.11 0.11 0.00 0.25 0.00 0.26 0.00 0.17 0.00 0.05 0.04 0.17 0.00 0.19 0.00 0.26 0.00 0.21 0.00 0.11 0.00 0.07 0.17 0.10 0.00 0.22 0.00 0.23 0.00 0.16 0.00 0.05 0.07 0.22 0.00 0.17 0.00 0.23 0.00 0.19 0.00 0.10 0.00 0.09 0.22 0.09 0.00 0.20 0.00 0.21 0.00 0.14 0.00 0.04 0.09 0.27 0.00 0.15 0.00 0.21 0.00 0.17 0.00 0.09 0.00 0.11 0.27 0.08 0.00 0.18 0.00 0.19 0.00 0.13 0.00 0.04 0.11 0.32 0.00 0.14 0.00 0.18 0.00 0.15 0.00 0.08 0.00 0.13 0.32 0.07 0.00 0.16 0.00 0.17 0.00 0.11 0.00 0.04 0.13 0.36 0.00 0.12 0.00 0.17 0.00 0.14 0.00 0.07 0.00 0.14 0.36 0.07 0.00 0.15 0.00 0.15 0.00 0.10 0.00 0.03 0.14 0.40 0.00 0.11 0.00 0.15 0.00 0.12 0.00 0.06 0.00 0.16 0.40 0.06 0.00 0.13 0.00 0.13 0.00 0.09 0.00 0.03 0.16 0.43 0.00 0.10 0.00 0.13 0.00 0.11 0.00 0.06 0.00 0.17 0.43 0.05 0.00 0.12 0.00 0.12 0.00 0.08 0.00 0.03 0.17 0.46 0.00 0.09 0.00 0.12 0.00 0.10 0.00 0.05 0.00 0.18 0.46 0.05 0.00 0.11 0.00 0.11 0.00 0.07 0.00 0.02 0.18 0.48 0.00 0.08 0.00 0.11 0.00 0.09 0.00 0.05 0.00 0.19 0.48 0.04 0.00 0.09 0.00 0.10 0.00 0.07 0.00 0.02 0.19 Tabla 1. El renglón i representa la distribución en el i-ésimo paso de la cadena de la ruina del jugador ó inclusive (qm − qm+1 ) = (q/p) (qm−1 − qm ) . Se sigue entonces que qm − qm+1 = −q1 (q/p) y por lo tanto qm = ( m q1 1−(q/p) 1−(q/p) q1 m si q 6= p m si q = p = 1/2 . 4. Ejemplos de utilización de la propiedad de Markov 29 Al utilizar la igualdad qn = 1 obtenemos finalmente ( 1−(q/p)m si q 6= p N . qm = 1−(q/p) m/N si q = p Cuando p = q = 1/2, podemos adicionalmente calcular vm = Em (T ). En efecto, por la propiedad de Markov y el hecho de que T ◦ θ1 = 1T ≥1 + T , vemos que v0 = vN = 0 y 2vm = 2 + vm+1 + vm−1 . La última ecuación se puede escribir, definiendo dm = vm −vm−1 como una igualdad matricial:      1 1 dm dm+1 , = −2 −2 0 1 por lo que   m    1 1 d1 dm+1 . = −2 −2 0 1 La potencia de la matriz resultante se puede calcular y nos devuelve dm+1 = d1 − 2m. Puesto que d1 = v1 y v0 = 0, vemos que vm = d1 + · · · + dm = v1 − v1 − 2 − · · · − v1 − 2 (m − 1) = mv1 − m (m − 1) . Al utilizar vN = 0, vemos que y que por lo tanto v1 = N − 1 vm = m (N − m) . Ejemplo 2.6 (Recurrencia y transitoriedad de una caminata aleatoria simple con barrera absorbente). El ejemplo anterior nos permite establecer la recurrencia o transitoriedad de una cadena de Markov con espacio de estados infinito que está ligada a la caminata aleatoria simple. Sea Pi,i+1 = p y Pi,i−1 = q si i ≥ 1 y P0,1 = 1. Sea X una cadena de Markov con transición P tal que comienza en i bajo Pi . Notemos que P se comporta como una caminata aleatoria si se encuentra en i ≥ 1, pero si llega a cero automáticamente pasa a 1. Puesto que esta cadena y la del ejemplo anterior tienen el mismo comportamiento hasta que llegan a cero o a n, vemos que si p 6= q: i Pi (Tn < T0 ) = 1 − (q/p) n. 1 − (q/p) Al tomar el lı́mite conforme n → ∞, obtenemos ( i 1 − (q/p) Pi (T0 = ∞) = 0 q<p . q≥p 4. Ejemplos de utilización de la propiedad de Markov 30 Ası́, la caminata aleatoria simple con barrera reflejante en cero es transitoria si p > 1/2 y recurrente si p ≤ 1/2. Ejemplo 2.7 (Cadenas de nacimiento y muerte). Para una sucesión de probabilidades de éxito pi ∈ (0, 1) (con qi = 1 − pi ), sea (Pi , i ≥ 0) la familia markoviana asociada a un proceso de nacimiento y muerte con matriz de transición deterninada por Pi,i+1 = pi . Haremos un análisis de primer paso para calcular hi = Pi (T0 < ∞). Notemos que h0 = 1. De nuevo, la propiedad de Markov nos permite obtener la siguiente relación de recurrencia para i ≥ 1: hi = pi hi+1 + qi hi−1 . De nueva cuenta, si escribimos di = hi − hi+1 , se tiene que pi di = qi di−1 y sabemos que d0 = 1 − h1 . Vemos que entonces di = γ i d0 donde γi = Ası́, podemos escribir pi · · · p1 . qi · · · q1 d0 (1 + γ1 + · · · + γi−1 ) = d0 + · · · + di−1 = 1 − hi . Ahora debemos hacer un análisis más preciso para determinar a la constante faltante dP 0. Si i γi = ∞, puesto que 1 − hi ∈ [0, 1], vemos que d0 = 1 − h1 = 0 y por lo tanto hP i = 1 para toda i ≥ 1. Si i γi < ∞, definamos h̃0 = 1 y h̃i = 1 − a (1 + γ0 + · · · + γi ) , por lo que h̃i ∈ [0, 1] y h̃i = pi h̃i+1 + qi h̃i−1 si i ≥ 1. Vemos que entonces   h̃i = Ei h̃Xn para toda n ≥ 1. Sin embargo:     Ei h̃Xn = Pi (T0 ≤ n) + Ei h̃Xn 1n<T0 ≥ Pi (T0 ≤ n) . Como lo anterior es válido para toda n, vemos que h̃i ≥ hi . Ası́, hi es la mı́nima solución no-negativa a la relación de recurrencia y por lo tanto P d0 = 1/(1 + i γi ). Vemos que por lo tanto hi < 1 para toda i ≥ 1. 4. Ejemplos de utilización de la propiedad de Markov 31 Ejemplo 2.8 (Tiempos de arribo para la caminata aleatoria simple unidimensional). Sea (Pk , k ∈ Z) la familia Markoviana asociada a la caminata aleatoria simple con matriz de transición Pi,i+1 = 1 − Pi,i−1 = p ∈ (0, 1). Sea T0 el primer arribo a cero dado por T0 = min {n : Xn = 0} . Nuestro objetivo será determinar a
φ(s) = E1 sT0 . Comenzando en 2, la caminata aleatoria simple debe pasar por uno para llegar a cero, y la trayectoria de 2 a 1 tiene la misma distribución que la de 1 a cero. Por lo tanto:
2 E2 sT0 = φ(s) . Por otra parte, la propiedad de Markov al instante 1 nos dice que

2 φ(s) = E1 sT0 = (1 − p)s + psE2 sT0 = (1 − p) s + psφ(s) . Ası́, puesto que φ(s) ∈ (0, 1) para s ∈ (0, 1), vemos que p 1 − 1 − 4p (1 − p) s2 . φ(s) = 2ps Esto tiene varias implicaciones. La primera es que E1 (T0 < ∞) = lim E1 s s→1 T0
1 − |1 − 2p| = = 2p ( 1 q p p < 1/2 . p ≥ 1/2 La segunda es el cálculo, para p ≤ 1/2 de la esperanza de T0 : al derivar la ecuación que satisface φ vemos que 0 = 1 − (2p − 1) φ′ (1−) Ejemplo 2.9 (El mı́nimo de caminatas aleatorias skip-free). Sea Pk , k ∈ Z la familia markoviana asociada a una caminata aleatoria cuyos saltos pertenecen a {−1, 0, 1, . . .}. Dichas caminatas se llaman sin saltos a la izquierda (skip-free to the left) aunque el punto es que el mı́nimo acumulativo tiene como imagen un intervalo de enteros. Sea −I = min Xn . n≥0 Obviamente I = 0 si P0 (X1 = −1) = 0. Supongamos por lo tanto que este no es el caso. Veamos que entonces bajo P0 , I es una variable aleatoria geométrica (aunque posiblemente degenerada en infinito). Recordemos que las variables geométricas están caracterizadas por la propiedad de pérdida de memoria, por lo que es suficiente verificar que P0 (I ≥ m + n) = P0 (I ≥ n) P0 (I ≥ m) . 4. Ejemplos de utilización de la propiedad de Markov 32 Sin embargo, notemos primero que si T−m = min {n ∈ N : Xn = −m} , entonces {I ≥ m} = {T−m < ∞}. Segundo, sobre el conjunto {T−m < ∞} se tiene que I = I ◦ θT−m . La propiedad de Markov fuerte nos permite afirmar que P0 (I ≥ m + n) = P0 (I ≥ m) P−m (I ≥ n + m) . Finalmente, puesto que P0 es la distribución de X + m bajo P−m , vemos que P−m (I ≥ n + m) = P0 (I ≥ n) . Ahora determinemos el parámetro de la distribución de I, al que denotaremos por p y escribimos q = 1 − p. Al aplicar la propiedad de Markov al tiempo 1, vemos que
q = P0 (I ≥ 1) = P0 (PX1 (I ≥ 1)) = E0 q 1+X1 . Si φ denota a la función generadora de la distribución de salto, vemos que 1 = φ(q) .
Por la desigualdad de Hölder, es fácil ver que la función φ e−λ es log-convexa. Por lo tanto la ecuación φ(s) = 1 sólo se satisface cuando s = 1 si φ′ (1−) = E0 (X1 ) ≤ 0 (lo cual dice que I es casi seguramente infinito) mientras que admite dos soluciones, digamos q̃ y 1 si E0 (X1 ) > 0. Ahora veremos que q̃ = q en este último caso. Basta mostrar que q ≤ q̃. Puesto que

Ek q̃ X1 = q̃ k E0 q̃1X = q̃ k φ(q̃) = q̃ k , la propiedad de Markov nos permite probar que
Ek q̃ Xn = q̃ k . Por lo tanto, al utilizar la propiedad de Markov n
X


q̃ = E0 q̃ 1+Xn = E0 1T−1 =k E−1 q̃ 1+Xn−k + E0 1T−1 >n q̃ 1+Xn k=0 = P0 (T1 ≤ n) + E0 1T−1 >n q̃ 1+Xn ≥ P0 (T1 ≤ n) .
Al tomar el lı́mite conforme n → ∞, vemos que q̃ ≥ q y por lo tanto, q es la solución mı́nima en [0, 1] a la ecuación φ(s) = 1. Ejemplo 2.10 (El principio de reflexión). Consideremos la familia markoviana (Pk )k∈Z asociada a la caminata aleatoria simple y sea Tm = min {n ∈ N : Xn = m} . 5. Medidas y distribuciones invariantes Definamos X̃n = ( 33 Xn 2m − Xn si n ≤ Tm . si n > Tm El principio de reflexión, que ahora demostaremos, afirma que la distribución de X̃ bajo P0 es precisamente P0 . En efecto, el principio de reflexión es consecuencia de la propiedad de Markov fuerte y el hecho de que la distribución de −X bajo P0 es también P0 . Como una aplicación, notemos que si k ≥ 0 entonces P0 (Tm ≤ n) = P0 (Tm ≤ n, Xn ≥ m) + P0 (Tm ≤ n, Xn < m) . Puesto que {Tm ≤ n} ⊂ {Xn ≥ m} y o n {Tm ≤ n, Xn ≤ m} = X̃n ≥ m , se deduce que P0 (Tm ≤ n) = P0 (Xn = m) + 2P0 (Xn > m) . 5. Medidas y distribuciones invariantes Comenzaremos con un sistema lineal importante que se obtiene de una aplicación de la propiedad de Markov fuerte. Ejemplo 2.11 (Cantidad esperada de visitas antes de la recurrencia). Consideremos de nuevo a los tiempos de primera visita Ty = min {n ∈ N : Xn = y} . Fijemos a x ∈ E y definamos νxx =1 y νyx = Ex Tx X i=0 1Xn =y ! para y 6= x. Teorema 2.3. Para una cadena irreducible y recurrente, ν x satisface νxx = 1, 0 < νyx < ∞ y X νyz Py,z = νzx . (3) z∈E Una medida invariante es un vector renglón νy , y ∈ E con entradas nonegativas tal que νP = ν o más explı́citamente, X νyz Py,z = νz . z∈E x Por lo tanto, el vector ν representa una construcción probabilı́stica de una medida invariante. 5. Medidas y distribuciones invariantes 34 Demostración. Claramente νxx = 1. Para ver que ν x satisface (3), utilizamos la propiedad de Markov de manera similar al análisis del primer paso. Supongamos primero que y 6= x. Entonces, puesto que Tx es finito Px casi seguramente   ! ∞ X X X 1Xn =y = Ex  1Xn =y  = Ex (1Xn =y 1Tx ≤n ) . Ex n<Tx n=1 n≤Tx El sumando n = 1 es fácil de calcular puesto que bajo Px , Tx ≥ 1: Ex (1X1 =y 1Tx ≥1 ) = Ex (1X1 =y ) = Px,y . Para los sumandos con n ≥ 2, aplicamos la propiedad de Markov al instante n − 1 (haciendo una especie de análisis de último paso), notando que {T ≤ n} ∈ Fn−1 : X
Ex 1Xn−1 =z 1Xn =y 1Tx ≤n Ex (1Xn =y 1Tx ≤n ) = z6=x X z6=x
Ex 1Xn−1 =z 1Tx ≤n Pz,y . Ası́, vemos que νyx X = Ex 1Xn =y n<Tx ! = Px,y + ∞ XX Ex (1Xn =y 1Tx >n ) Pz,y z6=x n=1 = νxx Px,y + X z6=x = X νzx Pz,y . Ex TX x −1 n=0 1Xn =y ! Pz,y z∈E νxx P νyx Py,x . Ahora veremos que 1 = = En efecto, basta descomponer respecto al valor de Tx y de XTx −1 , recordando que Tx es finito Px casi seguramente: 1= ∞ X Px (Tx = n) n=1 = Px,x + ∞ X X Px (Tx = n, Xn−1 = y) n=2 y6=x = Px,x + ∞ X X n=2 y6=x = Px,x + X y6=x Px (Tx ≥ n − 1, Xn−1 = y) Px,y νyx Px,y . 5. Medidas y distribuciones invariantes 35 m n Finalmente, si consideramos m tal que Px,y > 0 y n tal que Py,x > 0 entonces por una parte X m m νyx = νzx Pz,y ≥ Px,y >0 z y por otra 1 = νxx = X z implica que νyx < ∞. n n νzx Px,z ≥ νyx Px,y  Definición. Un estado x de una cadena de Markov es positivo recurrente si Ex (Tx ) < ∞. Denotamos por mx a dicha cantidad, a la que nos referiremos como tiempo medio de recurrencia de x. Notemos que, de acuerdo a nuestra definición de ν x , se tiene que ! X X X x 1Xn = y = Ex (Tx ) . ν y = Ex y y n<Tx Si x es positivo recurrente, podemos entonces definir al vector π x = νyx /Ex (Tx ) que satisface X X πyx Py,z = πzx . πyx ≥ 0, πyx = 1 y y
y∈E Una distribución invariante es un vector renglón que satisface las 3 condiciones anteriores. Una propiedad importante y útil de una distribución invariante es que si X0 tiene distribución π, entonces Xn tendrá distribución π para toda n. Por otra parte, en espacio de estados finito es fácil ver que existen distribuciones invariantes. Corolario 1. En una cadena irreducible con espacio de estados finito, todos los estados son positivo recurrentes. Demostración. Sea x cualquier estado de la cadena. Por el teorema anterior νyx < ∞ para toda y y como E es finito, entonces X νyx < ∞.  Ex (Tx ) = y∈E En caso de que exista un único vector de probabilidad invariante, se pueden calcular tiempos medios de recurrencia al resolver un sistema de ecuaciones. A continuación profundizaremos en esta idea. Teorema 2.4. Si ν es una medida invariante para una matriz de transición irreducible P y νx = 1 entonces ν ≥ ν x . Si además P recurrente entonces ν = ν x . Demostración. Al aplicar la invariancia de ν se obtiene X X X νy1 Py1 ,y = νx Px,y + νy = νy1 Py1 ,x = Px,y + νy1 Py1 ,y . y1 ∈E y1 6=x y1 6=x 5. Medidas y distribuciones invariantes 36 Al volverla a aplicar se obtiene X X νy2 Py2 ,y1 Py1 ,y Px,y1 Py1 ,y + νy = Px,y + y1 ,y2 6=x y1 6=x y al continuar repetidamente, vemos que X X νy = Px,y + Px,y1 Py1 ,y + · · · + y1 6=x + X y1 ,··· ,yn+1 6=x y1 ,...,yn 6=x Px,yn Pyn ,yn−1 · · · Py2 ,y1 Py1 ,y νyn+1 Pyn+1 ,yn · · · Py2 ,y1 Py1 ,x . Si y 6= x, encontramos la cota νy ≥ n X m=0 Px (Xm = y, Tx ≥ m) y el lado derecho converge conforme n → ∞ a νyx , por lo que νy ≥ νyx . Por otra parte, si la cadena es recurrente además de irreducible entonces ν x es invariante, por lo cual µ = ν − ν x es una medida invariante con µx = 0. Por irreducibilidad, para toda y ∈ E existe n ≥ 0 tal que P n y, x > 0 y entonces X n n , 0 = µx = µz Pz,x ≥ µy Px,y z por lo que µy = 0.  Teorema 2.5. Para una cadena irreducible las siguientes condiciones son equivalentes. (1) Todos los estados son positivo recurrentes (2) Algún estado es positivo recurrente (3) La cadena admite una distribución invariante. En este caso, la distribución invariante es única y asigna a x el recı́proco de su tiempo medio de recurrencia. Demostración. Primero probaremos que si la cadena admite una distribución invariante ν entonces todos los estados son positivo recurrentes. Para esto, notemos primero que νx > 0 para toda x ∈ E. En efecto, lo debe ser para alguna x y n al utilizar la irreducibilidad para encontrar n tal que Py,x > 0, vemos que X n n Pz,y νz ≥ νx Px,y > 0. νy = z∈E 5. Medidas y distribuciones invariantes 37 Si x ∈ E, entonces ν/νx es una medida invariante que asigna 1 a x, por lo cual ν/νx ≥ ν x . Puesto que ν es una distribución: X X νy 1 = < ∞. mx = Ex (Tx ) = νyx ≤ νx νx y y Ası́, todos los estados son positivo recurrentes. Al ser en particular recurrentes, hemos visto que νx = νxx = mx y como νx es una distribución entonces mx < ∞. Esto termina la demostración de las implicaciones y además nos produce la fórmula requerida para la distribución invariante.  Ejemplo 2.12 (La cadena de Ehrenfest). Sea E = {0, . . . , N } y Pi,i+1 = 1 − i/N = 1 − Pi,i−1 . Esta cadena es irreducible y con espacio de estados finito, por lo que todos los estados son positivo recurrentes y por lo tanto existe un único vector de probabilidad invariante π. Además, πx = 1/Ex (Tx ). El siguiente código nos permite obtener numéricamente estas cantidades. N <- 10 # Cantidad de bolas para la cadena de Ehrenfest P <- matrix (0 , N +1 , N +1) P [1 ,2] <-1 P [ N +1 , N ] <-1 for ( i in 2: N ) { P [i , i +1] <- 1 -(i -1) / N P [i ,i -1] <-(i -1) / N } I <- diag (1 , N +1 , N +1) Null <- function ( M ) # Calcula una base del kernel de la matriz M mediante factorizacion QR { tmp <- qr ( M ) set <- if ( tmp $ rank == 0) 1: ncol ( M ) else - (1: tmp $ rank ) qr . Q ( tmp , complete = TRUE ) [ , set , drop = FALSE ] } v <-t ( Null (I - P ) ) # El kernel de I - P tiene dimension 1 por ser la cadena irreducible y finita v <-v / sum ( v ) # Se obtiene el vector de probabilidad invariante m <-1 / v # Se obtienen los tiempos medios de recurrencia library ( xtable ) # Cargar el paquete que genera tablas en TeX print ( xtable ( v ) , type = " latex " , file = " E h r e n f e s t I n v a r i a n t D i s t . tex " ) # Genera el TeX de la tabla print ( xtable ( m ) , type = " latex " , file = " E h r e n f e s t M e a n R e c u r r e n c e . tex " ) # Genera el TeX de la tabla Listing 2.6. EhrenfestMatrix.R El resultado se puede visualizar en las tablas siguientes 1 1 0.00 2 3 4 0.01 0.04 0.12 5 0.21 6 0.25 7 8 9 0.21 0.12 0.04 10 0.01 11 0.00 6. Comportamiento a tiempos grandes 1 2 3 4 1 1024.00 102.40 22.76 8.53 5 4.88 38 6 4.06 7 8 9 4.88 8.53 22.76 10 102.40 11 1024.00 Por otra parte, para esta cadena particular, podemos calcular explı́citamente el vector de probabilidad
invariante y por lo tanto los tiempos medios de recurrencia. En efecto, sea πj = Nj 2−N y notemos que si j ∈ {1, . . . , N − 1} entonces X πi Pi,j = πj−1 Pj−1,j + πj+1 Pj+1,j i    N j+1 N −N N − j + 1 2 = 2−N + j−1 j+1 N N   1 N − j + 1 j+1 1 −N + = 2 N! (j − 1)! (N − j + 1)! N (j + 1)! (N − j − 1)! N   1 N ! 1 −N =2 + (j − 1)! (N − j − 1)! N (N − j) jN   1 N! = 2−N (j − 1)! (N − j − 1)! j (N − j)   N = 2−N j = πj .  Por lo tanto la distribución binomial de parámetros N y 1/2 es invariante para la cadena de Ehrenfest con espacio de estados (0, . . . , N ). Como ésta es irreducible, los tiempos medios de recurrencia son mi = 2N
. N i Si ahora utilizamos la fórmula de Stirling, notamos que conforme N → ∞ mN/2 ∼ √ 1 πN pero m0 = 2N . 6. Comportamiento a tiempos grandes El objetivo de esta sección es estudiar cómo se comporta una cadena de Markov recurrente cuando la dejamos correr un tiempo grande. Una aplicación de este tipo de ideas sirve para desplegar el resultado de una búsqueda en Google en orden descendiente de importancia de las páginas. Para describir cómo se define la importancia de una página, podemos pensar en la siguiente situación: se considera a la red como un conjunto de nodos (páginas web) y aristas (hipervı́nculos) que los conectan. Si Nv denota la cantidad de hipervı́nculos que hay en la página v, se 6. Comportamiento a tiempos grandes 39 puede definir la importancia de una página v como la cantidad Iv que satisface el sistema de ecuaciones: X X Iu , Iv = 1. Iv = Nu v La interpretación es que una página importante transmite su importancia a cada una de las ligas que está presente en la página y que la suma de las importancias es 1. Esta es la base del algoritmo PageRank de Google. Para interpretar lo anterior en términos de cadenas de Markov, consideremos la siguiente estrategia para navegar en esta red: al encontrarme en la página web, escojo un hipervı́nculo al azar y lo recorro. Entonces mi trayectoria por la red será una cadena de Markov con matriz de transición  1  si existe hipervı́nculo de u a v  Nu 1 Pu,v = Total de páginas no hay hipervı́nculos en u .   0 otro caso Notemos que la importancia de una página es simplemente una solución I a la ecuación I = IP . Por otra parte, mediante la estrategia de navegación postulada, podrı́amos pensar a la importancia de una página web como la fracción de tiempo asintótica que le dedico cuando recorro la web al azar. En esta sección veremos que de hecho, son la misma idea. En la página http://www.nd.edu/~networks/ resources.htm el los autores Albert, Jeong y Barabási proporcionan sus datos sobre la estructura de la red que recopilaron para un estudio de 1999. Por supuesto ha crecido mucho desde entonces, pero es un buen ejemplo de datos que se han hecho de acceso público. Un pequeño detalle es que los fundadores de Google se dieron cuenta de cuestiones numéricas que aparecen al tratar de encontrar la importancia de una página, por lo que modificaron la definición de importancia a: X Iu 1−α Iv = α + , Nu N donde α ∈ [0, 1] y N es la cantidad de páginas web que hay. La interpretación en términos de estrategia de navegación es que, con probabilidad α sigo la estrategia de navegación anterior y con probabilidad 1 − α salto a una página al azar. Ası́, se estará considerando la matriz de transición  α  si existe hipervı́nculo de u a v  Nu α P = total de1páginas no hay hipervı́nculos en u   1−α otro caso N El valor utilizado por Google es α = .85. La diferencia entre P y P α es que la segunda es irreducible y aperiódica mientras la primera no tiene por qué serlo. Como ejemplo concreto utilicemos la gráfica que, el 21 de Mayo del 2012, se encuentra en la página de Wikipedia de PageRank. Los vértices son serán etiquetados del 1 al 11 y la gráfica se presenta en la Figura 5. Con α = .85, la 6. Comportamiento a tiempos grandes 11 40 10 7 8 9 5 4 6 1 2 3 Figura 5. Ilustración sobre la definición de importancia de páginas web matriz de transición es  9.09 9.09 9.09 9.09 1.36 1.36 86.36 1.36  1.36 86.36 1.36 1.36  43.8 43.8 1.36 1.36  1.36 29.70 1.36 29.70  1.36 1.36 100∗Q =  1.36 43.8 1.36 43.8 1.36 1.36  1.36 43.8 1.36 1.36  1.36 43.8 1.36 1.36  1.36 1.36 1.36 1.36 1.36 1.36 1.36 1.36 9.09 9.09 9.09 9.09 9.09 9.09 1.36 1.36 1.36 1.36 1.36 1.36 1.36 1.36 1.36 1.36 1.36 1.36 1.36 1.36 1.36 1.36 1.36 1.36 1.36 29.70 1.36 1.36 1.36 1.36 43.8 1.36 1.36 1.36 1.36 1.36 43.8 1.36 1.36 1.36 1.36 1.36 43.8 1.36 1.36 1.36 1.36 1.36 43.8 1.36 1.36 1.36 1.36 1.36 86.3 1.36 1.36 1.36 1.36 1.36 86.3 1.36 1.36 1.36 1.36 1.36  9.09 1.36  1.36  1.36  1.36  1.36  1.36  1.36  1.36  1.36 1.36 El vector de probabilidad invariante, que nos da la importancia de cada página, es 100 ∗ π = (3.28, 38.44, 34.29, 3.91, 8.09, 3.91, 1.62, 1.62, 1.62, 1.62, 1.62) Como comentamos, otra manera de pensar a la importancia de un vértice es mediante la fracción de tiempo que paso en el. La fracción esperada de tiempo que paso en j comenzando en i después de n pasos es n 1X n Qi,j . n l=0 6. Comportamiento a tiempos grandes 41 Sin embargo, al calcular numéricamente la igual a  3.28 38.44 34.29 3.91 8.09 3.28 38.43 34.30 3.91 8.09  3.28 38.45 34.28 3.91 8.09  3.28 38.45 34.28 3.91 8.09  3.28 38.44 34.29 3.91 8.09  Q50 =  3.28 38.45 34.28 3.91 8.09 3.28 38.45 34.28 3.91 8.09  3.28 38.45 34.28 3.91 8.09  3.28 38.45 34.28 3.91 8.09  3.28 38.44 34.29 3.91 8.09 3.28 38.44 34.29 3.91 8.09 entrada 50 de la cadena, vemos que es 3.91 3.91 3.91 3.91 3.91 3.91 3.91 3.91 3.91 3.91 3.91 1.62 1.62 1.62 1.62 1.62 1.62 1.62 1.62 1.62 1.62 1.62 1.62 1.62 1.62 1.62 1.62 1.62 1.62 1.62 1.62 1.62 1.62 1.62 1.62 1.62 1.62 1.62 1.62 1.62 1.62 1.62 1.62 1.62 1.62 1.62 1.62 1.62 1.62 1.62 1.62 1.62 1.62 1.62 1.62  1.62 1.62  1.62  1.62  1.62  1.62  1.62  1.62  1.62  1.62 1.62 Vemos que desde el punto de vista numérico, Q50 es ya una matriz cuyos renglones son iguales al vector de probabilidad invariante. Esto no es ninguna coincidencia y uno de nuestros objetivos será explicarlo. Se utilizó el siguiente código para generar a la distribución invariante y a la potencia 50 de la matriz. N <- 11 A = matrix ( c ( 1 ,1 ,1 ,1 ,1 ,1 ,1 ,1 ,1 ,1 ,1 , 0 ,0 ,1 ,0 ,0 ,0 ,0 ,0 ,0 ,0 ,0 , 0 ,1 ,0 ,0 ,0 ,0 ,0 ,0 ,0 ,0 ,0 , 1 ,1 ,0 ,0 ,0 ,0 ,0 ,0 ,0 ,0 ,0 , 0 ,1 ,0 ,1 ,0 ,1 ,0 ,0 ,0 ,0 ,0 , 0 ,1 ,0 ,0 ,1 ,0 ,0 ,0 ,0 ,0 ,0 , 0 ,1 ,0 ,0 ,1 ,0 ,0 ,0 ,0 ,0 ,0 , 0 ,1 ,0 ,0 ,1 ,0 ,0 ,0 ,0 ,0 ,0 , 0 ,1 ,0 ,0 ,1 ,0 ,0 ,0 ,0 ,0 ,0 , 0 ,0 ,0 ,0 ,1 ,0 ,0 ,0 ,0 ,0 ,0 , 0 ,0 ,0 ,0 ,1 ,0 ,0 ,0 ,0 ,0 ,0) ,N ,N , byrow = TRUE ) grafica P = A / rowSums ( A ) a =.85 Pa = a * P +(1 - a ) / N # Define la matriz de adyacencia de la # Normaliza los renglones para obtener una matriz estocastica # Se produce la matriz irreducible l <- 50 aux <- Pa for ( i in 1: l ) { # Se calcula la potencia 50 de la matriz de transicion aux <- aux % * % Pa } I <- diag (1 ,N , N ) tmp <- qr (I - Pa ) # Se obtiene la factorizacion QR set <- if ( tmp $ rank == 0) 1: ncol ( M ) else - (1: tmp $ rank ) Null <- qr . Q ( tmp , complete = TRUE ) [ , set , drop = FALSE ] v <-t ( Null ) # El kernel de I - P tiene dimension 1 por ser la cadena irreducible y finita 6. Comportamiento a tiempos grandes 42 v <-v / sum ( v ) # Se obtiene el vector de probabilidad invariante library ( xtable ) # Cargar el paquete que genera tablas en TeX print ( xtable ( Pa ) , type = " latex " , file = " PageRank Matrix " ) # Genera el TeX de la tabla print ( xtable ( v ) , type = " latex " , file = " P a g e R a n k I n v a r i a n t D i s t r i b u t i o n . tex " ) print ( xtable ( aux ) , type = " latex " , file = " P a g e R a n k M a t r i x P o w e r . tex " ) Listing 2.7. PageRank.R Vale la pena comparar lo anterior con la cadena de Ehrenfest: con N = 10, a continuación se presentan las potencias 50 y 51 de la matriz de transición:   2 0 88 0 410 0 410 0 88 0 2 0 20 0 234 0 492 0 234 0 20 0   2 0 88 0 410 0 410 0 88 0 2   0 20 0 234 0 492 0 234 0 20 0   2 0 88 0 410 0 410 0 88 0 2   0 492 0 234 0 20 0 1000 ∗ P 50 =  0 20 0 234  2 0 88 0 410 0 410 0 88 0 2   0 20 0 234 0 492 0 234 0 20 0   2 0 88 0 410 0 410 0 88 0 2   0 20 0 234 0 492 0 234 0 20 0 2 0 88 0 410 0 410 0 88 0 2   0 20 0 234 0 492 0 234 0 20 0 2 0 88 0 410 0 410 0 88 0 2   0 20 0 234 0 492 0 234 0 20 0   2 0 88 0 410 0 410 0 88 0 2   0 20 0 234 0 492 0 234 0 20 0   0 410 0 410 0 88 0 2 1000 ∗ P 51 =   2 0 88 0 20 0 234 0 492 0 234 0 20 0   2 0 88 0 410 0 410 0 88 0 2   0 20 0 234 0 492 0 234 0 20 0   2 0 88 0 410 0 410 0 88 0 2 0 20 0 234 0 492 0 234 0 20 0 Casi lo mismo se repite al calcular potencias superiores. En este caso, paracerı́a haber convergencia de las potencias pares y de las potencias impartes a un lı́mite distinto. La diferencia entre la matriz de transición de PageRank y la de la cadena de Ehrenfest se encuentra en el siguiente concepto. Definición. Sea P una matriz de transición. Definimos al periodo de un  n estado x como el máximo común divisor de n ≥ 0 : Px,x > 0 . Decimos que la cadena es aperiódica si el periodo de todo estado es 1. kd Si x tiene periodo d, entonces Px,x > 0 salvo para una cantidad finita de ı́ndices kd+(d−1) kd+1 k mientras que Px,x = · · · = Px,x = 0. 6. Comportamiento a tiempos grandes 43 Es fácil verificar que el periodo es el mismo para todos los elementos de una clase de comunicación. Entonces, para verificar que una cadena irreducible es aperiódica, basta verificar que un algún estado tiene periodo 1. Cuando el espacio de estados es finito, una cadena es irreducible y aperiódica N si y sólo si existe N > 0 tal que Px,y > 0 para todo x, y ∈ E. Definición. Decimos que una cadena finita es regular si existe una potencia de la matriz de transición con todas las entradas positivas. Puesto que la matriz PageRank con α < 1 tiene todas las entradas positivas, entonces es regular. Sin embargo, la cadena de Ehrenfest tiene periodo 2, por lo que no es regular. Esta es la gran diferencia entre ambos tipos de cadenas. El resultado más importante que relaciona a las potencias de la matriz de transición con el vector de probabilidad invariante es el siguiente. Teorema 2.6 (Teorema fundamental de convergencia). Si P es la matriz de transición de una cadena de Markov irreducible, aperiódica y positivo recurrente y P n π es la única distribución invariante de P entonces y Px,y − πy → 0 conforme n → ∞. Utilizando el hecho de que si una sucesión converge también converge la sucesión de sus promedios, vemos n 1x=y + Px,y + · · · + Px,y = πy . lim n→∞ n Por lo tanto la entrada y del vector de probabilidad también se puede interpretar como la fracción esperada asintótica de veces en que la cadena se encuentra en el estado y, independientemente de dónde comenzó. Hay dos tipos de prueba muy distintos para este resultado. El primero se enfoca en espacio de estado finito y utiliza el hecho de que la matriz de transición de la cadena será regular, ası́ como un argumento de punto fijo. El segundo se basa en una idea probabilı́stica conocida como acoplamiento, en el que se prueba que para dos cadenas de Markov independientes que satisfacen las hipótesis del teorema toman el mismo valor en un tiempo aleatorio y finito. A continuación exploraremos la prueba del teorema fundamental de convergencia en el caso de espacio de estados finito. Comenzaremos con el caso en que no sólo P es regular sino que además P tiene todas las entradas positivas. Definimos a ρ = min Px,y . x,y∈E Puesto que el espacio de estados es finito, el mı́nimo anterior está definido y es positivo. Sea y un vector columna. Probaremos a continuación que las entradas de P n y se van acercando entre si conforme n crece. Para esto, sean mn el mı́nimo de las entradas de P n y y M n el máximo. Primero mostraremos que m0 ≤ m1 ≤ M1 ≤ M0 y que M1 − m1 ≤ (1 − 2ρ) (M0 − m0 ). La idea, que se puede justificar con cálculo 6. Comportamiento a tiempos grandes 44 P diferencial, es que para maximizar el promedio j Pi,j yj sujeto a la restricción de que maxj yj = M0 y minj yj = m0 formamos un vector con todas sus coordenas iguales a M0 salvo una, que será igual a m0 y que colocamos con el menor peso posible, es decir, en la coordenada j tal que Pi,j = minj Pi,j ≥ ρ. Ası́, obtendremos X Pi,j yj ≤ M0 (1 − Pi,j ) + m0 Pi,j ≤ M0 (1 − ρ) + m0 ρ. j Con un argumento similar se obtiene X Pi,j yj ≥ m0 (1 − ρ) + M0 ρ j y por lo tanto M1 − m1 ≤ (1 − 2ρ) (M0 − m0 ) como se anunció. Al aplicar lo anterior con P n−1 y en vez de y, vemos que m0 ≤ m1 ≤ · · · ≤ mn ≤ Mn ≤ · · · ≤ M1 ≤ M0 y que n (Mn − mn ) ≤ (1 − 2ρ) (M0 − m0 ) → 0. Por lo tanto, P n y converge a un vector con todas sus entradas iguales. Si aplicamos lo anterior al vector ej que tiene 1 en la entrada j y cero en cualquier otra entradas, vemos que (P n ej )i (la i-ésima entrada de P n ej ) es igual a n Pi,j que esta cantidad converge a uj para cualquier i. Notemos además que m1 > 0 y por lo tanto uj > 0 para toda j. Ası́, la matriz P n converge a una matriz cuyos renglones son todos iguales al vector cuya entrada j es uj . Finalmente, notamos que u es una distribución invariante. Es una distribución pues, como el espacio de estados es finito, entonces X X X n n uj . lim Pi,j = Pi,j = 1 = lim n→∞ j j n→∞ j Por otro lado, puesto que P n+1 = P n P , vemos que X X n+1 n uj = lim Pi,j = lim (P n P )i,j = lim Pi,k Pk,j = uk Pk,j . n→∞ n→∞ n→∞ k k Finalmente, aprovechamos la ocasión para dar otra prueba de por qué en este caso, hay a lo más una distribución invariante. Si v es una distribución invariante entonces vP n = v y por otro lado X X n v i uj = uj . vi Pi,j → vj = (vP n )j = i Por lo tanto u = v. i 7. Cadenas absorbentes 45 2000 0 1000 Frequency 3000 4000 Histogram of apariciones Figura 6. Histograma de una muestra del tiempo de aparición de la palabra 101 7. Cadenas absorbentes Consideremos el siguiente problema: se tiene un abecedario fijo (digamos para fijar ideas igual a {0, 1}) con cuyos elementos conformamos palabras (como podrı́a ser 1001). Si se extraen las letras del abecedario con igual probabilidad y con independencia entre las extracciones, ¿Cuánto se debe esperar para que aparezca la palabra escogida? Más generalmente: ¿con qué probabilidad aparece la palabra por primera vez al observar n letras? Un caso particular de este problema es cuando la palabra tiene sólo una letra (digamos p); el tiempo de aparición de la palabra es entonces geométrico de parámetro p. El siguiente código toma una muestra de tamaño 10000 de la distribución del tiempo de aparición de la palabra 101 cuando las letras 1 y 0 aparencen con probabilidades 2/3 y 1/3. Se obtiene el histograma de la muestra (ver la Figura 6) ası́ como la media empı́rica, que resultó ser de 8.26. # Obtener una muestra de tama~ n o 10000 del tiempo de aparición de la palabra 101. # 1 aparece con probabilidad 2 / 3 p <-2 / 3 # Probabilidad de que salga 1 N <- 10000 # Tama \ ~ no de muestra llevo <-0 # Cu \ ’ antas simulaciones llevo apariciones <-c () # En que tiempos aparecen las palabras while ( llevo < N ) { # Mientras no lleve N simulaciones llevo <- llevo +1 # Realizo una mas n <-0 # Inicializando a $ n $ en cero p a l a b r a A l e a t o r i a <-c () # Con un nuevo repositorio de letras aux <-0 # Y una variable que indica si ya he visto la palabra while ( aux ==0) { # Mientras no haya visto la palabra n <-n +1 # Espero una unidad de tiempo p a l a b r a A l e a t o r i a <-c ( palabraAleatoria , rbinom (1 , size =1 , prob =2 / 3) ) # Obtengo una letra nueva 7. Cadenas absorbentes 46 if (n >2) { if ( p al a b r a A l e a t o r i a [n -2]==1 & p a l a b r a A l e a t o r i a [n -1]==0 & p a l a b r a A l e a t o r i a [ n ]==1) { aux <-1 }} } # Si he visto la palabra , la variable aux nos sacar \ ’ a del ciclo apariciones <-c ( apariciones , n ) # Agrego el nuevo tiempo de aparici \ ’ on } hist ( apariciones ) mean ( apariciones ) Listing 2.8. Palabra.R Consideremos otro problema similar: para el laberinto de la Figura 1, si la rata comienza su trayectoria en el cuarto i, ¿Cuál es la cantidad de tiempo esperada para que llegue a la celda 9? Estas preguntas pueden responderse con la teorı́a de cadenas de Markov absorbentes. Ya hemos discutido por qué el segundo problema está relacionado con cadenas de Markov, en particular con una cadena de Markov en la que el estado 9 es absorbente. El primero parecerı́a no estarlo. Sin embargo, asociemos una cadena de Markov a nuestro ejemplo particular. El espacio de estados es {0, 1, 2, 3, 4} donde el estado i representa que han aparecido i letras correctas de la palabra 1001. El estado 4 será absorbente y nos señalará que ya ha aparecido por primera vez la palabra. Si p representa la probabilidad de que la enésima letra sea 1 y q = 1 − p, la matriz de transición de la cadena será   q p 0 0 0 p 0 q 0 0    P = p 0 0 q 0  q 0 0 0 p 0 0 0 0 1 Si por ejemplo tuvieramos la palabra  q p q 0  0 0  q 0  q 0 0 0 11011, la matriz de transición serı́a  0 0 0 0 p 0 0 0  p q 0 0  0 0 p 0  0 0 0 p 0 0 0 1 puesto que si tengo 2 letras correctas es por que estoy viendo la palabra 11; si sale 1 seguiré observando 11 mientras que si sale 0, observaré 110. n Ahora interpretemos una entrada particular de la matriz de transición: P0,4 , para la palabra 1001. Es la probabilidad de que la cadena se encuentre en el estado al tiempo n, que equivale a que la cadena se haya absorbido en 4 antes del tiempo n, que significa a su vez que haya aparecido la palabra a lo más en la enésima extracción. Pasemos a la abstracción común de ambos problemas. 7. Cadenas absorbentes 47 Definición. Una cadena absorbente es una cadena absorbente con espacio de estados finito, que tiene al menos un estado absorbente y que desde cualquier estado se puede acceder a un estado absorbente. Los estados de una cadena absorbente los podemos dividir en no-absorbentes (T , pues son transitorios) y absorbentes (A). Si los enumeramos del 1 al n y ponemos al final a los absorbentes, la matriz de transición tomará la forma   Q R P = , 0 I donde Q es una matriz de tamaño m×m (m < n) en donde se encuentran las probabilidades de transición entre estados no-absorbentes, R es una matriz de tamaño m × n correspondiente a las probabilidades de transición de estados no-absorbentes a absorbentes, mientras que I es la matriz identidad de tamaño (n − m) × (n − m). n representa la probabilidad Si i es no-absorbente y j es absorbente, la entrada Pi,j de que comenzando en iPla cadena se haya absorbido por el estado j antes del instante n, mientras que j≤m Qni,j representa la probabilidad de que comenzando en i, la cadena todavı́a no se haya absorbido en el instante n. Es usual escribir X Qni,j = Qn 1 j≤m donde 1 es un vector columna de longitud m y entradas iguales a 1. Proposición 2.11. Sea ( ∞ T = min {n ≥ 0 : Xn ∈ A} {n ≥ 0 : Xn ∈ A} = ∅ en otro caso el tiempo de absorción de una cadena absorbente. Entonces Pi (T < ∞) = 1 para toda i ∈ E y limn→∞ Qn = 0. En el enunciado anterior, 0 es una matriz de m × m cuyas entradas son iguales a cero. Demostración. En la Proposición 2.10 hemos probado que Pi (T < ∞) = 1 puesto que T es una clase abierta y finita. Por otra parte, si i, j ≤ m: 0 = Pi (T = ∞) = lim Pi (T ≥ n) ≥ lim Pi (, Xn = j, T ≥ n) = lim Qni,j . n→∞ n→∞ n→∞  Ahora mostraremos un resultado que nos permite encontrar la cantidad esperada de visitas a los diversos estados absorbentes. Proposición 2.12. La matriz I − Q es invertible. Si M = (I − Q) ! ∞ X 2 3 1Xn =j . Mi,j = I + Q + Q + Q + · · · = Ei n=0 −1 entonces 7. Cadenas absorbentes Si ti = Ei ( P∞ n=0 48 1Xn ∈T ) entonces t = M 1. Si Bi,j = Pi (XT = j) entonces B = M R. Demostración. Debemos mostrar que si (I − Q) x = 0 entonces x = 0. Sin embargo, si x = Qx entonces x = Qn x → 0, por lo que x = 0. Por otra parte, puesto que
(I − Q) 1 + Q + Q2 + · · · + Qn = I − Qn+1 , al multiplicar por M = (I − Q) −1 vemos que 1 + Q + Q2 + · · · + Qn = M I − Qn+1 y al tomar el lı́mite conforme n → ∞, vemos que
Mi,j = I + Q + Q2 + Q3 + · · · . Al utilizar el teorema de convergencia dominada, vemos que X n Qni,j = X Ex (1Xn =j ) = Ex n X n 1Xn =j ! . Notemos que XX j∈T 1Xn =j = n X 1Xn ∈T n es la cantidad de veces que X está en el conjunto T , que es igual al tiempo de absorción. Al sacar esperanza, vemos que X Mi,j = M 1. ti = Ei (T ) = j Finalmente, notemos que Pi (XT = j) = XX Pi (XT = j, T = n, Xn−1 = z) n j∈T = XX n−1 Pi,k Pk,j n k∈T = X Mi,k Pk,j . k∈T  7. Cadenas absorbentes 49 Apliquemos ahora la proposición anterior al tiempo de aparición de la palabra 101 cuando p = 2/3. La matriz de transición y la matriz Q asociadas son:     1/3 2/3 0 0 1/3 2/3 0  0  2/3 1/3 0  2/3 1/3 . y Q= 0 P = 1/3 0 0 2/3 1/3 0 0 0 0 0 1 El siguiente código nos permite calcular numéricamente la matriz fundamental, la matriz B y el vector t. P = matrix (0 ,4 ,4) P [ c (1 ,3 ,10) ] <-1 / 3 P [ c (5 ,6 ,15) ] <-2 / 3 P [4 ,4] <-1 Q = P [1:3 ,1:3] R = P [1:3 ,4] I = diag (1 ,3 ,3) M = solve (I - Q ) uno = matrix (1 ,3 ,1) t = M % * % uno B=M%*%R Listing 2.9. PalabraFundamental.R Se obtiene  2.25 M = 0.75 0.75 4.50 4.50 1.50  1.50 1.50 1.50 y   8.25 t = 6.26 , 3.75 mientras que la matriz B tiene obviamente las entradas iguales a 1. Notemos que el tiempo esperado para la aparición de la palabra coincide (a nivel numérico) con el obtenido por simulación con el código Palabra.R. En el caso de la rata y el laberinto de la Figura 1, se utilizó el código siguiente para determinar la matriz fundamental M y al vector t. P = t ( matrix ( c (0 ,1 / 2 ,0 ,1 / 2 ,0 ,0 ,0 ,0 ,0 ,1 / 3 ,0 ,1 / 3 ,0 ,1 / 3 ,0 ,0 ,0 ,0 ,0 ,1 / 2 ,0 ,0 ,0 ,1 / 2 ,0 ,0 ,0 ,1 / 3 ,0 ,0 ,0 ,1 / 3 ,0 ,1 / 3 ,0 ,0 ,0 ,1 / 4 ,0 ,1 / 4 ,0 ,1 / 4 ,0 ,1 / 4 ,0 ,0 ,0 ,1 / 3 ,0 ,1 / 3 ,0 ,0 ,0 ,1 / 3 ,0 ,0 ,0 ,1 / 2 ,0 ,0 ,0 ,1 / 2 ,0 ,0 ,0 ,0 ,0 ,1 / 3 ,0 ,1 / 3 ,0 ,1 / 3 ,0 ,0 ,0 ,0 ,0 ,1 / 2 ,0 ,1 / 2 ,0) ,9) ) P [9 ,6] <-0 P [9 ,8] <-0 P [9 ,9] <-1 # Genera la matriz de transici \ ’ on ( absorbente ) para la rata en un laberinto Q = P [1:8 ,1:8] I = diag (1 ,8 ,8) M = solve (I - Q ) uno = matrix (1 ,8 ,1) t = M % * % uno Listing 2.10. LaberintoFundamental.R 7. Cadenas absorbentes Se obtuvieron  3.00 2.00  1.50  2.00 M = 1.50  1.00  1.50 1.00 50 los siguientes resultados numéricos: 3.00 3.62 2.62 2.37 2.25 1.62 1.87 1.37 1.50 1.75 2.50 1.25 1.25 1.25 1.00 0.75 3.00 2.37 1.87 3.62 2.25 1.37 2.62 1.62 3.00 3.00 2.50 3.00 3.50 2.00 2.50 2.00 1.50 1.62 1.87 1.37 1.50 2.12 1.12 0.87 1.50 1.25 1.00 1.75 1.25 0.75 2.50 1.25  1.50 1.37  1.12  1.62  1.50  0.87  1.87 2.12 y   18.00 17.00   15.00   17.00   t=  15.00 11.00   15.00 11.00 CAPÍTULO 3 Procesos de renovación Retomaremos el ejemplo de contéos aleatorios del Capı́tulo 1 (Ejemplo 1.3). Consideremos variables aleatorias independientes e idénticamente distribuidas S1 , S2 , . . . con valores estrı́ctamente postivos. Podemos pensar en que Si es el tiempo de vida de un componente crucial para el funcionamiento de cierto sistema y que al fallar se debe reemplazar. Los tiempo de reemplazo serán T0 = 0, T1 = S 1 , T2 = S1 + S2 , . . . . Por otra parte la cantidad de componentes que han fallado hasta el tiempo t serı́a Nt = min {n : Tn > t} = ∞ X n1Tn ≤t<Tn+1 . n=0 A este modelo general se le llama modelo de renovación. La sucesión S será la sucesión de tiempos de vida, la sucesión T la de tiempos de renovación y la sucesión N será el proceso de contéo asociado. Hay un par de procesos adicionales que son útiles e interesantes: el proceso de tiempo residual (hasta el próximo reemplazo) es Rt = TNt +1 − t, el proceso de edad es At = t − T N t . Su suma es el proceso de tiempos totales Lt = SNt +1 . En la Figura 1 se ilustran las definiciones asociadas a los procesos de renovación. Imaginemos ahora que en vez de cambiar al componente crucial en cuanto falla, se realiza una revisión diaria en la que se ve si se cambia o no. Entonces es más conveniente medir el tiempo en dı́as en vez de continuamente. En términos del modelo, se puede imponer simplemente que el tiempo de vida Si sea una variable aleatoria con valores en {1, 2, . . .} y las definiciones tienen sentido como las hemos puesto, salvo que Rn , An y Ln , con n ∈ N, toman valores en N y determinan a los procesos R, A y L. En este caso hablamos de un proceso de renovación aritmético. 51 52 4 3 2 Lt 1 0 At T1 Rt t Nt T2 T3 T4 Figura 1. Ilustración de las definiciones de proceso de renovación En este capı́tulo se hará una introducción a los procesos de renovación. Se enunciarán los resultados tanto en el caso general como en el aritmético, pero las justificaciones se harán en general en el caso aritmético. Ejemplo 3.1. Imaginemos la trayectoria de la rata en el laberinto de la Figura 1, que hemos modelado mediante una cadena de Markov. Supongamos que inicialmente le damos de comer a la rata en la casilla central y que, al terminar, la rata se va a dar un paseo aleatorio. Cuando regresa a la casilla central, encontrará la comida de nuevo dispuesta. En este caso, nuestros tiempos entre sucesos Si serán la cantidad de pasos entre dos comidas de la rata y nuestros tiempos de reemplazo (o de renovación) serán los instantes de las visitas sucesivas de la rata a la casilla central. La variable Rn se puede interpretar como la cantidad de tiempo que le falta a la rata, después de n pasos, para volver a comer. La variable An es la cantidad de tiempo que lleva la rata sin comer al paso n mientras que Ln es la cantidad total de pasos que pasará la rata sin comer desde la última vez que comió anterior al paso n, hasta la siguiente vez que lo hará. No es completamente trivial verificar que la situación descrita corresponde a un fenómeno de renovación, pues no hemos discutido por qué los tiempos entre las visitas sucesivas de la rata a la casilla central conforman una sucesión de variables independientes e idénticamente distribuidas. Sin embargo, la propiedad de Markov fuerte nos permite ver que ası́ es. 53 El tipo de preguntas a las que responde la teorı́a de renovación en este contexto son las siguientes: Al instante n: ¿cuántas veces ha comido la rata? ¿ Qué pasa conforme n → ∞? ¿ Qué pasa en promedio? ¿ Cuál es la distribución del tiempo que le falta para volver a comer (Rn )? ¿ Qué le pasa a esta distribución conforme n → ∞? ¿Cuál es la probabilidad de que al paso n la rata coma? ¿Qué pasa con dicha probabilidad conforme n → ∞? El ejemplo anterior es caso particular de uno mucho más general: si 0 = T0 < T1 < · · · son los instantes sucesivos en que una cadena de Markov recurrente visita a su estado inicial, que fijamos igual a x, entonces las variables Ti − Ti−1 son independientes e idénticamente distribuidas por lo que conforman un fenómeno de renovación. Si X0 = y 6= x, entonces la distribución de S1 es distinta a la de S2 , S3 , . . ., aunque aún ası́ son iid. En este caso hablamos de un proceso de renovación demorado. Ejemplo 3.2 (El proceso Bernoulli). Se trata de un proceso de renovación aritmético en el cual los tiempos entre sucesos tienen distribución geométrica: P(Si = k) = p (1 − p) k−1 para k = 1, 2, . . .. En este caso particular se pueden calcular las distribuciones de los procesos asociados al fenómeno de renovación. Además, admite una interpretación adicional: sean B1 , B2 , . . . variables aleatorias Bernoulli de parámetro p y sean T0 = 0 y Tn+1 = min {i > Tn : Bi = 1} . Entonces T es un proceso Bernoulli en el sentido de que la sucesión (Tn+1 − Tn ) es iid con distribución geométrica (concentrada en {1, 2, . . .}) de parámetro p. Comencemos con la distribución de Tn : de acuerdo a nuestra interpretación, Tn es el tiempo en el que ocurre el enésimo éxito de la sucesión Bernoulli B, por lo que Tn es binomial negativa de parámetros n y p con valores en {n, n + 1, · · · }, es decir:   k−1 k k−n P(Tn = k) = p (1 − p) . k−n Al proceso de contéo asociado también se le puede calcular la distribución exacta: notemos que Nn es la cantidad de unos en la sucesión Pn B1 , . . . , Bn , y como estas variables toman los valores cero y uno, pues Nn = i=1 Bi . Ası́, Nn tiene distribución binomial de parámetros n y p. Esto además implica que E(Nn ) = p/n y por lo tanto E(Nn /n) → p conforme n → ∞, que es un caso particular del teorema de renovación elemental que demostraremos más adelante. Respecto al proceso de tiempos resiguales, notemos que {Rn = k} = {Bn+1 = 0, . . . , Bn+k−1 = 0, Bn+k = 1} de donde concluimos que P(Rn = k) = p (1 − p) k−1 k = 1, 2, . . . . 54 En otras palabras, Rn es geométrica de parámetro p con valores en {1, 2, . . .}. Un argumento similar funciona para el proceso de edad pues {An = k} = {Bn = 0, Bn−1 = 0, . . . , Bn−k+1 = 0, Bn−k = 1} por lo que P(An = k) = p (1 − p) k k = 0, 1, . . . . En otras palabras, An es geométrica de parámetro p con valores 0, 1, 2, . . .. De hecho, las variables An y Rn son independientes pues al ser las Bi independientes vemos que P(An = j, Rn = k) = P(Bn = 0, . . . , Bn−j−1 = 0, Bn−j = 1, Bn+1 = 0, . . . , Bn+k−1 = 0, Bn+k = 1) = P(Bn = 0, . . . , Bn−j−1 = 0, Bn−j = 1) P(Bn+1 = 0, . . . , Bn+k−1 = 0, Bn+k = 1) = P(An = j) P(Rn = k) Finalmente, la distribución del proceso de tiempos totales se puede entonces calcular; recordemos que la suma de dos geométricas independientes es binomial negativa. Sólo hay que tener cuidado pues este cálculo asume que ambas variables geométricas toman valores en el mismo conjunto. El resultado es que P(Ln = m) = mp2 (1 − p) m−1 . m = 1, 2, . . . A continuación exploraremos otra liga entre cadenas de Markov y procesos de renovación. Proposición 3.1. En procesos de renovación aritméticos, el proceso de tiempos residuales es una cadena de Markov. Si M es el supremo del soporte de la distribución de S1 , la matriz de transición de R es irreducible en {i ∈ N : i ≤ M }. R es aperiódica si {n ≥ 1 : P(S1 = n) > 0} 6⊂ hN para toda h ≥ 2. Un proceso de renovación aritmético aperiódico es aquel para el cual se satisface la anterior condición para que R sea aperiódica. Demostración. Consideremos a la función ( i−1 i>1 f (i, r) = . r i=1 Al definir a la sucesión R̃ medante R̃0 = 1 y   R̃n+1 = f R̃n , Sn+1 , vemos que R̃ es una cadena de Markov con matriz de transición ( 1 i > 1, j = i − 1 Pi,j = . P(S1 = j) i = 1 55 Sin embargo, puesto que la sucesión S es iid, R y R̃ tienen la mismas distribuciones finito-dimensionales en el sentido siguiente:   P(R0 = i0 , . . . , Rn = in ) = P R̃0 = i0 , . . . , R̃n = in . En efecto, si i0 , . . . , in ∈ {1, 2, . . .} e il = 1 si y sólo si l ∈ I ⊂ {0, . . . , n} e il = il−1 − 1 si l 6∈ I y digamos que I = {l1 , . . . , lm } entonces     P R̃0 = i0 , . . . , R̃n = in = P Silk = ilk , k ≤ m = P(Sk = ilk , k ≤ m) = P(R0 = i0 , . . . , Rn = in ) . Esto prueba que R es una cadena de Markov con matriz de transición P . Si M es el supremo del soporte de la distribución de S1 y M < ∞, entonces P(S1 = M ) > 0 y P(S1 = n) = 0 para toda n ≥ M . Entonces de 1 se accede a M , a partir del cualse accede a M − 1, . . . , 1, por lo que {0, . . . , M } es una clase de comunicación, que de hecho es cerrada pues P(S1 = n) = 0 para toda n ≥ M . Esto hace que la cadena en {1, . . . , M } sea irreducible. Si M = ∞, entonces para toda M existe n ≥ M tal que P(S1 = n) > 0, por lo que 0 se comunica con M vı́a n. Es decir, la cadena es irreducible en {0, 1, . . .}. Si {n ≥ 1 : P(S1 = n) > 0} 6⊂ hN para h ≥ 2 entonces existen n1 , n2 , k1 y k2 naturales tal que n1 k1 = 1 + n2 k2 y tales que P(S1 = n1 ) , P(S1 = n2 ) > 0. Vemos entonces que es posible ir de cero a cero en n1 k1 = n2 k2 + 1 pasos y en n2 k2 pasos, por lo que el periodo de 1 es 1 y por lo tanto la cadena es aperiódica.  La respuesta a cuántas veces a comido la rata al paso n se puede interpretar en términos de Nn en el caso de procesos de renovación aritméticos o de Nt en el caso general. Los siguientes dos resultados nos permiten dar una posible respuesta: para tiempos grandes, la variable Nt se puede predecir determinı́sticamente. Proposición 3.2 (Ley fuerte de los grandes números). Si µ = E(Si ) < ∞ entonces Nt /t → 1/µ casi seguramente Demostración. Por la ley fuerte de los grandes números, sabemos que Tn /n → µ casi seguramente. Notemos que si t ∈ [Tn , Tn+1 ) entonces Tn+1 n + 1 t Tn n ≤ ≤ . n n+1 Nt n+1 n Los extremos de la desigualdad convergen al mismo lı́mite, µ, conforme n → ∞ de manera casi segura. Por lo tanto Nt /t → 1/µ.  Nuestro siguiente resultado involucra a la llamada función de renovación. Es la función m : [0, ∞) → [0, ∞) dada por m(t) = E(Nt ). Proposición 3.3 (Teorema de renovación elemental). Si µ = E(S1 ) < ∞ entonces m(t) /t → 1/µ. 56 La demostración utilizará un resultado sobre camintas aleatorias conocido como la identidad de Wald y que se prueba como una aplicación de la teorı́a de martingalas. Se trata de una generalización de la fórmula E(Tn ) = nE(T1 ) al caso en que n ya no es una cantidad fija sino aleatoria. En nuestro contexto, veremos que (4) E(TNt +1 ) = E(Nt + 1) E(T1 ) = E(Nt + 1) µ. La demostración en este caso particular es la siguiente: notamos que {Nt + 1 = n} = {Nt = n − 1} = {Tn−1 ≤ t < Tn } . Se sigue que {Nt + 1 > n} pertenece a la σ-álgebra generada por S1 , . . . Sn y es por lo tanto independiente de Sn+1 . Al aplicar el teorema de Tonelli (dos veces): ! X 1i=Nt +1 Ti E(TNt +1 ) = E i   = E XX 1i=Nt +1 Sj  = E XX 1i=Nt +1 Sj  = E X   = X j =µ i j j j≤i i≥j   1Nt +1≥j Sj  E(Sj ) P(Nt + 1 ≥ j) X j P(Nt + 1 ≥ j) = µE(Nt + 1) . Prueba del teorema de renovación elemental. Ocupémonos primero de obtener una cota inferior. A partir de la identidad de Wald, vemos que t ≤ E(TNt +1 ) ≤ E(Nt + 1) µ. Por lo tanto 1 E(Nt ) 1 − ≤ . µ t t Notemos que lim t→∞ 1 1 1 − = . µ t µ 1. La medida de renovación en el caso aritmético 57 Para obtener una cota superior, definimos a S̃i = Si ∧b. Con la notación obvia, notemos que T̃i ≤ Ti y que por lo tanto Nt ≤ Ñt . Además,     E S̃1 ∧ b E(Nt + 1) = E T̃Ñt +1 ≤ t + b, por lo que Al utilizar b = por lo que √ t+b E(Nt ) ≤ . t tE(S1 ∧ b) t, el teorema de convergencia monótona nos dice que  √ E S1 ∧ t → µ, √ t+ t 1 √
= . lim t→∞ tE S1 ∧ µ t Podemos entonces deducir que lim t→∞ 1 E(Nt ) = . t µ  1. La medida de renovación en el caso aritmético Continuaremos el estudio de procesos de renovación al introducir a la medida de renovación. En realidad se introducirá a partir de la función u dada por X un = P(∃m, Tm = n) = P(Tm = n) . m En términos de la rata en el laberinto, un es la probabilidad de que la rata coma al instante n. El comportamiento asintótico de un se conoce y es el objeto del teorema de renovación para procesos de renovación aritméticos de Erdős, Feller y Pollard. Teorema 3.1. Para un proceso de renovación aritmético, aperiódico y con media finita: 1 lim un → . n→∞ µ Demostración. Primero probemos que el proceso de tiempos residuales R tiene una distribución invariante. En efecto, sean πi = P(S1 ≥ i) /µ 2. La medida de renovación en el caso continuo 58 y P la matriz de transición de R. Entonces X πi Pi,j (πP )j = i = πj+1 Pj+1,j + π1 P1,j 1 P(S1 ≥ j + 1) Pj+1,j + P(S1 = j) µ µ P(S1 ≥ j) = = πj . µ Al ser la cadena irreducible y aperiódica y con una distribución invariante, es positivo recurrente y se puede aplicar el teorema fundamental de convergencia para concluir que 1 n lim Pi,1 = . n→∞ µ n Por otra parte, podemos calcular explı́citamente P1,1 pues = n P1,1 = P(Rn = 1) = P(n = Tm para alguna m) = un . Ası́, vemos que un → 1/µ conforme n → ∞.  2. La medida de renovación en el caso continuo Ahora presentaremos el resultado análogo al teorema de renovación aritméticos en el caso de procesos no aritméticos. No se abordará la prueba pues depende de resultados análogos al teorema fundamental de convergencia pero para cadenas de Markov a tiempo discreto y espacio continuo (en este caso, el espacio de estados del proceso de tiempos residuales es (0, ∞)). Sea S un proceso de renovación no-aritmético, es decir una sucesión de variables aleatorias independientes e idénticamente distribuidas con valores en (0, ∞) y tal que no toman valores en δN. Se puede pensar que Si es una variable aleatoria con densidad. Sea T la sucesión de tiempos entre sucesos y N el proceso de contéo asociado. Definición. La medida de renovación es la medida U en [0, ∞) tal que U ([0, t]) = E(Nt ) = ∞ X n=1 P(Tn ≤ t) . Se utilizará la notación Ut para U ([0, t]). En general, para procesos de renovación no aritméticos, la cantidad P(∃nTn = t) es cero para casi toda t. Es por eso que se utiliza la medida de renovación en vez de la función de renovación del caso aritmético. En efecto, si [a, b] es un intervalo muy pequeño, puede que haya un punto Tn en dicho intervalo, pero serı́a raro que hubiera más de uno. Por lo tanto U ([a, b]) ≈ P(∃n, Tn ∈ [a, b]). El siguiente teorema puede entonces interepretarse como es una extensión del teorema de renovación para procesos aritméticos. 3. La ecuación de renovación 59 Teorema 3.2 (Teorema de renovación de Blackwell). Para todo h > 0: lim Ut+h − Ut = t→∞ h . µ 3. La ecuación de renovación La ecuación de renovación aparece en el cálculo de algunas cantidades de la teorı́a de renovación. En general se trata de un análisis al primero instante de renovación al notar que el proceso T ′ dado por 0, T2 − T1 , T3 − T1 , · · · es un proceso de renovación idéntico en distribución a T0 , T1 , T2 , · · · e independiente de T1 . Ejemplos de su utilidad es que nos permite estudiar el comportamiento asintótico de las distribuciones de los procesos de tiempo residual, de edad ó de tiempos totales. Sin embargo, también se puede entender como una relación de recurrencia para cantidades de interés en teorı́a de renovación. Nos concentraremos en el caso aritmético, al ser técnicamente menos complicado, al comenzar con algunos ejemplos. Escribiremos pk = p(k) = P(S1 = k) . Ejemplo 3.3. Sea uk la densidad de renovación para un proceso de renovación aritmético: uk = P(∃n ≥ 0 tal que Tn = k) . Escribiremos a dicha cantidad como u(k) cuando la tipografı́a lo requiera. Notemos que u0 = 1, mientras que para k ≥ 1, claramente P(T0 = k) = 0 y por lo tanto: uk = P(∃n ≥ 1 tal que Tn = k) = k X j=1 = k X j=1 = k X j=1 = k X P(S1 = j, ∃n ≥ 0 tal que Tn − j = k − j) P(S1 = j, ∃n ≥ 0 tal que Tn′ = k − j) pj P(∃n ≥ 0 tal que Tn = k − j) pj uk−j j=1 = E(u(k − S1 )) , donde por supuesto se utiliza el hecho que por definición uj = 0 si j < 0. La ecuación uk = δ0 (k) + E(u(k − S1 )) 3. La ecuación de renovación 60 que también se escribe uk = δ0 (k) + X pj uk−j es la llamada ecuación de renovación para la densidad de renovación. Ejemplo 3.4. Sea Ln , n ≥ 0 el proceso de tiempos totales de un proceso de renovación aritmético. Para r ≥ 0, sea z(n) = P(Ln = r) para n ≥ 0 y z(n) = 0 si n < 0. Notemos que si n ≤ r entonces z(n) = P(S1 = r) = pr . Por otra parte, si n > r, nos fijamos en la primera renovación T1 = S1 , que si Ln = r es necesariamente es más chica que n. Entonces se obtiene: X X z(n) = P(S1 = j, Ln−j (T ′ ) = r) = pj z(n − j) . j≤n j≤n La ecuación de renovación resultante es z(n) = 1n≤r pr + X j pj z(n − j) . Ejemplo 3.5. Sea Rn el proceso de tiempos residuales de un proceso de renovación aritmético. Sea z(n) = P(Rn = r) si n ≥ 0 y z(n) = 0 si n < 0. El cálculo de Rn se puede dividir en dos casos, dependiendo de si S1 ≤ n o no: z(n) = P(Rn = r) = pn+r + n X j=1 P(S1 = r) z(n − r) . En general, la ecuación de renovación es la siguiente: dada una función b (tal que bn = 0 si n < 0) se busca una función z tal que z(n) = 0 si n < 0 y X z(n) = b(n) + pk z(n − j) . j≤n La solución se puede encontrar al iterar la ecuación de renovación para escribir: z(n) = b(n) + E(z(n − S1 )) = b(n) + E(b(n − S1 )) + E(b(n − S1 − S2 )) = b(n) + E(b(n − S1 )) + E(b(n − S1 − S2 )) + E(b(n − S1 − S2 − S3 )) = · · · . Puesto que b(n − Tm ) = 0 si m ≥ n, vemos que X z(n) = E(b(n − Tm )) m y al sumar sobre los posibles valores de Tm , se obtiene XX X z(n) = b(n − x) P(Tn = x) = ux bn−x . m x x 3. La ecuación de renovación 61 Una vez establecida la relación entre soluciones a la ecuación de renovación y la densidad de renovación, se tienen los elementos clave para probar el siguiente resultado. Teorema 3.3 (Teorema clave de renovación en el caso discreto). Si z es solución a la ecuación de renovación X z(n) = b(n) + pk z(n − j) . j≤n yPb es sumable entonces z está dada por z(n) = x b(x) /µ conforme n → ∞. P x b(x) u(n − x) y z(n) → Demostración. Sólo hace falta verificar el comportamiento asintótico de z. Sin embargo, si b es sumable, puesto que un es una probabilidad, se tiene la cota X X b(x) u(n − x) ≤ b(x) < ∞ x x y por lo tanto el teorema de convergencia dominada y el teorema de renovación de Erdős-Feller-Pollard (EFP) nos dicen que X X lim z(n) = lim b(x) u(n − x) = b(x) /µ. n→∞ n→∞ x x  Veamos ahora algunas aplicaciones del teorema de renovación clave. En el caso del tiempo total asintótico, vemos que X lim P(Rn = r) = pr+x /µ = P(S1 ≥ r) /µ, n→∞ x aunque en realidad esto ya lo sabı́amos y lo utilizamos en la prueba del teorema de renovación de EFP. Un ejemplo más interesante es el del tiempo total: X lim P(Ln = r) = px /µ = rpr /µ. n→∞ x≤r CAPÍTULO 4 Procesos de Poisson Imaginemos la siguiente situación: al tiempo cero (inicial) una serie de personas contratan un seguro por el que pagan una prima. La compañia de seguros tienen entonces la obligación de dar una compensación económica a los contratantes en caso de que suceda un cierto evento llamado siniestro. En muchos casos, el monto que debe pagar la compañı́a por el siniestro también es aleatorio. El modelo matemático se realiza en dos etapas: primero nos enfocamos en los instantes en que suceden los siniestros y luego en sus montos. Si 0 = T0 < T1 < T2 < · · · son los instantes en que acaecen los siniestros, se introduce al proceso de contéo asociado N = (Nt , t ≥ 0), donde Nt nos dice cuántos siniestros han sucedido al instante t, mismo que está caracterizado por satisfacer las identidades X Nt = n1Tn ≤t<Tn+1 , n {Nt = n} = {Tn ≤ t < Tn+1 } y {Nt ≥ n} = {Tn ≤ t} . En general, un proceso de contéo es un proceso con trayectorias constantes por pedazos que toma valores enteros y va incrementando de uno en uno. A Si = Ti − Ti−1 se le llama i-ésimo tiempo interarribo. Es más natural imponer supuestos sobre N que sobre los tiempos Ti . Los supuestos se pueden traducir como sigue: Incrementos estacionarios: La distribución de Nt+s − Nt sólo depende de s (en otras palabras, si los siniestros se presentan homogéneamente en el tiempo, se tiene la igualdad en distribución; en particular no hay fenómenos estacionales). Incrementos independientes: Si 0 = t0 ≤ t1 ≤ · · · ≤ tn , las variables Nt1 − Nt0 , . . . , Ntn − Ntn−1 son independientes. (Ası́, el que ocurran siniestros por la mañana no afecta lo que ocurre por la tarde. Definición. Un proceso de Poisson es un proceso de contéo con incrementos independientes y estacionarios. 62 63 Esta definición captura la idea de que el riesgo de que ocurra un siniestro se distribuye homogéneamente tanto en el tiempo como en la población. El primer objetivo de este capı́tulo es mostrar como una tal definición tiene implicaciones prácticas y teóricas que hacen que podamos hacer cálculos con el proceso de Poisson que no es posible hacer explı́citamente con otros modelos. En efecto, hay una definición alternativa de proceso Poisson que automáticamente nos deja hacer cálculos. Teorema 4.1. Un proceso de contéo N es de Poisson si y sólo si existe λ > 0 tal que los tiempos interarribo son variables exponenciales independientes de parámetro λ. Al parámetro λ se le conoce como intensidad del proceso de Poisson. Antes de probar el teorema, veámos algunas de sus consecuencias. La primera es que Nt es el proceso de contéo asociado a un proceso de renovación con tiempos interarribo exponenciales. Por lo tanto, se satisface la ley fuerte de los grandes números en particular. Además, se pueden calcular explı́citamente las distribuciones al tiempo t del proceso de edad, de tiempos residuales y de tiempos totales, aún en el contexto de procesos a tiempo continuo. Una segunda consecuencia es que se conoce la distribución exacta del proceso de contéo. Proposición 4.1. Sea N un proceso de Poisson de parámetro λ. Entonces Nt tiene distribución Poisson de parámetro λt. Demostración. Puesto que los tiempos interarribo son exponenciales independientes de parámetro λ, entonces Tn tiene distribución Γ de parámetros λ y n y es independiente de Sn+1 . Ası́, vemos que la densidad conjunta de (Tn , Sn+1 ) está dada por λn tn−1 e−λt e−λs . fTn ,Sn+1 (t, s) = n − 1! λ Por lo tanto: P(Nt = n) = P(Tn ≤ t < Tn + Sn+1 ) Z tZ ∞ fTn ,Sn+1 (t1 , s) ds dt1 = 0 = Z tZ 0 = Z t t−t1 ∞ λn t1n−1 e−λt1 −λs λe . ds dt1 n − 1! t−t1 λn tn−1 e−λt1 −λ(t−t1 ) 1 n − 1! n (λt) = e−λt n! e 0  64 Como corolario, podemos calcular las distribuciones finito dimensionales de N : sean 0 = t0 < t1 < · · · < tm y n0 = 0 ≤ n1 ≤ · · · ≤ nm por lo que P(Nt1 = n1 , . . . , Ntm = nm ) = P Nt1 − Nt0 = n1 − n0 , . . . , Ntm − Ntm−1 = nm − nm−1 y al utilizar la independencia de incrementos = m Y i=1 P Nti − Nti−1 = ni − ni−1

ası́ como su estacionariedad = = m Y i=1 m Y P Nti −ti−1 = ni − ni−1
e−λ(ti −ti−1 ) (λ (ti − ti−1 )) (ni − ni−1 )! i=1 ni −ni−1 . Ya con esto, podemos hacer un primer código para simular al proceso de Poisson. De hecho se trata simplemente del código Poisson.R que ya se habı́a introducido: se trata de generar variables exponenciales e irlas sumando hasta que se sobrepase el tiempo t: lambda =10 # Intensidad xi = rexp (1 , lambda ) # xi representa el tiempo del primer siniestro T = c (0 , xi ) # El vector T irá acumulando los tiempos en que van ocurriendo siniestros N =0 # N nos dirá cuantos siniestros ocurren hasta el tiempo 1 while ( xi <1) { # Mientras no se haya sobrepasado el instante 1 N <-N +1 # Aumentamos el número de siniestros xi <- xi + rexp (1 , lambda ) # Vemos el tiempo en el que ocurre el siguiente siniestro T = c (T , xi ) # Aumentamos un evento temporal } plot (T , c (1:( N +2) ) ) Listing 4.1. Poisson2.R Para ver otro esquema de simulación para el proceso de Poisson, en el que el ciclo while se substituye por un ciclo for, haremos algunos cálculos adicionales. La primera pregunta es: si sabemos que han ocurrido n siniestros en [0, t], ¿en dónde han caido estos saltos? Para responder a dicha pregunta, calcularemos la densidad 65 conjunta de T1 , . . . , Tn dado que Nt = n. Puesto que P(T1 ∈ dt1 , . . . , Tn ∈ dtn , Nt = n) = P(T1 ∈ dt1 , . . . , Tn ∈ dtn , tn ≤ t ≤ tn + Sn+1 ) = n Y λe−λ(ti −ti−1 ) e−λ(t−tn ) i=1 vemos que la densidad condicional buscada es: fT1 ,...,Tn |Nt =n (t1 , . . . , tn ) = n Y λe−λ(ti −ti−1 ) e−λ(t−tn ) i=1 n! n (λt) e−λt n! . tn Se concluye que condicionalmente a Nt = n, las variables T1 , . . . , Tn los valores ordenados de n uniformes independientes en [0, t]. Por lo tanto, hay otra manera de simular la trayectoria de un proceso de Poisson en el intervalo de tiempo [0, t]: = l <- 10 % Intensidad t <-3 % Umbral de tiempo N <- rpois (1 , l * t ) % Cantidad de siniestros u <- runif ( N ) * t % Tiempos en que ocurren ( desordenados ) u <- sort ( u ) % Tiempos de siniestros ordenados plot ( u ) Listing 4.2. PoissonConUniformes.R Antes de continuar con el modelo completo de reclamaciones en una compañı́a de seguros, veremos algunos procesos estocásticos importantes asociados al proceso de Poisson. Se definirá a Ft = σ(Xs : s ≤ t) . Ejemplo 4.1 (El proceso de Poisson compensado). El proceso de Poisson compensado es el proceso Mt = Nt − λt. Este proceso satisface la siguiente propiedad: si s ≤ t E( Mt | Fs ) = Ms . En efecto, al notar que Mt = Mt − Ms + Ms vemos que E( Mt | Fs ) = E( Mt − Ms | Fs ) + E( Ms | Fs ) = E(Mt − Ms ) + Ms = Ms puesto que Mt − Ms tiene media cero y es independiente de Fs . La propiedad anterior nos dice que M es una martingala. Las martingalas son fundamentales en el desarrollo de la probabilidad moderna. Ejemplo 4.2 (La martingala exponencial). Consideremos al proceso −q Et = e−qNt +tλ(1−e ) . 66 Este proceso también es una martingala. Para verificarlo, primero calculemos ∞ n
X −q −q e−λt (λt) = e−λt eλte = e−tλ(1−e ) . e−qn E e−qNt = n! n=0 Luego utilizamos la independencia entre Nt − Ns y Fs para obtener:   −q E( Et | Fs ) = etλ(1−e ) e−qNs E e−q(Nt −Ns ) Fs −q −q = etλ(1−e ) e−qNs e−(t−s)λ(1−e ) = Es . Una consecuencia del cálculo de la transformada de Laplace de Nt es que podemos calcular la esperanza y la varianza de Nt : si definimos
−q ϕt (q) = E e−qNt = e−tλ(1−e ) , y derivamos bajo el signo integral cuando q > 0, se obtiene

∂ϕt (q) ϕt (q) −λte−q = = −E Nt e−qNt . ∂q Luego el teorema de convergencia monótona nos permite ver que λt = E(Nt ) . Al derivar dos veces se tiene que para q > 0:
2
∂ϕt (q) ϕt (q) λte−q + ϕt (q) λte−q = = E Nt2 e−qNt , ∂q por lo que el teorema de convergencia monótona nos dice que
2 E Nt2 = (λt) + λt. Finalmente, se obtiene Var(Nt ) = λt. Continuaremos ahora con el modelo de las reclamaciones en una compañı́a de seguros: es natural que los montos de los siniestros son independientes entre si y de los tiempos en que acontecen. Por lo tanto se plantéa el siguiente modelo completo: (1) Los tiempos en que acontecen los siniestros tienen como proceso de contéo asociado a un proceso de Poisson N . Sea λ su intensidad. (2) Los montos de los siniestros son variables aleatorias independientes e idénticamente distribuidas ξ1 , ξ2 , . . .. Sea Sn = ξ1 + · · · + ξn la sucesión de sumas parciales asociada a los montos de los siniestros. Entonces, al tiempo t, el monto total que ha sido reclamado a la compañı́a de seguros es X X t = SN t = ξi . i≤Nt 67 Al proceso X se le conoce con el nombre de proceso de Poisson compuesto. A la distribución de ξi se le conoce como distribución de salto de X y al parámetro de N se le conoce como intensidad de X. El siguiente código permite simular al proceso de Poisson compuesto cuando la distribución de salto es la de 100,000 veces una Beta de parámetros 1/2 y 1/2. La idea es modelar una situación en la que es más probable que los montos de los siniestros sean o muy pequeños o muy grandes (hasta el tope de 100,000 que podrı́a ser la suma asegurada). El resultado se puede apreciar en la Figura 1 lambda <- 10 # Intensidad s <- rexp (1 , lambda ) # s representa el tiempo del primer siniestro x <- 100000 * rbeta (1 ,1 / 2 ,1 / 2) # x representa el monto del primer siniestro T <-c (0 , s ) # El vector T irá acumulando los tiempos en que van ocurriendo siniestros X <-c (0 , x ) # X irá acumulando los montos de los siniestros N <-0 # N nos dirá cuantos siniestros ocurren hasta el tiempo 1 while (s <1) { # Mientras no se haya sobrepasado el instante 1 N <-N +1 # Aumentamos el número de siniestros s <-s + rexp (1 , lambda ) # Vemos el tiempo en el que ocurre el siguiente siniestro x <- 100000 * rbeta (1 ,1 / 2 ,1 / 2) # Calculamos su monto T <-c (T , s ) # Aumentamos un evento temporal X <-c (X , tail (X ,1) + x ) # Agregamos el monto de la reclamación } plot (T , X ) Listing 4.3. PoissonCompuesto.R Se pueden hacer cálculos paralelos a los del proceso de Poisson en el caso Poisson compuesto. Por ejemplo, calculemos la media y varianza de Xt : si µ = E(S1 ) y σ 2 = Var(S1 ) son finitas entonces E Xt2
= ∞ X E Sn2 1(Nt =n) n=0
= ∞ X n=0

nσ 2 + n2 µ2 P(Nt = n) = E σ 2 Nt + µ2 Nt2 < ∞, ya que las variables Poisson tienen momentos de cualquier orden. La media está dada por E(Xt ) = ∞ X n=0 ∞
X nµP(Nt = n) = λtµ = E(Nt ) E(S1 ) . E Sn2 1(Nt =n) = n=0 Ası́, también se tiene que    2 Var(Xt ) = E E (SNt − λµt) Nt    2 2 = E E (SNt − Nt µ) + (Nt µ − λµt) Nt   2 = E Nt σ 2 + µ2 (Nt − λt) = σ 2 λt + µ2 λt. 150000 0 50000 100000 X 200000 250000 68 0.0 0.2 0.4 0.6 0.8 1.0 T Figura 1. Trayectoria de un proceso de Poisson compuesto Las igualdades anteriores también se pueden interpretar como un caso particular de las identidades de Wald. Continuaremos ahora con un análisis del Teorema 4.1. Primero probaremos que un proceso de contéo con tiempos interarribo independientes y exponenciales de parámetro λ tiene incrementos independientes y estacionarios. Sean S1 , S2 , . . . exponenciales independientes de parámetro λ, T las sumas parciales asociadas dadas por T0 = 0, Tn = S1 + · · · + Sn y N = (Nt , t ≥ 0) el proceso de contéo asociado dado por Nt = X n1Tn ≤t<Tn+1 . Consideremos al proceso N t dado por Nst = Nt+s − Nt . Este proceso es un proceso de contéo. Sus tiempos interarribo son S1t = TNt +1 − t, S2t = SNt +2 , S3t = SNt +3 , . . . . 69 Al descomponer respecto al valor de Nt vemos que
P S1t > t1 , S2t > t2 , . . . , Snt > tn X P(Nt = m, Sm+1 − t > t1 , Sm+2 > t2 , . . . , Sm+n > tn ) = m = X m = X m P(Tm ≤ t < Tm + Sm+1 , Tm + Sm+1 − t > t1 , Sm+2 > t2 , . . . , Sm+n > tn ) P(Tm ≤ t < Tm + Sm+1 , Tm + Sm+1 − t > t1 ) P(Sm+2 > t2 , . . . , Sm+n > tn ) . Al condicionar sobre el valor de Tm , su independencia con Sm+1 y utilizar la propiedad de pérdida de memoria para S1 , vemos que P(Tm ≤ t < Tm + Sm+1 , Sm+1 − t > t1 ) = E E 1Tm ≤t<Tm +Sm+1 1Tm +Sm+1 −t>t1 T1   = E 1Tm ≤t e−λ(t−Tn ) e−λt1 .

Sin embargo, al utilizar t1 = 0, nos damos cuenta de que   P(Nt = m) = P(Tm ≤ t < Tm + Sm+1 ) = E 1Tm ≤t e−λ(t−Tn ) , por lo cual se deduce que
X P(Nt = m) e−λt1 e−λt2 · · · e−λtn . P S1t > t1 , S2t > t2 , . . . , Snt > tn = m Se deduce que los tiempos de salto de N t son exponenciales independientes de parámetro λ. Se puede concluir que N t tiene las mismas distribuciones finitodimensionales que N , es decir, que si 0 = t0 ≤ t1 ≤ · · · entonces
P(Nt1 = m1 , . . . , Ntn = mn ) = P Ntt1 = m1 , . . . , Nttn = mn . (Esto se verifica al expresar los eventos anteriores en términos de los tiempos interarribo.) En particular, N t tiene incrementos estacionarios. Ahora veamos que N t tiene incrementos independientes. Para esto, consideremos al evento A = {Ns1 = k1 , . . . , Nsm = km } , donde s1 , . . . , sn ≤ t. Al descomponer sobre el valor de Nt se obtiene
P A, Nt = m, S1t > t1 , S2t > t2 , . . . , Snt > tn X = P(A, Tm ≤ t < Tm + Sm+1 , Tm + Sm+1 − t > t1 ) P(Sm+2 > t2 , . . . , Sm+n > tn ) . m Sin embargo, como el conjunto A ∩ {Nt = m} se escribe en términos de S1 , . . . , Sm y de {Nt = m}, podemos repetir el argumento anterior para verificar que P(A, Tm ≤ t < Tm + Sm+1 , Tm + Sm+1 − t > t1 ) = P(A ∩ {Nt = m}) e−λt1 . 70 Por lo tanto, P A, Nt = m, S1t > t1 , S2t > t2 , . . . , Snt > tn

= P(A, Nt = m) P S1t > t1 , S2t > t2 , . . . , Snt > tn . Lo anterior implica que Nst es independiente de Ns1 , . . . , Nsm si s1 , . . . , sm ≤ t, y por inducción se puede probar entonces que N tiene incrementos independientes. Ası́, N es un proceso de contéo con incrementos independientes y estacionarios. Consideremos ahora un proceso de contéo con incrementos independientes y estacionarios y veamos que se trata de un proceso de renovación con tiempos interarribo exponenciales. De la propiedad de incrementos independientes y estacionarios, podemos inmediatamente verificar que T1 tiene una distribución exponencial. En efecto, notemos que P(T1 > t + s) = P(Nt = 0, Nt+s − Nt = 0) = P(Nt = 0) P(Ns = 0) = P(T1 > t) P(T1 > s) Ası́, o T1 es infinito casi seguramente, lo cual dice que N es idénticamente cero, ó es exponencial de algún parámetro λ > 0. Ahora sólo debemos ver la razón por la que los tiempos interarribo son independientes y con la misma distribución que el primero. Para esto, analizaremos dos implicaciones de nuestras hipótesis conocidas como la propiedad de Markov y de Markov fuerte. Proposición 4.2 (Propiedad de Markov del proceso de Poisson). El proceso N t dado por Nst = Nt+s − Nt es un proceso de Poisson con las mismas distribuciones finito-dimensionales que N y es independiente de Ft . Demostración. Lo primero es notar que N t es un proceso de Lévy y de contéo. Para ver que tiene las mismas distribuciones finito-dimensionales que N , notamos que   Ntt1 , Ntt2 − Ntt1 , . . . , Nttn − Nttn−1
= Nt+t1 − Nt , Nt+t2 − Nt+t1 , . . . , Nt+tn − Nt+tn−1
d = Nt1 , Nt2 − Nt1 , . . . , Ntn − Ntn−1 . Al considerar f (x1 , . . . , xn ) = (x1 , x1 + x2 , . . . , x1 + · · · + xn ), vemos que entonces  
Ntt1 , Ntt2 , . . . , Nttn = f Ntt1 , Ntt2 − Ntt1 , . . . , Nttn − Nttn−1
d = f Nt1 , Nt2 − Nt1 , . . . , Ntn − Ntn−1 = (Nt1 , Nt2 , . . . , Ntn ) . 71 Ahora sólo se debe probar que N t es independiente de Ft . Sin embargo, si 0 = s0 ≤ s1 ≤ sm ≤ t = t0 ≤ t1 ≤ tn entonces las variables
Ns1 − Ns0 , . . . , Nsm − Nsm−1 y 
 Nt1 − Nt0 , . . . , Ntm − Ntm−1 = Ntt1 , . . . , Nttm −t − Nttm−1 −t son independientes. Eso implica que (Ns1 , . . . , Nsm ) y  Ntt1 , . . . , Nttm −t − Nttm−1 −t  son independientes y que por lo tanto, si definimos a los π-sistemas C1 = {Ns1 ∈ A1 , . . . , Nsm ∈ Am } y  C2 = Ntt1 ∈ B1 , . . . , Nttn ∈ Bn se sigue que M1 = {A ∈ Ft : A ⊥ B para todo B ∈ C2 } contiene a C1 . Puesto que M1 es un λ-sistema y σ(C1 ) = Fs entonces Fs ⊂ M1 . Por lo tanto el λ-sistema M2 = {B ∈ Ft : A ⊥ B para todo A ∈ Fs } contiene al π-sistema C2 . Esto implica que σ(Nst : s ≥ 0) ⊂ M2 y que por lo tanto N t es independiente de Ft .  Ahora, haremos una extensión adicional de la propiedad anterior. Proposición 4.3 (Propiedad de Markov fuerte). Si T : Ω → [0, ∞) es una variable aleatoria tal que {T ≤ t} ∈ Ft para toda t ≥ 0 y FT = {A ∈ F : A ∩ {T ≤ t} ∈ Ft ∀t} entonces el proceso N T dado por NtT = NT +t − NT es un proceso de Lévy y de contéo con las mismas distribuciones finito-dimensionales que N e independiente de FT . Veamos cómo se aplica la propiedad de Markov fuerte: puesto que los tiempos interarribo de N Tn son Tn+1 −Tn , Tn+2 − Tn+1 , . . . y son independientes de FTn , se deduce que Tn+1 − Tn es exponencial e independiente de T1 , T2 − T1 , . . . , Tn − Tn−1 . Esto termina la prueba de que un proceso de contéo y de Lévy es el proceso de contéo asociado a un proceso de renovación exponencial. 72 Prueba de la propiedad de Markov fuerte. Para cada n ≥ 1, definimos T n igual a (k + 1)/2n si T ∈ [k/2n , (k + 1)/2n ). Formalmente, podemos escribir T n = ⌈2n T ⌉/2n . Entonces T n es una sucesión de variables aleatorias que decrecen hacia T . Además, notemos que {T n = k + 1/2n } = {k/2n ≤ T < (k + 1) /2n } = {T < (k + 1) /2n } \ {T < k/2n } ∈ F(k+1)/2n . Por el mismo argumento, si A ∈ FT entonces A ∩ {T n = k/2n } ∈ F(k+1)/2n . Ası́, por la propiedad de Markov, vemos que   n m P A, T n = k + 1/2n , NtT1 = k1 , . . . , NtTn = km = P(A, T n = k + 1/2n ) P(Nt1 = k1 , . . . , Ntm = km ) . Al sumar sobre k, vemos que   n n P A, NtT1 = k1 , . . . , NtTm = km = P(A) P(Nt1 = k1 , . . . , Ntm = km ) . Como N es continuo por la derecha y T n decrece a T , vemos que confome n → ∞:
T T P A, Nt1 = k1 , . . . , Ntm = km =P(A)P(Nt1 = k1 , . . . , Ntm = km ). Se concluye que N T es un proceso de contéo y de Lévy con las mismas distribuciones finitodimensionales que N . Por lo tanto, sus tiempos interarribo tienen la mismas distribuciones conjuntas, pues si S̃ son los tiempos interarribo de N T entonces   P S̃1 > s1 , S̃2 > s2 , . . . , S̃n > sn   = P T̃1 > s1 , T̃2 > s1 + s2 , . . . , T̃n > s1 + · · · + sn
= P NsT1 ≤ 0, NsT1 +s2 ≤ 1, . . . , NsT1 +···+sn ≤ n − 1 = P(Ns1 ≤ 0, Ns1 +s2 ≤ 1, . . . , Ns1 +···+sn ≤ n − 1) = P(S1 > s1 , S2 > s2 , . . . , Sn > sn ) .  Proposición 4.4. El proceso de Poisson compuesto es un proceso de Lévy y satisface la propiedad de Markov siguiente: si FtX = σ(Xr : r ≤ t) entonces el proceso X t dado por Xst = Xt+s − Xt es un proceso de Poisson compuesto con la misma intensidad y distribución de salto que X y es independiente de FtX . Demostración. El proceso X comienza en cero y tiene trayectorias que son constantes por pedazos, en particular, continuas por la derecha y con lı́mites por la izquierda. 73 Para ver que X tiene incrementos independientes y estacionarios, notemos que si 0 = t0 ≤ t1 ≤ · · · ≤ tn y Ai ∈ BR para i = 1, . . . , n, entonces
P Xti − Xti−1 ∈ Ai , i = 1, . . . , n X
= P Ski − Ski−1 ∈ Ai , Nti = ki , i = 1, . . . , n 0=k0 ≤k1 ≤···≤kn = X n Y 0=k0 ≤k1 ≤···≤kn i=1
P Ski − Ski−1 ∈ Ai P(Nti = ki , i = 1, . . . , n) , puesto que las variables {Si : i = 0, 1, · · · } y {Nt : t ≥ 0} son independientes. Si k2 ≥ k1 , entonces Sk2 − Sk1 tiene la misma distribución que Sk1 −k2 , entonces n X Y
P Ski − Ski−1 ∈ Ai P(Nti = ki , i = 1, . . . , n) = = = 0=k0 ≤k1 ≤···≤kn i=1 n X Y 0=k0 ≤k1 ≤···≤kn i=1 n X Y 0≤j1 ,...,jn i=1 n Y i=1 P Ski −ki−1 ∈ Ai n
Y i=1 P Sji ∈ Ai , Nti −ti−1 = ji P Nti − Nti−1 = ki − ki−1


P Xti −ti−1 ∈ Ai . Si utilizamos el caso n = 2, con A1 = R, podemos concluir que Xt2 − Xt1 y Xt2 −t1 tienen la misma distribución, de donde la igualdad n

Y P Xti − Xti−1 ∈ Ai , i = 1, . . . , n = P Xti −ti−1 ∈ Ai i=1 nos permite concluir la independencia de los incrementos y por lo tanto, que el proceso Poisson compuesto es un proceso de Lévy. La propiedad de Markov es válida más generalmente para procesos de Lévy y tiene una demostración muy similar a la del proceso de Poisson.  Veamos una aplicación adicional del Teorema 4.1. Proposición 4.5. Si N 1 y N 2 son dos procesos de Poisson independientes de intensidad λ1 y λ2 entonces N 1 + N 2 es un proceso de Poisson de intensidad λ1 + λ2 . En efecto, N 1 + N 2 es un proceso de contéo y de Lévy y su primer salto, que es igual a S11 ∧ S21 , tiene distribución exponencial de parámetro λ1 + λ2 . CAPÍTULO 5 Procesos de Markov constantes por pedazos En este capı́tulo analizaremos una clase de procesos estocásticos que generaliza al proceso de Poisson y tiene similitudes con las cadenas de Markov: los procesos de Markov a tiempo continuo con valores en un conjunto a lo más numerable E. Comenzaremos con un estudio del proceso de Poisson, al expresar su propiedad de Markov de forma que se parezca a la de las cadenas de Markov. Luego, se introducirá un segundo ejemplo, el de los procesos de nacimiento puro. Finalmente, comenzaremos el estudio de procesos con trayectorias constantes por pedazos y daremos una descripción probabilı́stica de estos al introducir parámetros que determinan a un proceso de Markov: la distribución inicial y la matriz infinitesimal. Luego, estudiaremos sus probabilidades de transición mediante unas ecuaciones diferenciales que satisfacen y que están ligadas con su matriz infinitesimal: las ecuaciones backward y forward de Kolmogorov. Finalmente, veremos cómo las ecuaciones backward nos permiten estudiar a las distribuciones invariantes de los procesos de Markov a tiempo continuo. 1. El proceso de Poisson como proceso de Markov El proceso de Poisson toma valores en los naturales. Como lo hemos definido, siempre comienza en cero, a diferencia de las cadenas de Markov que podemos comenzar en cualquier parte de su espacio de estados. Recordemos que si N es un proceso de Poisson de parámetro λ, al tener incrementos independientes y estacionarios, podemos escribir, para = t0 < t1 < t2 < · · · < tn P(Nt1 = k1 , Nt2 = k2 , . . . , Ntn = kn )
= P(Nt1 = k1 ) P(Nt2 −t2 = k2 − k1 ) · · · P Ntn − Ntn−1 = kn − kn−1 . Por lo tanto, si definimos n Pt (i, j) = P(Nt = j − i) = e−λt (λt) , n! vemos que P(Nt1 = k1 , Nt2 = k2 , . . . , Ntn = kn ) = Pt1 (0, k1 ) Pt2 −t1 (k1 , k2 ) · Ptn (kn−1 , kn ) . La expresión anterior ya es muy paralela a la que vemos en cadenas de Markov; la entrada i, j de la matriz (infinita) Pt se puede interpretar como la probabilidad de 74 1. El proceso de Poisson como proceso de Markov 75 ir de i a j en t unidades de tiempo. Esto además nos dice cómo podrı́amos definir a un proceso estocástico que fuera como el proceso de Poisson pero que comenzara en k: simplemente substituimos Pt1 (0, k1 ) por Pt1 (k, k1 ). Sin embargo, no tenemos un proceso estocástico que satisfaga lo anterior; demos una construcción de él. Sea N un proceso de Poisson de parámetro λ. Si t > 0, hemos visto que el proceso N t dado por Nst = Nt+s − Nt es un proceso de Poisson independiente de Nt . Por lo tanto: P(Nt+s1 = k1 , . . . , Nt+sn = kn |Nt = k) P(Nt+s1 = k1 , . . . , Nt+sn = kn , Nt = k) P(Nt = k) P(Ns1 = k1 − k, . . . , Nsn = kn − k) P(Nt = k) = P(Nt = k) = P(Ns1 + k = k1 , . . . , Nsn + k = kn ) . = Vemos entonces que, la distribución condicional del proceso de Poisson en tiempos posteriores a t condicionalmente a Nt = k es la misma que la del proceso k + N . Ası́, es natural definir al proceso de Poisson que comienza en k como k + N , pues este proceso tiene trayectorias como las del proceso de Poisson (aumenta de uno en uno en ciertos tiempos aleatorios) pero comienza en k. Además: P(Nt+s1 = k1 , . . . , Nt+sn = kn |Nt = k) = Pt1 (k, k1 ) Pt2 −t1 (k1 , k2 ) · Ptn (kn−1 , kn ) Recordemos que en el caso de cadenas de Markov, se satisfacen las ecuaciones de Chapman-Kolmogorov. En este contexto, las ecuaciones se pueden expresar como sigue: X Pt+s (i, k) = Ps (i, j) Pt (j, k) . j Ejercicio 5.1. Al utilizar el teorema del biniomio, pruebe directamente que la ecuación anterior se satisface. Dé además un argumento probabilı́stico, basado en condicionar con lo que sucede al tiempo s, para probar dicha ecuación. La diferencia con las cadenas de Markov y el proceso de Poisson es que, al ser el segundo un proceso de Markov a tiempo continuo, en lugar de tener una sola matriz cuyas potencias nos permiten describir al proceso, tenemos toda una colección de matrices (Pt , t ≥ 0). Uno de los objetivos de este capı́tulo es mostrar como, con una definición adecuada, podemos generar a todas las matrices (Pt , t ≥ 0) mediante una sola matriz, la llamada matriz infinitesimal o matriz de tasas de transición. La idea es la siguiente: si x es un número, podemos interpolar a la sucesión xn , n ∈ N mediante la función exponencial: xn = en log x . El lado derecho ya tiene sentido si n ≥ 0 y no sólo si n ∈ R. En cierto sentido, si Pt = etQ , podrı́amos interpretar a Q como el logaritmo de P1 . Para ver cómo podrı́amos obtener a Q, notemos que para 1. El proceso de Poisson como proceso de Markov 76 obtener a log x si conocemos a f (t) = et log x , simplemente calculamos log x = f ′ (0) y además f ′ (t) = f ′ (0) f (t). En el caso de nuestra matriz Pt , vemos que j−i d (λt) d Pt (i, j) = e−λt dt dt (j − i)! h i 1 j−i = e−λt (j − i) tj−i−1 λj−i 1j≥i+1 − λe−λt (λt) (j − i)! por lo que al evaluar en cero se obtiene:   −λ j = i d Pt (i, j) = λ j =i+1 .  dt t=0  0 j 6= i, i + 1 Ejercicio 5.2. Sea   −λ Q(i, j) = λ   0 j=i j =i+1 . j 6= i, i + 1 Pruebe directamente que d (5) Pt (i, j) = QPt (i, j) = Pt Q(i, j) , dt donde QPt es el producto de las matrices Q y Pt . A las ecuaciones (diferenciales) (5) se les conoce como ecuaciones de Kolmogorov. Esbocemos una interpretación y deducción probabilı́stica de dichas ecuaciones. Para calcular Pt+h (0, j), podemos descomponer respecto del valor de Nh para obtener Pt+h (i, j) = P(Nt+h = j) = P(Nh = i) P(Nt = j − i) . Por otra parte, al utilizar explı́citamente la distribución Poisson, vemos que P(Nh = 1) P(Nh ≥ 2) 1 − P(Nh = 0) = λ lim =λ lim = 0. lim h→0 h→0 h→0 h h h La interpretación es que es muy probable que haya cero saltos en un intervalo de longitud h, hay probabilidad proporcional al tamaño del intervalo de que haya un sólo salto, y muy improbable que haya más de dos saltos. Se concluye que Ph (i, i) Pt (i, k) − Pt (i, k) ∂Pt (i, k) = lim h→0 ∂t h Ph (i, i + 1) Pt (i + 1, k) + lim h→0 h = Pt (i + 1, k) λ − Pt (i, j) λ = Pt Q(i, k) A esta ecuación se le denomina ecuación hacia atrás de Kolmogorov (también llamada, aún en español, ecuación backward) y fué obtenida al descomponer a 2. El proceso de nacimiento puro 77 Nt+h respecto del valor de Nh . Ésta implica automáticamente la ecuación hacia adelante de Kolmogorov ∂Pt (i, k) = QPt (i, k) ∂t al utilizar la relación Pt (i + 1, j) = Pt (i, j − 1). 2. El proceso de nacimiento puro El proceso de nacimiento puro es una generalización del proceso de Poisson que puede presentar un fenómeno interesante: la explosión. Daremos una definición de dicho proceso paralela a la del proceso de Poisson. Definición. Sean q0 , q1 , . . . ∈ (0, ∞). Un proceso de nacimiento puro es el proceso de contéo cuyos tiempos interarribo conforman una sucesión de variables exponenciales independientes S1 , S2 , . . . de parámetros λ0 , λ1 , . . .. La interpretación del proceso de nacimiento puro es que λi representa la tasa a la que un nuevo individuo se agrega a una población cuando esta tiene i elementos. Ası́, el proceso de nacimiento puro se construye al definir a los tiempos de salto ∞ X Si 0 = T 0 T n = S 1 + · · · + Sn T ∞ = n=0 y al proceso de contéo asociado Nt = X n1Tn ≤t<Tn+1 . n∈N Claramente, cuando qi = λ para toda i el proceso que hemos construido es el proceso de Poisson. Sin embargo, cuando tenemos qi dependientes de i se presenta un nuevo fenómeno. En efecto, en el caso del proceso de Poisson, sabemos que Tn /n → 1/λ = E(S1 ) > 0 por la ley fuerte de los grandes números y por lo tanto T∞ = ∞ casi seguramente. Sin embargo, en el caso del proceso de nacimiento puro puede suceder que T∞ < ∞, en cuyo caso: lim Nt = ∞. t→T∞ − Decimos entonces que ha ocurrido la explosión (en este caso demográfica) del proceso de nacimiento puro y es natural definir X Nt = ∞1T∞ ≤t + n1Tn ≤t<Tn+1 n∈N pues con la definición anterior el proceso se vuelve cero despues de que la población haya explotado. Puesto que el proceso de nacimiento puro admite una construcción sencilla, se puede dar un criterio explı́cito para la explosión, el cual además prueba que la explosión se da o con probabilidad cero o con probabilidad uno. 2. El proceso de nacimiento puro 78 P P Proposición 5.1. T∞ = ∞ casi seguramente si E(T∞ ) = i 1/qi = ∞. Si i 1/qi < ∞ entonces T∞ < ∞ casi seguramente. P Demostración. Si i 1/qi < ∞ entonces podemos aplicar el teorema de convergencia monótona a la sucesión Tn para deducir que X X E(T∞ ) = lim E(Tn ) = lim 1/qm = 1/qn < ∞, n n n m≤n por lo que T∞ < ∞ casi Pseguramente. Por otra parte, si i 1/qi = ∞, entonces podemos aplicar el teorema de convergencia acotada para ver que n Y

E e−T∞ = lim E e−Tn = lim n n qi . 1 + qi i=1 Ahora dividiremos en dos partes nuestro análisis: si existe µ > 0 tal que qi ≤ µ para una cantidad infinita de ı́ndices i entonces n  n Y qi 1 →0 ≤ 1 + qi 1 + 1/µ i=1 conforme n → ∞. Si por otra parte qi → ∞ entonces utilizamos el hecho de que qi log (1 + 1/qi ) → 1 conforme n → ∞ y por lo tanto n Y Pn Pn qi = e− i=1 log(1+1/qi ) ≤ e−C i=1 1/qi → 0. 1 + qi i=1 Vemos que en cualquier caso
E e−T∞ = 0, lo cual nos dice que T∞ = ∞ casi seguramente.  Para analizar la propiedad de Markov del proceso de nacimiento puro, necesitamos definirlo cuando comienza en i ∈ N. Una definición posible es que se trata del proceso (NTi +t , t ≥ 0). Este proceso va tomando los valores sucesivos i, i + 1, i + 2, . . . y los tiempos que permanece en cada uno de los estados son Si , Si+1 , . . . ,. Pues este proceso tiene la misma estructura probabilı́stica que i más un proceso de contéo cuyos tiempos interarribo son exponenciales independientes de parámetros λi , λi+1 , . . ., o sea, que una definición equivalente es que un proceso de nacimiento puro con parámetros λ0 , λ1 , . . . que comienza en i es i más un proceso de nacimiento puro con parámetros λi , λi+1 , . . .. Con estos preliminares podemos enunciar la propiedad de Markov del proceso de nacimiento puro. Proposición 5.2. Sean N un proceso de nacimiento puro de parámetros λ0 , λ1 , . . . y s ≥ 0. Condicionalmente a Ns = i, el proceso estocástico (Nt+s − i, t ≥ 0) es un proceso de nacimiento puro de parámetros λi , λi+1 , . . . y es condicionalmente independiente de Fs = σ(Nr : r ≤ t) dado Ns = i. 2. El proceso de nacimiento puro 79 En otras palabras, condicionalmente a Ns = i, el proceso Nt+s , t ≥ 0 es un proceso de nacimiento puro que comienza en i. Demostración. El proceso Nt+s −i, t ≥ 0 es un proceso de contéo cuyos tiempos interarribo son TNs +1 − s, SNs +2 , SNs +3 , . . .. Denotémoslos como S0s , S1s , . . .. Debemos ver que condicionalmente a Ns = i, estos tiempos interarribo son exponenciales, independientes y de parámetros λi , λi+1 , . . .. En efecto: P(Ns = i, S0s > s0 , · · · , Sns > sn ) = P(Ti ≤ s < Ti + Si , s0 + s < Ti + Si , s1 < Si+1 , . . . , sn < Si+n ) . Al condicionar por Ti , utilizar la independencia de las variables Ti, Si , . . . , Si+n y la propiedad de pérdida de memoria de la distribución exponencial vemos que P(Ti ≤ s < Ti + Si , s0 + s < Ti + Si , s1 < Si+1 , . . . , sn < Si+n ) .
= E 1Ti ≤s e−λi (s − Ti ) e−λi s1 · · · e−λi+n sn . Sin embargo, al poner s0 = · · · = sn = 0, vemos que
E 1Ti ≤s e−λi (s − Ti ) = P(Ns = i) y por lo tanto, al condicional por Ns = i, vemos que N s es un proceso de nacimiento puro que comienza en i. Falta demostrar que N s es condicionalmente independiente de Fs dado que Ns = i, esto es, que para todo A ∈ Fs , se tiene que
P Nts1 = k1 , . . . , Ntsn = kn , A |Ns = i
= P Nts1 = k1 , . . . , Ntsn = kn |Ns = i P(A |Ns = i) . Por clases monótonas, basta verificarlo cuando A = {Ns1 = j1 , . . . , Nsm = jm } con s1 ≤ · · · ≤ sm ≤ t y j1 ≤ · · · ≤ jm ≤ i. Para esto, seguimos el mismo razonamiento anterior al notar que A ∩ {Ns = i} = {Tj1 ≤ s1 < Tj1 +1 , · · · , Tjm ≤ s < Tjm +1 } ∩ {Ti ≤ s < Ti+1 } .  Ejercicio 5.3. Haga un programa en R que simule al proceso de nacimiento puro que comienza en 1 si λi = iλ para algún λ > 0. ¿Este proceso explota? Haga otro programa que simule el caso λi = λi2 y diga si ocurre explosión o no. La propiedad de Markov se puede interpretar de manera similar a la del proceso de Poisson. Dados los parámetros q = (λ0 , λ1 , . . .), sea Pt (i, j) la probabilidad de que un proceso de nacimiento puro de parámetro q que comienza en i se encuentre en j al tiempo t. Entonces: P(Ns+t1 = k1 , . . . , Ns+tn = kn |Ns = k) = Pt1 (k, k1 ) Pt2 −t1 (k1 , k2 ) · · · Ptn −tn−1 (kn−1 , kn ) . 3. Matrices infinitesimales y construcción de procesos de Markov 80 Aquı́ es más difı́cil obtener una expresión explı́cita para Pt (i, j). Ésto se puede lograr cuando por ejemplo λi = iλ para alguna λ > 0 (que es el caso del Proceso de Yule, básico en la teorı́a de procesos de ramificación a tiempo continuo). Sin embargo, podemos argumentar por qué son válidas las ecuaciones hacia atrás de Kolmogorov: Si h es pequeño y nuestro proceso de nacimiento puro comienza en i, será raro que haya más de dos saltos en un intervalo de tiempo pequeño. Lo que se afirma es que, exactamente como en el caso de procesos de Poisson, P Ph (i, i + 1) 1 − Ph (i, i) j≥2 Ph (i, j) = λi lim = λi lim = 0. lim h→0 h→0 h→0 h h h Por lo tanto, las ecuaciones hacia atrás de Kolmogorov serı́an ∂ Pt (i, j) = λi Pt (i + 1, j) − λi Pt (i, j) . ∂t Vemos entonces que la matriz infinitesimal o de tasas de transición estarı́a dada por   −λi i = j Q(i, j) = λi j =i+1   0 en otro caso y que las ecuaciones hacia atrás se pueden escribir como ∂ Pt (i, j) = QPt (i, j) . ∂t 3. Matrices infinitesimales y construcción de procesos de Markov Comenzaremos por describir a un tercer ejemplo de proceso de Markov a tiempo continuo. Este ejemplo se basa en calcular el mı́nimo de variables exponenciales independientes. Comenzaremos con un ejemplo con espacio de estados E finito. Para cada x, y ∈ E con x 6= y sea λx,y ≥ 0. Sean Ei,x,y , ≥ 1, x, y ∈ E, x 6= y variables aleatorias exponenciales independientes, tal queEi,x,y tiene parámetro λx,y . A λx,y lo interpretaremos como la tasa a la que dejo el estado x para acceder al estado y y supondremos que para cada x ∈ E existe y ∈ E distinta de x tal que λx,y > 0. Si x ∈ E, definamos a T1 = min E1,x,y y6=x donde T1 = E1,x,X1 . Recursivamente, definiremos a Tn+1 = min En+1,Xn ,y y6=Xn donde Tn+1 = E1,Xn ,Xn+1 . Definamos finalmente, al proceso a tiempo continuo de interés: Zt = Xn si t ∈ [Tn , Tn+1 ). 3. Matrices infinitesimales y construcción de procesos de Markov 81 Hagamos un ejemplo en R para ilustrar la definición. Comencemos con E = {1, 2}, por lo que necesitamos definir a dos cantidades positivas λ1,2 (que denotaremos por λ1 y digamos que es igual a 2) y λ2,1 (denotada λ2 y digamos igual a 1). Si queremos comenzar nuestra construcción con x = 1, podemos utilizar el siguiente código. # Simulemos una cadena de Markov en tiempo continuo en $ E ={1 ,2} $ # l (1) es la tasa a la que dejo $ 1 $ ( para ir a 2) y viceversa l = c (2 ,1) ; # x es el estado actual , comenzamos en x =1; X=c(x); # T contendr \ ’ a los tiempos en que cambia la cadena T =0 # n es la cantidad de pasos n =20; for ( i in 1: n ) { # Cambio de estado x =3 - x ; # Agrego el nuevo estado X = c (X , x ) # Simulo la exponencial de la tasa adecuada y la agrego a mis tiempos T = c (T , tail (T ,1) + rexp (1 , l [ x ]) ) ; } # Grafico la trayectoria plot (T ,X , type = " s " ) Listing 5.1. QCadena2Estados.R Una exponencial de parámetro λ tiene la misma distribución que una exponencial de parámetro 1 dividida entre λ, por lo cual en el ciclo es equivalente utilizar el comando T=c(T,tail(T,1)+rexp(1)/l[x]);. Un ejemplo un poco más complicado lo obtenemos si E = {1, 2, 3} y ponemos las tasas de transición λ1,2 = 2, λ2,1 = λ2,3 = 1, λ3,1 = 1/3 y todas las demás igual a cero. En este caso, nuestro primer impulso podrı́a ser el utilizar el siguiente código: # Simulemos una cadena de Markov en tiempo continuo en $ E ={1 ,2 ,3} $ # La matriz de tasas de transici \ ’ on L = matrix (0 ,3 ,3) L [1 ,2]=2 L [2 ,1]=1 L [2 ,3]=1 L [3 ,1]=1 / 3 # Estado inicial x =1 # Vector de estados X=c(x) # Vector de tiempos T = c (0) # N \ ’ umero de pasos n =20 3. Matrices infinitesimales y construcción de procesos de Markov 82 for ( i in 1:20) { # Genero las exponenciales e = rexp (3) / L [x ,] # El nuevo tiempo de salto t = min ( e ) T = c (T , tail (T ,1) + t ) # El nuevo estado x = which . min ( e ) X = c (X , x ) } plot (T ,X , type = " s " ) Listing 5.2. QCadena3Estados.R Hay una manera de realizar la simulación anterior con menos variables aleatorias. En vez de las 3 exponenciales, podemos utilizar una exponencial y una uniforme. La idea es utilizar la proposición siguiente para ver que el mı́nimo de variables exponenciales es exponencial. Proposición 5.3. Sean T1 , . . . , Tn variables exponenciales independientes de parámetros respectivos λ1 , . . . , λn . Sea T = mini Ti . Entonces, con probabilidad 1, existe un único ı́ndice K tal que T = TK . T y K son independientes, T tiene distribución exponencial de parámetro λ = λ1 + · · · + λn y P(K = k) = λk /λ. Demostración. Calculemos: P(T ≥ t, Tj > Tk para toda j 6= k) = P(Tk ≥ t, Tj > Tk para toda j 6= k)   Y e−λk Tk  = E1Tk ≥t j6=k = Z ∞ λk e−λk s t λk −λt = e . λ Y e−λj s ds j6=k Al evaluar en t = 0 vemos que λk λ por lo que al utilizar que los eventos anteriores son ajenos conforme variamos a k y su unión es la probabilidad de que el mı́nimo de T1 , . . . , Tn se alcance para un sólo ı́ndice, vemos que X λk λ = = 1. P(Existe un único k ∈ {1, . . . , n} tal que T = Tk ) = λ λ P(Tj > Tk para toda j 6= k) = k Al sumar sobre k, vemos que P(T ≥ t) = e−λt , 3. Matrices infinitesimales y construcción de procesos de Markov 83 lo cual implica que T es exponencial de parámetro λ Puesto que con probabilidad 1 hay un único ı́ndice (aleatorio) en el que se alcanza el mı́nimo, digamos K tal que TK = T , vemos que P(K = k, T ≥ t) = λk −λt e = P(K = k) P(T ≥ t) . λ  Ası́, nuestro proceso a tiempo continuo se puede simular también mediante el siguiente código. # Simulemos una cadena de Markov en tiempo continuo en $ E ={1 ,2 ,3} $ # La matriz de tasas de transici \ ’ on L = matrix (0 ,3 ,3) L [1 ,2]=2 L [2 ,1]=1 L [2 ,3]=1 L [3 ,1]=1 / 3 # Tasas totales de salto l = rowSums ( L ) # Matriz de transici \ ’ on P = L / rowSums ( L ) # Estado inicial x =1 # Vector de estados X=c(x) # Vector de tiempos T = c (0) # N \ ’ umero de pasos n =20 for ( i in 1:20) { # Genero la exponencial e = rexp (1 , l [ x ]) # El nuevo tiempo de salto T = c (T , tail (T ,1) + e ) # El nuevo estado x = sample (3 ,1 , prob = P [x ,]) X = c (X , x ) } plot (T ,X , type = " s " ) Listing 5.3. QCadena3EstadosConMinimo.R Continuemos con el análisis de la cadena general. A partir del código, vemos que los parámetros se pueden reducir a una matriz de transición P con ceros en la diagonal y un vector l de tasas totales de salto. Explı́citamente, X λx,y . l(x) = λx,y y Px,y = l(x) y Además, es posible construir nuestra cadena mediante una sucesión de variables exponenciales estándar independientes ξ1 , ξ2 , . . . entre sı́ y de una sucesión iid U1 , U2 , . . . de variables uniformes. En efecto, podemos definir a X0 = x, T0 = x 3. Matrices infinitesimales y construcción de procesos de Markov 84 y a T1 = S1 = ξ1 /l(X0 ). Luego, utilizamos a la variable uniforme U1 para escoger a un elemento aleatorio X1 de E tal que P(X1 = y) = PX0 ,y . Finalmente, continuamos este procedimiento de manera recursiva: si ya tenemos definidos a X0 , . . . , Xn y a T0 , . . . , Tn en términos de S1 , . . . , Sn y U1 , . . . , Un entonces definimos Sn+1 = ξn+1 /l(Xn ), Tn+1 = Tn + Sn+1 y utilizamos a las variables Un+1 y Xn para construir a Xn+1 de tal manera que P(Xn+1 = y |Xn = x) = Px,y . Finalmente, recordemos que Zt = Xn si Tn ≤ t < Tn+1 . Resulta ser que Z es un proceso de Markov a tiempo continuo. Formalmente, sea Pt (x, y) la probabilidad de que Zt = y cuando Z0 = x y probemos que Px0 (Xt1 = x1 , . . . , Xtn = yn ) = Pt1 −t0 (x0 , x1 ) Pt2 −t1 (x1 , x2 ) · · · Ptn −tn−1 (xn−1 , xn ) donde Px es la medida de probabilidad que rige a Z cuando comienza en x y 0 = t0 < t1 < · · · < tn . El argumento es parecido a cuando probamos que el proceso de Poisson (o más generalmente el de nacimiento y muerte) es un proceso de Markov, se basa en notar que el proceso posterior a s, condicionalmente a Zs = x tiene la misma construcción probabilı́stica que Z comenzando en x. Esto es, que sus tiempos y lugares de salto se pueden obtener a partir de una sucesión de variables uniformes y exponenciales estándar independientes. Concentrémonos mejor en las aplicaciones de este hecho. La primera aplicación es a la obtención de las ecuaciones de Chapman-Kolmogorov. En efecto, vemos que Pt+s (x, z) = Px (Zt+s = z) X = Px (Zs = y, Zt+s = z) y = X Ps (x, y) Pt (y, z) . y Seguido de esto, obtendremos las ecuaciones de Kolmogorov como en el caso de procesos de procesos de Poisson. Se afirma que lim 1 − Px (Nh = 0) = l(x) , h lim Px (Nh ≥ 2) = 0. h h→0 lim h→0 Px (Nh = 1) , Zt = y = l(x) Px,y h y h→0 En efecto, puesto que bajo Px el primer salto de N ocurre en una variable exponencial de parámetro l(x), entonces 1 − e−l(x)h 1 − Px (Nh = 0) = lim = l(x) . h→0 h→0 h h lim 3. Matrices infinitesimales y construcción de procesos de Markov 85 Por otra parte, vemos que {Nh ≥ 2} = {S1 + S2 ≤ h}   ξ2 ξ1 + ≤h = l(X0 ) l(X1 )   ξ1 ξ2 ⊂ + ≤h . λ λ El último evento corresponde a que un proceso de Poisson de parámetro λ, digamos N λ tenga más de dos saltos en el intervalo [0, t], que ya hemos estimado. Se concluye que
Px Nhλ ≥ 2 Px (Nh ≥ 2) 0 ≤ lim sup ≤ lim = 0. h→0 h h h→0 Finalmente, veamos cómo ocuparnos de Px (Nh , Xh = y) = 1. Puesto que U1 es independiente de ξ1 , ξ2 , entonces   ξ1 1 ξ1 ξ2 1 Px,y Px (Nh = 1, Zh = y) = Px ≤h< + h h l(x) l(x) l(y)   ξ1 ξ1 ξ2 1 ≤h< + = Px,y Px h l(x) l(x) λ Z h 1 l(x) e−l(x)s e−l(y)(h−s) ds = Px,y h 0 Z 1 h −s[l(y)−l(x)] e ds = Px,y l(x) e−l(y)h h 0 →h→0 l(x) Px,y . Podemos entonces deducir que lim h→0 Ph (x, x) − 1 = −l(x) h y para y 6= x lim h→0 Ph (x, y) = Px,y l(x) . h Llámemósle Q a la matriz Qx,y = ( −l(x) l(x) Px,y y=x . y 6= x A esta matriz se le conoce como matriz infinitesimal del proceso de Markov Z. 4. Descripción probabilı́stica de procesos de Markov constantes por pedazos 86 Al utilizar la finitud de E para intercambiar derivadas y sumas y las ecuaciones de Chapman-Kolmogorov se obtiene Pt+h (x, z) − Pt (x, z) d Pt (x, z) = lim h→0 dt h Ph (z, z) − 1 X 1 + Pt (x, y) lim Ph (y, z) = lim Pt (x, z) h→0 h h→0 h y6=z X = Pt (x, y) l(y) Py,z − Pt (x, y) l(y) y6=z = (Pt Q)x,y . Estas son las ecuaciones forward de Kolmogorov. De igual manera, al utilizar la finitud de E se obtienen las ecuaciones backward de Kolmogorov d Pt (x, z) = (QPt )x,y . dt Como ya se ha comentado, la anterior ecuación diferencial matricial es parecida a la ecuación diferencial f ′ (t) = λf (t) (con f (t) = Pt y λ = Q) cuya única solución es f (t) = eλt . Aún en el caso matricial se le puede dar sentido a esta ecuación diferencial al definir la exponencial de la matriz tQ por medio de X tn Q n etQ = . n! n Entonces, se deduce que Pt = etQ . La pregunta que sigue es si es posible realizar la construcción que acabamos de explorar en espacio de estados numerable. La respuesta es básicamente que sı́, salvo que se debe prestar atención al fenómeno de explosión. Otra pregunta importante es si cualquier cadena de Markov a tiempo continuo se puede construir como lo acabamos de hacer. La respuesta es básicamente afirmativa de nuevo cuidandonos de la explosión. La siguiente sección aborda este segundo cuestionamiento. 4. Descripción probabilı́stica de procesos de Markov constantes por pedazos En esta sección definiremos a las cadenas de Markov a tiempo continuo y analizaremos su estructura y su comportamiento a tiempos grandes. Se realizará por lo tanto un estudio paralelo al de las cadenas de Markov a tiempo discreto. Sea E un conjunto a lo más numerable, ∆ 6∈ E algún punto que denominaremos cementerio y a donde mandaremos las trayectorias de un proceso de Markov cuando explote, y sea C el conjunto de funciones f : [0, ∞) → E ∪ {∆} que satisfacen: 4. Descripción probabilı́stica de procesos de Markov constantes por pedazos 87 (1) si f (t) ∈ E entonces existe δ > 0 tal que f (s) = f (t) para s ∈ [t, t + δ], (2) si f (t) = ∆ entonces f (s) = ∆ para toda s ≥ t y (3) si t1 < t2 < · · · , tn → t < ∞ y f (tn+1 ) 6= f (tn ) entonces f (t) = ∆. Definición. Una cadena de Markov a tiempo continuo con espacio de estados E es un proceso estocástico X = (Xt , t ≥ 0) tal que (1) Xt toma valores en E ∪ {∆} (2) para todo ω ∈ Ω, la trayectoria t 7→ Xt (ω) es un elemento del conjunto de funciones C, (3) existe una colección de matrices estocásticas (Pt , t ≥ 0) indexadas por los elementos de E tales que si 0 = t0 < t1 < · · · < tn entonces: P(X0 = x, Xt1 = x1 , . . . , Xtn = xn ) = P(X0 = x) Pt1 −t0 (x0 , x1 ) Pt2 −t1 (x1 , x2 ) · · · Ptn −tn−1 (xn−1 , xn ) . Denotaremos por Px a P condicionada por X0 = x. A la colección de matrices estocásticas (Pt , t ≥ 0) les llamamos probabilidades de transición. Comenzamos por verificar una propiedad familiar de las probabilidades de transición. Proposición 5.4. Las probabilidades de transición satisfacen las ecuaciones de Chapman-Kolmogorov Pt+s = Pt Ps . Demostración. Al descomponer Px (Xt+s = z) por el valor que toma Xt se obtiene X X Px (Xt+s = z) = Px (Xt = y, Xt+s = z) = Pt (x, y) Ps (x, y) , y y de lo cual se sigue que Pt+s es el producto de las matrices Pt y Ps .  Proposición 5.5 (Propiedad de Markov). Si X es una cadena de Markov a tiempo continuo y Xst = Xt+s , entonces X t también es una cadena de Markov a tiempo continuo con las mismas probabilidades de trancisión que X. X t es independiente de Xs , s ≤ t condicionalmente a Xt . Demostración. Sean t1 < t2 < · · · . Entonces Px X0t = x0 , Xtt1 = x1 , . . . , Xttn = xn
= Px (Xt+0 = x0 , Xt+t1 = x1 , . . . , Xt+tn = xn ) Pt (x, x0 ) Pt1 (x0 .x1 ) · · · Ptn −tn−1 xn−1 , xn . Puesto que Pt (x, x0 ) = Px (X0t = x0 ), vemos que X t es una cadena de Markov a tiempo continuo con las mismas probabilidades de transición que X. 4. Descripción probabilı́stica de procesos de Markov constantes por pedazos 88 Para ver la independencia, hacemos un cálculo similar: sean 0 < s1 < · · · < sm ≤ t y 0 < t1 < · · · < tn . Al utilizar la definición de cadena de Markov, obtenemos
Px Xs1 = x1 , . . . , Xsn = xn , X0t = y0 , Xtt1 = y1 , . . . , Xttn = yn = Ps1 (x, x1 ) Ps2 −s1 (x1 , x2 ) · · · Psn −sn−1 (xn−1 , xn ) Pt−sn (xn , y0 ) Pt1 (x1 , y0 ) Pt2 −t1 (y1 , y2 ) · · · Ptn −tn−1 (yn−1 , yn ) . Al condicional por Xt = y0 vemos que Px Xs1 = x1 , . . . , Xsn = xn , X0t = y0 , Xtt1 = y1 , . . . , Xttn = yn Xt = y0
Px Xs1 = x1 , . . . , Xsn = xn X0t = y0 Py0 (Xt1 = y1 , . . . , Xtn = yn ) por lo que X t es independiente de Xs , s ≤ t condicionalmente a Xt .
 Consideremos a los tiempos aleatorios T0 = 0, Tn+1 = inf {t ≥ Tn : Xt 6= XTn } y ζ = lim Tn n→∞ con la convención inf ∅ = ∞. Proposición 5.6. Tn y ζ son tiempos de paro respecto de la filtración canónica. Existen tres categorı́as para las trayectorias en términos de estos tiempos aleatorios: Absorción: Cuando existe n tal que Tn < ∞ = Tn+1 , en cuyo caso Xt = XTn para toda t ≥ Tn , Explosión: cuando ζ < ∞ y Movimiento perpetuo: cuando Tn < ∞ para toda n y ζ = ∞. Proposición 5.7. Bajo Px , T1 es exponencial de parámetro c(x) ∈ [0, ∞). Si c(x) > 0 entonces las variables XT1 y T1 son independientes. Demostración. Al utilizar la propiedad de Markov, vemos que Px (T1 > t + s) = Px (T1 > s, Xs = x, T1 (X s ) > t) = Px (T1 > s) Px (T1 > t) y por lo tanto, bajo Px , T1 es exponencial. Por otra parte Px (1T1 >t XT1 = y)   = Px 1T1 >t , Xt = x, XTt 1 (X t ) = y = Px (1T1 >t ) Px (XT1 = y) por lo que T1 es independiente de X ◦ θT1 .  4. Descripción probabilı́stica de procesos de Markov constantes por pedazos 89 A c(x) la interpretamos como la tasa a la que dejamos el estado x. Definamos ahora ( 0 c(x) = 0 Px,y = y α(x, y) = c(x) Px,y . Px (XT1 = y) c(x) 6= 0 A α se le conoce como la matriz de tasas de transición y la interpretación de α(x, y) es la tasa a la que dejamos el estado x para pasar al estado y. La matriz de tasas de transición es el parámetro que nos permitirá caracterizar a la familia Markoviana. Para verificar por qué, es necesario extender la propiedad de Markov. Teorema 5.1 (Propiedad de Markov fuerte). Sea T un tiempo de paro finito. Entonces el proceso X T dado por XtT = XT +t es una cadena de Markov con las mismas probabilidades de transición que X que es independiente de Xs , s ≤ t condicionalmente a XT y a T ≤ t. La prueba es un tanto más complicada del nivel que se quiere para este curso y no se dará. El teorema anterior nos permitirá caracterizar a la cadena de Markov en tiempo continuo en términos de la matriz de tasas de transición α, o equivalentemente, de c y P . Sea Z el proceso estocástico a tiempo discreto definido por Z n = XT n si Tn < ∞. En el caso absorbente, definimos Zn+m = Zn para toda m ≥ 1 si Tn < ∞ = Tn+1 . Teorema 5.2. El proceso Z es una cadena de Markov de matriz de transición P que comienza en x bajo Px . Si c(x) > 0 para toda x ∈ E, condicionalmente a Z, las variables S1 , S2 , . . . con Si = Ti − Ti−1 son independientes y exponenciales de parámetros c(Z0 ) , c(Z1 ) , . . .. Demostración. Al utilizar el lema de clases de Dynkin, vemos que es suficiente verificar que (6) Px (Z1 = x1 , . . . , Zn = xn , S1 > t1 , . . . , Sn > tn ) = Px,x1 · · · Pxn−1 ,xn e−c(x)t1 e−c(x1 )t2 · · · e−c(xn−1 )tn . Esto se sigue por inducción al utilizar la propiedad de Markov fuerte. La base inductiva es la Proposición 5.7. Por otra parte, si suponemos válida la ecuación (6) vemos que al aplicar la propiedad de Markov fuerte al instante Tn y la Proposición 5.7 se sigue que Px (Z1 = x1 , . . . , Zn+1 = xn+1 , S1 > t1 , . . . , Sn+1 > tn+1 ) = Px (Z1 = x1 , . . . , Zn = xn , S1 > t1 , . . . , Sn > tn ) e−c(xn )tn+1 Pxn ,xn+1 puesto que Zn+1 es el primer estado al que salta X Tn y Sn+1 es el tiempo que tarda en realizar X Tn su primer salto.  5. Las ecuaciones backward y forward de Kolmogorov 90 P Dada una función α : E × E → [0, ∞) tal que c(x) = y α(x, y) < ∞, podemos definir a Px,y = α(x, y) /c(x) cuando c(x) > 0 y a Px,y = δx,y cuando c(x) = 0 y preguntarnos cuándo existe una cadena de Markov a tiempo continuo cuya matriz de tasas de transición sea α. Nos abocaremos ahora a verificar que se puede construir una cadena de Markov cuya matriz de tasas de transición sea α. En efecto, sean S̃1 , S̃2 , . . . variables aleatorias exponenciales de parámetro 1 y Z una cadena de Markov con matriz de transición P que comienza en X. Ahora definamos T0 = 0, y Tn+1 = Tn + S̃n /c(Zn ) . Consideremos al proceso X definido mediante X̃t = Zn si t ∈ [Tn , Tn+1 ). Definimos a Px como P condicionada por Z0 = x. Se afirma que Px es una familia Markoviana cuya matriz de tasas de transición es α. Por definición, bajo la medida de probabilidad Px es válida la ecuación (6). Teorema 5.3. La colección (Px )x∈E es una familia Markoviana con matriz de tasas de transición α. De hecho, hemos verificado este teorema en la sección anterior y la prueba en este caso en el que el espacio de estados es posiblemente infinito es muy similar. 5. Las ecuaciones backward y forward de Kolmogorov Tal como la propiedad de Markov y de Markov fuerte nos llevan a relaciones de recurrencia para probabilidades que deseamos calcular, en tiempo continuo nos llevan a ecuaciones diferenciales. Una de ellas es la ecuación backward de Kolmogorov. Sea (Px , x ∈ E) una familia markoviana con probabilidades de transición Pt (x, y) = Px (Xt = y). Estas probabilidades de transición satisfacen las ecuaciones de Chapman-Kolmogorov X Pt+s (x, z) = Ps (x, y) Pt (y, z) . y Teorema 5.4 (Ecuaciones backward de Kolmogorov). Para cualquier x, y ∈ E, las probabilidades de transición satisfacen la ecuacion backward de Kolmogorov X ∂ α(x, z) Pt (x, z) − Pt (x, y) . Pt (x, y) = ∂t z∈E Dada la matriz de tasas de transición, definiremos a la matriz infinitesimal Q mediante: ( α(x, y) x 6= y Qx,y = . −c(x) x = y 5. Las ecuaciones backward y forward de Kolmogorov 91 Entonces la ecuación backward de Kolmogorov se puede escribir como la ecuación diferencial para la matriz Pt ∂ Pt = QPt . ∂t Esto explica la conexión con ecuaciones diferenciales: las probabilidades de transición de una familia markoviana satisfacen una ecuación diferencial. Veremos que en el caso de espacio de estados finito, la teorı́a clásica de ecuaciones diferenciales lineales nos permite verificar que existe una única solución para la ecuación backward de Kolmogorov y por lo tanto nos da una manera de obtener, a veces explı́citamente pero inclusive también numéricamente, a las probabilidades de transición de la familia Markoviana. Demostración. Heuriı́sticamente, la prueba es una aplicación de la propiedad de Markov fuerte. Sin embargo, necesitamos una verisón que también depende del tiempo. Especı́ficamente, notemos que si s ≤ t Pt (x, y) = Px (Xt = y) = Ex (Pt−s (Xs ) y) . Podemos por lo tanto pensar que por la propiedad de Markov fuerte aplicada al tiempo de paro σ = t ∧ T1 se satisface (7) Pt (x, y) = Ex (Pt−σ (Xσ , y)) para t > 0. Esto es en efecto cierto pero no se ha demostrado y se sigue del hecho de que podemos aproximar al tiempo de paro σ por σ n = ⌈σ2n ⌉/2n y tomar el lı́mite conforme n → ∞ para verificar (7). Veamos cómo se aplica dicha ecuación. De (7) se deduce que: Z tX −c(x)t e−c(x)s α(x, y) Pt−s f (y) ds. Pt f (x) = Ex (Pt−σ f (Xσ )) = f (x) e + 0 y Al multiplicar por ec(x)t de ambos lados se obtiene Z t X ec(x)t Pt f (x) = f (x) + ec(x)s α(x, y) Ps f (y) ds. 0 y Finalmente, la expresión del lado derecho muestra que el lado izquierdo es derivable y por la continuidad del integrando vemos que X ∂ Pt f (x) + c(x) Pt f (x) = αx, yPt f (y) .  ∂t y Ahora recordemos que en espacio de estados finito, digamos de cardinalidad n, estas ecuaciones backward admiten una solución en términos de la matriz infinitesimal Q y esto nos permitirá introducir a las ecuaciones forward. Cuando I es finito, debemos resolver el sistema de ecuaciones ∂ Pt = QPt . ∂t 5. Las ecuaciones backward y forward de Kolmogorov 92 Este es un sistema lineal y si Q fuera de tamaño 1 × 1, tendrı́a como única solución a la función exponencial. Lo mismo sucede en el caso n×n, si definimos a la matriz etQ = ∞ n n X t Q n! n=0 para cualquier t ∈ R.La convergencia de la serie se sigue pues la sucesión N X tn Q n . n! n=0 es de Cauchy cuando se utiliza la norma kQk = max {kQxk : x ∈ Rn , kxk = 1} . En efecto, puesto que esta norma es submultiplicativa, se sigue que: sup k n≥m n X tk Q k k=0 k! − m k k X t Q k=0 k! k≤ ∞ X k=m+1 |t| k kQkk →m→∞ 0. k! tQ Ahora veremos que e , t ≥ 0 es la única solución a las ecuaciónes backward y forward de Kolmogorov: ∂ tQ ∂ tQ e = QetQ y e = etQ Q. ∂t ∂t Además, satisfacen las ecuaciones de Chapman-Kolmogorov e(s+t)Q = esQ etQ . En efecto, Chapman-Kolmogorov se sigue de la definición de la exponencial de una matriz. Por otra parte, podemos derivar término a término la serie de potencias (cuyo radio de convergencia, con la norma matricial, es infinito) para obtener ( ∞ ∞ k k+1 X QetQ t Q ∂ tQ X tk−1 Qk e = = = , ∂t k! k! etQ Q k=1 k=0 lo cual muestra que se satisfacen las ecuaciones backward y forward. Además e0Q = Id. Para la unicidad de la solución a estas últimas, supongamos que Pt , t ≥ 0 es una colección de matrices en Rn que satisface las ecuaciones backward (el argumento para las forward es similar) y tal que P0 = Id. Notemos que la inversa de etQ es e−tQ . Entonces ∂ −tQ e Pt = −Qe−tQ Pt + e−tQ QPt = −Qe−tQ Pt + Qe−tQ Pt = 0. ∂t Por lo tanto e−tQ Pt es constante y como la constante es P0 = Id, vemos que Pt = etQ . 6. Distribuciones invariantes 93 6. Distribuciones invariantes Ahora pasaremos al estudio de las distribuciones invariantes para familias markovianas. La liga entre el tiempo continuo y discreto nos lo proporciona el siguiente resultado que se sigue de las ecuaciones backward de Kolmogorov. Definición. Decimos que una distribución ν en E (identificada con la colección numérica νx = ν({x})) es invariante para una familia Markoviana si X νx Pt (x, y) = νy . x En Potras palabras, la distribución ν es invariante si la distribución de Xt bajo Pν = x νx Px es igual a la de X0 . P Teorema 5.5. Una medida de probabilidad ν tal que x νx c(x) < ∞ es invariante para X si y sólo si cν = (cx νx , x ∈ E)) es invariante para la cadena asociada. Demostración. Por la ecuación backward de Kolmogorov y el teorema de Fubini-Tonelli se sigue que Z t XX X X (8) νx Pt (x, z) = νx P0 (x, z) + νx α(x, y) [Ps (y, z) − P (x, z)] ds. x x 0 y x Ası́, ν es invariante si y sólo si la integral del lado derecho es igual a 0 para cualquier t. Escribamos a α(x, y) = c(x) P P (x, y), donde P es la matriz de transición de la cadena asociada. Puesto que x cx νx y t 7→ Pt (x, y) es continua, podemos aplicar el teorema de convergencia dominada para concluir que el integrando en el lado derecho de (8) es continuo. Por lo tanto, ν es invariante si y sólo si XX 0= νx α(x, y) [P0 (y, z) − P0 (x, z)] y = x XX y x c(x) νx Px,y [1y=z − 1x=z ] = X x c(x) νx Px,z − cz νz En otras palabras, cν es invariante para la cadena asociada.  Recordemos que en el teorema fundamental de convergencia para cadenas de Markov (en tiempo discreto) la periodicidad juega un rol importante. Ahora veremos que en tiempo continuo, en cierto sentido el proceso ya es periódico. Proposición 5.8. Pt (x, y) > 0 para alguna t > 0 si y sólo si Pt (x, y) > 0 para toda t > 0. En particular Pt (x, x) > 0 para toda t ≥ 0. Demostración. El caso particular es simple: Pt (x, x) ≥ Px (T1 > t) > 0. 6. Distribuciones invariantes 94 Por otra parte, si y 6= x y para la cadena asociada se accede de x a y entonces existen x0 , . . . , xn ∈ E tales que x0 = x, xn = y y xk+1 6= xk para los cuales Px,x1 · · · Pxn−1 ,y > 0. En particular, se tiene que c(xi ) > 0 para i < n. Si S1 , . . . , Sn+1 son exponenciales de parámetro 1 independientes entonces   X Sk X Sk  P (x0 , x1 ) · · · Pxn−1 ,y > 0. Pt (x, y) ≥ P ≤t< c(xk−1 ) c(xk−1 ) k≤n k≤n+1 (Sólo se debe tener cuidado si c(y) = 0.) Finalmente, si de x no se accede a y para la cadena asociada Z entonces Px (Xt 6= y para toda t ≥ 0) = Px (Zn 6= y para toda n ≥ 0) = 1.  Una familia markoviana es irreducible si Px (Xt = y) > 0 para toda t > 0 y toda y ∈ E. Proposición 5.9. Si la cadena asociada a una familia markoviana irreducible es recurrente entonces no hay explosión. Lo anterior nos dice que los conjuntos {t ≥ 0 : Xt = y} y {n ∈ N : Zn = y} o son ambos acotados o ambos no acotados para familias markovianas irreducibles. En el primer caso hablamos de transitoriedad y en el segundo de recurrencia. Demostración. Sólo hay que notar que si Px (Zn = x i.o. ) entonces o Zn se absorbe en x (que sucede si y sólo si c(x) > 0 y no es compatible con la irreducibilidad de la cadena) ó c(x) > 0 y X 1 ≥ ∞/c(x) = ∞ c(Zn ) n Px -casi seguramente, en cuyo caso, al condicionar con Z, vemos que no hay explosión. (Recordemos que si τi son exponencialesPindependientes de parámetro λi P entonces τi = ∞ casi seguramente si y sólo si 1/λi = ∞.)  Teorema 5.6. Si (Px ) es una familia markoviana irreducible entonces son equivalentes: (1) Existe una única distribución invariante ν para la familia que satisface νx > 0 para toda x ∈ E y para cualquier distribución inicial µ: X lim |Pν (Xt = y) − νy | = 0. t→∞ x (2) Para alguna h > 0, la sucesión de variables aleatorias (Xnh , n ∈ N) es una cadena de Markov positivo recurrente. En caso contrario, no existe ninguna distribución invariante y Px (Xt = y) → 0 conforme t → ∞. 6. Distribuciones invariantes 95 Demostración. Sólo demostraremos la equivalencia. (La prueba completa se puede verificar en el libro de Kallenberg.) Sea h > 0. Notemos que (Xnh , n ≥ 0) es una cadena de Markov con matriz de transición Ph (x, y) , x, y ∈ E. En efecto, vemos que Px (Xh = x1 , . . . , Xnh = xn ) = Ph (x, x1 ) Ph (x1 , x2 ) · · · Ph (xn−1 , xn ) . Si para alguna h, dicha cadena de Markov es positivo recurrente, entonces al ser irreducible y aperiódica, existe una única distribución invariante νh . Por otra parte, la cadena de Markov Xnh/2 n ≥ 0 debe también ser positivo recurrente pues su tiempo de primer retorno está acotado por dos veces el tiempo de primer retorno de Xnh , n ≥ 0, el cual es integrable. Ası́, existe una única distribución invariante para Xnh/2 , digamos νh/2 pero como ésta también es invariante para Xnh , vemos que νh/2 = νh . Escribamos por lo tanto ν = ν h . Generalizando, vemos que para cualquier racional no-negativo q, la distribución de Xqh bajo Pν es ν y, al aproximar a cualquier t > 0 por la derecha por reales de la forma qh, vemos que ν es invariante para la familia markoviana. Para mostrar la convergencia en variación, notemos que, de acuerdo al teorema fundamental de convergencia para cadenas de Markov, se tiene que X |Pnh (x, y) − νy | → 0 x conforme n → ∞. Por lo tanto, al escribir a t (de manera única) en la forma nh + r con n ∈ N y 0 ≤ r < h, las ecuaciones de Chapman-Kolmogorov y la invariancia de ν nos dicen que X XX |Pt (x, y) − νy | ≤ |Pnh (x, z) − νz |Pr (z, y) → 0. x x y Por lo tanto, el teorema de convergencia dominada nos permite afirmar que X lim |Pν (Xt = y) − νy | = 0. t→∞ x Por otra parte, si existe una distribución invariante ν para la familia markoviana, entonces ν es una distribución invariante para Xnh , lo que implica que esta es positivo recurrente para cualquier h > 0.  Finalmente, pasamos a la relación entre el comportamiento asintótico de la probabilidad de transición y los tiempos medios de recurrencia. Sea T̃ y = min {t > T1 : Xt = y} . Teorema 5.7. Si y no es absorbente entonces lim Pt (x, y) = t→∞ Px (T y < ∞)  . c(y) Ey T̃ y 6. Distribuciones invariantes 96 Demostración. Sólo podremos probarlo en el caso transitorio y positivo recurrente. En el caso nulo recurrente, tendremos la convergencia en el sentido de Cesàro. Primero nos concentraremos en el caso x = y. Si y es transitorio entonces Ey (T̃ y ) = ∞ y por lo tanto el enunciado es válido. Si por otra parte y es positivo recurrente y nos concentramos en su clase de comunicación, esta será irreducible y sabemos que Pt (x, y) converge a νy donde ν es la distribución invariante única (en la clase de comunicación de y). Ası́, los tiempos medios de ocupación Z 1 t 1Xs =y ds Lt = t 0 satisfacen: Z Z 1 t 1 t Px (Xs = y) ds = Ps (x, y) ds →t→∞ νy . Ex (Lt ) = t 0 t 0 y y Por otra parte, si T̃ny = T̃ y + T̃n−1 (X T ) representa al tiempo del enésimo retorno de la cadena al estado y, la propiedad de Markov fuerte nos dice que T̃ny es una caminata aleatoria. Como T̃ y (X T1 ) se puede acotar en términos del tiempo de visita a y por la cadena XT1 +nh , n ≥ 0, que es finito por ser positivo recurrente, vemos que Ex (T̃ y ) < ∞, por lo que podemos aplicar la ley fuerte de los grandes números y deducir que bajo Py se tiene que T̃ny /n → Ey (T̃ny ). Por esto, observamos que LT̃ny ξ1 + · · · + ξn 1 → y = y c(y) E T̃n T̃n x (Ty ) donde ξi = T1 ◦ θT̃ny son variables exponenciales de parámetro c(y) (a la cuales también les aplicamos la ley fuerte de los grandes números). Finalmente, por convergencia dominada vemos que 1 Ex (Lt ) → , c(y) Ex (Ty ) lo cual prueba el resultado en este caso.