# Cadenas de Markov

Las cadenas de Markov se han vuelto un tema de gran desarrollo y aplicación en muy distintos ambitos científicos. Esto debido a que se trata de un proceso estocástico (no determinista, i.e. probabilista) que cumple con una gran propiedad a la cual se le conoce como *propiedad de Markov*. __Esta es que el proceso *carece de memoria*__. 

Al usar las CM (Cadenas de Markov), en general uno trabajo con *conjuntos de estados*. Tal cual lo hicimos en el Notebook pasado. La **carencia de memoria** se refiere a que un conjunto de estados al tiempo $t$,  $\mathbb{S}_{t}$ depende *únicamente* del conjunto de estados en el paso de tiempo anterior, i.e. $\mathbb{S}_{t-1}$.

De hecho (¡oh casualidad!), los caminantes aleatorios son un tipo particular de CM, ya que la posición al tiempo $t+1$ en una caminata depende sólo de la posición al tiempo $t$ (no de posiciones a tiempos anteriores). 

Supongamos ahora que $\mu_n$ es el estado del proceso al tiempo $n$. Por ahora nos interesaremos en procesos de Markov tal que $P(\mu_{n+1} = j \, | \, \mu_n = i) = p_{ij}$, independiente del tiempo. Lo que tenemos aquí es una *probabilidad condicional* la cual nos dice la probabilidad de que el sistema pase del estado $i$ al estado $j$.

Generalizando un poco, podemos llegar fácilmente a que $p_{ij}$ son en realidad entradas de una matriz -cuya importancia es *muy* grande-.

## Matriz de transición

Empecemos con ejemplo para ilustrar una matriz de transición. 

Sea un caminante aleatorio que se mueve en una arreglo de 6 celdas, $X = 0, \dotsc, 5$ con fronteras reflejantes y con una probabilidad $p$ de dar un paso a la derecha, $q$ de dar el paso a la izquierda y $r$ de quedarse en el mismo lugar.

No es dificil llegar a que la matriz de transición es:

$$
(p_{ij}) = 
\begin{pmatrix}
r + q & p & 0 & 0 & 0 & 0 \\
    q & r & p & 0 & 0 & 0 \\
    0 & q & r & p & 0 & 0 \\
    0 & 0 & q & r & p & 0 \\
    0 & 0 & 0 & q & r & p \\
    0 & 0 & 0 & 0 & q & r + p 
\end{pmatrix}
$$

De aquí notamos una primera propiedad muy importante que cualquier matriz de transición ha de cumplir. Si sumamos los elementos de cada fila, hemos entonces de encontrar que el resultado es igual a $1$. Dicho de otro modo:

$$ \sum_{j}p_{ij}=1$$

También podemos notar otra propiedad de este caso, y que tomaremos como general para nuestros propósitos. la matriz de transición $\mathbf{P}$ es **invariante en el tiempo**. 

Hagamos un segundo ejemplo. Supongamos una cadena que contiene 4 estados. Desde cada estado se puede saltar a cualquier otro estado con la misma probabilidad.

La matriz de transición en este caso será

$$
(p_{ij}) = 
\begin{pmatrix}
r & p & p & p \\
p & r & p & p \\
p & p & r & p \\
p & p & p & r 
\end{pmatrix}
$$

De aquí podemos llegar a otra propiedad sobre los elementos de la matriz de transición: $r \geq 0$. Esto podría ser evidente, pero tiene consecuencias importantes.

Evidentemente también para cualquier $i, \ j$: $0 \leq p_{ij} \leq 1$.

**[1]** Como primer ejercicio, el lector tendrá que hacer un `kernel` en el cual se calcule un proceso a partir de una condición inicial, con una matriz de transición y una duración arbitraria. 

Consideramos que para este ejercicio, el lector tendrá ya todas las herramientas, sino en la cabeza, en todos los notebooks anteriores.

## Enumeración exacta

Regresemos al método de enumeración exacta para cambiar un poco la perspectiva del problema y enfocarnos más en la *distribución de probabilidad*.

Supongamos que el estado empieza con una distribución $\mu_0 = \{ P_{0}(i) : i = 0, \dotsc, L \}$.

La distribución al paso $t=1$ estará dada por 
$$ \mu_{1} = (p_{ij}) \ \mu_{0} $$

¿Qué pasará en el paso $t=2$? Tendremos de nuevo que $\mu_{2} = (p_{ij}) \ \mu_{1}$. Sin embargo según la ecuación siguiente eso implica entonces que

$$\mu_{2} = (p_{ij})^2 \ \mu_0$$

El paso siguiente sería saltar a $t = n$, cuyo resultado es:

$$\mu_{n} = (p_{ij})^n \ \mu_0$$

**[2]** Implementa estas deducciones al código ya escrito y simula de nuevo el proceso.

**[3]** ¿Qué pasa con $\mu_t$ cuando $t \to +\infty$?

## Ergodicidad y balance

En la pregunta pasada se pidió buscar el comportamiento de $\mu_t$ cuando $t \to +\infty$. De esta manera se logró observar como $\mu_t$ convergía hacia una cierta distribución estable. Esta es una propiedad de suma profundidad, y se da únicamente con ciertas caracteristicas de las cuales hablaremos más adelante.

A esta distribución a la que se converge se le llama *distribución estacionaria* y la denotaremos como $\mu_{\infty}$.  Una derivación un poco más formal sería la siguiente:

