7.3 Algoritmo hacia adelante-hacia atrás (forward-backward)

Utilizamos el algoritmo forward-backward para calcular \(\hat{\gamma}_{ij}\) y \(\hat{\delta}_j\). Recordemos que

\[ \begin{aligned} \hat{\delta}_{j}&=P(S_{t}=j|x)=\frac{P(x, S_t=j)}{P(x)}\\ \hat{\gamma}_{ij}&=P(S_{t-1}=i,S_t=j|x)=\frac{P(x,S_{t-1}=i,S_t=j)}{P(x)} \end{aligned} \]

El cálculo de cada una de las probabilidades de arriba es computacionalmente intensivo, por ejemplo para calcular \(P(x)\):

\[P(x)=\sum_{\mathcal{S}} p_{s_1}(x_1)p_{s_1,s_2}p_{s_2}(x_2)\cdots p_{s_{T-1},s_T}p_{s_T}(x_T)\]

donde \(\mathcal{S}\) son las combinaciones de posibles estados (\(M^T\) posibilidades) por tanto esta aproximación no es factible. Es por esto que surge la necesidad de un algoritmo más eficiente.

El algoritmo hacia adelante-hacia atrás usa el principio de programación dinámica(recursión inteligente) para calcular \(\hat{\gamma}_{ij}\) y \(\hat{\delta}_j\) en tiempo lineal (\(M^2T\)), consta de dos pasos y explota las independencias condicionales del modelo.

Probabilidad hacia adelante

Definimos la probabilidad hacia adelante \(\alpha_i(t)\) como la probabilidad conjunta de observar las primeras \(t\) observaciones \(x^j\) (\(j=1,...,t\)) y siendo \(i\) el estado al tiempo \(t\):

\[\alpha_i(t)=P(X_1=x_1,...,X_T=x_t,S_t=i)\]

La probabilidad se puede evaluar de manera recursiva siguiendo la fórmula:

  1. \(\alpha_i(1) = \pi_k p_i(x_1)\) para \(i=1,...,M\)

  2. \(\alpha_i(t) = p_i(x_t)\sum_{j=1}^M \alpha_j(t-1)p_{j,i}\) para \(t=2,...,T\) e \(i=1,...,M\).

Prueba:

La idea clave es usar \((S_t,X_t) \perp (X_1,...,X_{t-1})|S_{t-1}\) \[ \begin{aligned} \alpha_i(t)&=P(x_1,...,x_t,S_t=i)\\ &=\sum_{j=1}^M P(x_1,...,x_t,S_t=i, S_{t-1}=j)\\ &=\sum_{j=1}^M P(x_1,...,x_{t-1},S_{t-1}=j)P(x_t|S_t=i)P(S_t=1|S_{t-1}=j)\\ &=\sum_{j=1}^M \alpha_j(t-1)p_i(x_t)p_{j,i}\\ \end{aligned} \]

Probabilidad hacia atrás

Definimos la probabilidad hacia atrás \(\beta_i(t)\) como la probabilidad condicional de las observaciones posteriores al tiempo \(t\) (\(x_{t+1},...,x_T\)) dado que el estado al tiempo \(t\) es \(i\).

\[\beta_i(t)=P(x_{t+1},...,x_T|S_t=i)\] para \(t=1,...T-1\).

La recursión de la probabilidad hacia atrás se evalúa como:

  1. \(\beta_i(T)=1\), para \(i=1,...,M\).

  2. \(\beta_i(t)=\sum_{i=1}^M p_{i,j}p_j(x_{t+1})\beta_i(t+1)\) para \(t = 1,...,T-1\).

Prueba:

La idea clave es usar \(X_{t+1} \perp (X_{t+2},...,X_T)|S_{t+1}\)