Puesto que $\mu_t = (p_{ij})^t \ \mu_{0}$, entonces

$$
\begin{align}
\lim_{t \to \infty}\mu_{t}=\lim_{t \to \infty}(p_{ij})^t \ \mu_{0}
\end{align}
$$

Sin embargo $$\lim_{t \to \infty} (p_{ij})^t \ \mu_{0} = (p_{ij}) \lim_{t \to \infty} \mu_t$$ Y al hacer $$\lim_{t \to \infty} \mu_t = \mu_{\infty}$$ obtenemos entonces que

$$ \mu_{\infty} = (p_{ij}) \ \mu_{\infty} \ \ \ (1)$$

Esta última ecuación nos indica una de las propiedades más importantes de la *distribución estacionaria*:
$\mu_{\infty}$ es un __*eigenvector*__ de $(p_{ij})$ cuyo __*eigenvalor*__ $\lambda = 1$.

Otra información importante obtenida de (1) es que una vez que la *distribución estacionaria* es alcanzada, esta ya no cambiará (de ahí su nombre). Es decir que si tomamos como condición inicial $\mu_{\infty}$, la distribución será la misma en cualquier paso de tiempo.

Se puede demostrar (ver las referencias, especialmente el libro de Sethna) que esta distribución existe y es única si la cadena cumple con dos condiciones: es *irreducible* y *aperiódica*

La segunda condición se refiere a que en el proceso estocástico no existen órbitas en las cuales la cadena pueda caer. 

La condición de irreducibilidad se traduce como el hecho de que dada una cadena y sus estados, uno puede llegar de un estado dado a cualquier otro con un número finito de pasos. Algo así como que no hay camino que no pueda ser cruzado en un número finito de pasos.

A la suma de estas dos condiciones se le conoce como **ergodicidad**.

El hecho de que el proceso estocástico no sea *aperiódico* puede lograrse pidiendo la condición de __*balance*__. Con ella obtenemos que

$$ \mu \mathbf{P}_{\mu \nu} = \nu \mathbf{P}_{\nu \mu} $$

donde $\mu$ y $\nu$ son dos estados y $\mathbf{P}_{\mu \nu}$ y $\mathbf{P}_{\nu \mu}$ es la misma matriz de transición (el sub índice sólo denota el sentido de la transición. Por ejemplo $\mathbf{P}_{\mu \nu}$ denota la transición $\mu \to \nu$)

Con esta condición podemos decir que la probabilidad de ir de un estado a otro es la misma que hacer el camino de regreso, y por lo tanto los dos caminos son en promedio igual de frecuentes. Se puede demostrar que con este comportamiento las órbitas desaparecen. 

Para profundizar más sobre la formalidad y propiedades recomendamos leer el libro de Newman & Barkema listado en las referencias.

**[4]** Dada una matriz de transición, calcula la distribución estacionaria relacionada a la matriz. Puedes buscar distintos métodos para calcularla (tal vez encuentres algo en la librería cuBLAS).

**[5]** ¿Qué tan rápido converge una cadena de Markov hacia su distribución estacionaria?

## Conclusiones hasta el momento

Esperamos que hasta este momento el lector no sólo se sienta ya cómodo con PyCUDA, sino que también sus habilidades a la hora de escribir código tanto en PyCUDA como en CUDA C hayan aumentado. 

Este es el último notebook "introductorio" a algunas herramientas y tópicos de física estadística y a partir del siguiente notebook empezaremos a introducir el *modelo de Ising*.

Este modelo tiene grandes ventajas, como profundizar en algunos temas de física estadística, pero sobre todo el conocer y aprender los métodos MCMC (Markov Chain Monte Carlo) y el método de Metropolis-Hastings. 

¡Así que con grandes ansias los esperamos en los siguiente notebooks!

## Referencias

+ Wiki de [Markov Chain](https://en.wikipedia.org/wiki/Markov_chain)
+ Newman & Barkema, Monte Carlo Methods in Statistical Physics
+ [Sethna, Statistical Mechanics: Entropy, Order Parameters, and Complexity](http://pages.physics.cornell.edu/~sethna/StatMech/EntropyOrderParametersComplexity.pdf)
+ Documentación de [cuBLAS](http://docs.nvidia.com/cuda/cublas/)


```python

```

2a7e1e0d3091237827657ca66abeb2e5bfa77472
ipynb
Jupyter Notebook
Parte 2 - PyCUDA y aplicaciones/05 - Cadenas de Markov.ipynb
brincolab/Servicio_social
ac3cb224a2c1934e84c490a40420cf8bbad30235
MIT
2016-03-20T00:45:31.000Z
2020-11-25T23:54:22.000Z
Parte 2 - PyCUDA y aplicaciones/05 - Cadenas de Markov.ipynb
Sebzero77/Cuda_y_PyCuda
ac3cb224a2c1934e84c490a40420cf8bbad30235
MIT
2015-08-25T19:36:00.000Z
2018-04-20T07:10:23.000Z
Parte 2 - PyCUDA y aplicaciones/05 - Cadenas de Markov.ipynb
Sebzero77/Cuda_y_PyCuda
ac3cb224a2c1934e84c490a40420cf8bbad30235
MIT
2015-02-11T17:55:55.000Z
2020-11-25T23:54:34.000Z
Qwen/Qwen-72B
1. YES
2. YES