\[ \begin{aligned} \beta_i(t)&=P(x_{t+1},...,x_T|S_t=i) \\ &=\sum_{j=1}^M P(x_{t+1},...,x_T,S_{t+1}=j |S_t=i)\\ &=\sum_{j=1}^M P(S_{t+1}=j|S_t=i) P(x_{t+1},...,x_T|S_{t+1}=j)\\ &=\sum_{j=1}^M p_{i,j}P(x_{t+1},...,x_T|S_{t+1}=j)\\ &=\sum_{j=1}^M p_{i,j}P(x_{t+1}|S_{t+1}=j)P(x_{t+2},...,x_T|s_{t+1}=j)\\ &=\sum_{j=1}^M p_{i,j}p_j(x_{t+1})\beta_j(t+1) \end{aligned} \]

Escribimos \(\delta\) y \(\gamma\)

Ahora vemos como escrbir \(\delta_j\) y \(\gamma_{i,j}\) usando las probabilidades hacia adelante y hacia atrás:

\[ \begin{aligned} \hat{\delta}_{j}(t)&=P(S_{t}=j|x)\\ &=\frac{P(x,S_{t}=j)}{P(x)}\\ &=\frac{\alpha_j(t)\beta_j(t)}{\sum_{i=1}^M \alpha_i(T)} \end{aligned} \]

Prueba: \[ \begin{aligned} P(x_1,...,x_T,S_t=j)&=P(x_1,...,x_t,S_t=j)P(x_{t+1},...,x_T|S_t=j)\\ &=\alpha_j(t)\beta_j(t) \end{aligned} \]

Para el denominador notemos que:

\[ \begin{aligned} P(x) &= \sum_i^M P(x,S_{t}=i)\\ &= \sum_{i=1}^M \alpha_i(t)\beta_i(t) \end{aligned} \] esto se cumple para cualquier \(t\), así que si tomamos \(t=T\):

\[P(x) = \sum_{i=1}^M \alpha_i(T)\]

En el caso de \(\gamma_{i,j}\) tenemos:

\[ \begin{aligned} \hat{\gamma}_{ij}&=P(S_{t-1}=i,S_t=j|x)\\ &=\frac{P(x,S_{t-1}=i,S_t=j)}{P(x)}\\ &=\frac{\alpha_i(t-1)\beta_j(t)p_{i,j}p_j(x_{t})}{P(x)} \end{aligned} \]

Prueba: \[ \begin{aligned} P(x_1,...,x_T,S_{t-1}=i,S_t=j)&=P(x_1,...,x_{t-1},S_{t-1}=i)P(x_{t+1},...,X_T|S_t=j)P(S_t=j|S_{t-1}=i)P(x_t|S_t=j)\\ &=\alpha_i(t-1)\beta_j(t)p_{i,j}p_j(x_{t}) \end{aligned} \]

7.3.1 Resumen de algoritmo de estimación

Entonces, el algoritmo de estimación itera de la siguiente manera:

  1. Comenzamos con valores inciales \(\hat{\delta}_j\) y \(\hat{\gamma}_{ij}\).

  2. Actualizamos los parámetros \[\hat{p}_{ij}=\frac{\hat{\gamma}_{ij}}{\sum_l \hat{\gamma}_{il}}\] y los correspondientes a las densidades \(p_j(x)\).

  3. Utilizando el conjunto de parámetros actuales (\(\hat{p}_{ij}\) y los correspondientes a \(p_j(x)\)) calculamos \(\hat{\delta}_j\) y \(\hat{\gamma}_{ij}\) a través del algoritmo hacia adelante-hacia atrás.

  4. Iteramos entre 2 y 3.

7.3.2 Algoritmo de Viterbi

En muchas de las aplicaciones de HMM nos interesa hacer inferencia de la secuencia de estados \(\{S_1,...,S_T\}\), en este caso el criterio de optimización es:

\[MAP(S|x)=argmax_s P(S|x) = argmax_s P(S,x)\]

Aqui estamos buscando el camino más probable. Si consideramos un algoritmo de fuerza bruta, esto es, realizamos búsqueda exahustiva sobre todas las posibles secuencias tendríamos que considerar \(M^T\) casos. Es por ello que nuevamente recurrimos a un algoritmo de programación dinámica.