⛹🏻 😓 📮 Desmontamos el algoritmo EM en pequeños ladrillos 🔘 👃🏽 ➰

En este artículo, como probablemente ya haya adivinado, hablaremos sobre el algoritmo EM del dispositivo. En primer lugar, el artículo puede ser de interés para quienes se unen lentamente a la comunidad de datasainists. El material presentado en el artículo será más útil para aquellos que recientemente comenzaron a tomar el tercer curso "Búsqueda de estructuras de datos" dentro de la especialización "Aprendizaje automático y análisis de datos" de MIPT y Yandex.

El material presentado en el artículo, en cierto sentido, es una adición a la primera semana de capacitación en el curso antes mencionado, es decir, le permite responder algunas preguntas importantes con respecto al principio de funcionamiento del algoritmo EM. Para una mejor comprensión del material, es recomendable que nuestro estimado lector pueda realizar operaciones con matrices (multiplicación de matrices, encontrar el determinante de una matriz y una matriz inversa), comprender los conceptos básicos de la teoría de la probabilidad y el matstat y, por supuesto, tener al menos una comprensión básica de los algoritmos de agrupamiento básicos y comprender qué la agrupación tiene lugar en el aprendizaje automático. Aunque, por supuesto, sin este conocimiento, puede leer el artículo, algo y sí, sin duda será claro :)

Además, según la antigua tradición, el artículo no contendrá una investigación teórica profunda, sino que estará lleno de ejemplos simples y comprensibles. Cada ejemplo posterior explicará el funcionamiento del algoritmo EM un poco más profundo que el anterior, lo que finalmente nos lleva directamente al análisis del algoritmo mismo. Para cada ejemplo, se escribirá el código. Todo el código está escrito en Python 2.7, y por esto me disculpo de antemano. Resultó que ahora estoy usando esta versión en particular, pero después de cambiar a python 3, intentaré cambiar el código en el artículo.

Conozcamos el bosquejo del artículo: 1) Consideremos en términos generales cómo funciona el algoritmo EM.

2) Repetimos los puntos principales del teorema de Bayes (fórmula):

P_{x_{i}} (B_{j}) = \frac{P (B_{j}) \cdot P_{B_{j}} (x_{i})}{P (x_{i})}

$P_{x_i}(B_j) = \frac{P(B_j) \centerdot P_{B_j}(x_i)}{P(x_i)}$

3) Considere el primer ejemplo del teorema bayesiano, en el que todos los elementos de la fórmula (o más bien probabilidad

P (B_{j}), P_{B_{j}} (x_{i})

$P(B_j), P_{B_j}(x_i)$ ) después del signo igual a la derecha son conocidos. En este caso, solo necesitamos entender correctamente qué números sustituir y luego realizar operaciones básicas.

4) Considere el segundo ejemplo del teorema de Bayes. Para resolver el problema, necesitaremos calcular adicionalmente la densidad de probabilidad de un determinado evento (

x_{i}

$x_i$ ) sujeto a la hipótesis (

B_{j}

$B_j$ ) -

ρ_{B_{j}} (x_{i})

$\rho_{B_j}(x_i)$ . Para los cálculos, se nos darán parámetros de una variable aleatoria, a saber, la expectativa matemática:

μ

$\mu$ y desviación estándar -

σ

$\sigma$ . Los cálculos de densidad de probabilidad se llevarán a cabo de acuerdo con la siguiente fórmula:

ρ_{B_{j}} (x_{i}) = \frac{1}{σ_{j} \cdot \sqrt{2 π}} \cdot \exp (- \frac{(x_{i} - μ_{j})^{2}}{2 σ_{j}^{2}})

$\rho_{B_j}(x_i) = \frac{1}{\sigma_j \centerdot \sqrt{2 \pi}} \centerdot \exp(- \frac{(x_i-\mu_j)^2}{2 \sigma_j^2})$

El uso de la fórmula anterior supone que la variable aleatoria tiene una dimensión unidimensional y se distribuye normalmente (en otras palabras, la distribución gaussiana o gaussiana-Laplace es una distribución de probabilidad).

5) Considere el tercer ejemplo, que difiere del anterior solo en que considera el caso de una distribución normal multidimensional (en nuestra versión bidimensional) de una variable aleatoria. En este caso, el cálculo de densidad se realiza de acuerdo con la siguiente fórmula:

ρ_{B_{j}} (x_{i}) = \frac{1}{(2 π)^{k / 2} \cdot ∣ Σ_{j} ∣^{1 / 2}} \cdot \exp (- \frac{1}{2} \cdot (x_{i} - μ_{j})^{T} Σ_{j}^{- 1} (x_{i} - μ_{j}))

$\rho_{B_j}(x_i) = \frac{1}{(2\pi)^{k/2} \centerdot \mid \Sigma_j\mid ^{1/2}} \centerdot \exp(- \frac{1}{2} \centerdot (x_i-\mu_j)^T \Sigma_j^{-1} (x_i-\mu_j))$

6) Finalmente, modificamos el tercer ejemplo de tal manera que demuestre claramente el funcionamiento del algoritmo EM.

7) En conclusión, comparamos la calidad del algoritmo que hemos deducido con la calidad del algoritmo EM, que está integrado en la biblioteca sklearn (clase sklearn.mixture.GaussianMixture)

Si, en lugar de las fórmulas especificadas en el esquema del artículo, vio caracteres japoneses, entonces no se asuste - Después del teorema de Bayes (el segundo párrafo del bosquejo del artículo), analizaremos tales ejemplos que probablemente entenderá al menos intuitivamente. ¡Bueno, vamos!

Algoritmo EM en general

En el curso, en función del cual se escribió el artículo, el algoritmo EM se presenta como uno de los métodos de agrupación. En otras palabras, este es un método de aprendizaje automático sin un maestro, cuando no sabemos las respuestas verdaderas de antemano. En nuestro artículo, también consideraremos el algoritmo como parte de uno de los métodos de agrupación de una variable aleatoria con una distribución normal. Sin embargo, debe tenerse en cuenta que el algoritmo tiene una aplicación más amplia. El algoritmo EM (maximización de expectativas en inglés) es un método general para encontrar estimaciones de la función de probabilidad en modelos con variables ocultas, que a partir de una mezcla de distribuciones le permite construir distribuciones de probabilidad complejas (aproximadas).

El algoritmo EM en problemas de agrupamiento se utiliza como un algoritmo iterativo, que realiza dos pasos en cada iteración:

E-step. En el primer paso E, de alguna manera, por ejemplo, al azar, seleccionamos variables ocultas, en nuestro caso será una expectativa matemática:

μ

$\mu$ y desviación estándar -

σ

$\sigma$ . Usando las variables seleccionadas, calculamos la probabilidad de asignar cada objeto a un grupo particular. Los pasos E posteriores utilizan las variables ocultas definidas en los pasos M.

M-step. En el paso M, de acuerdo con las probabilidades de asignar cada objeto a un grupo particular obtenido en el paso E, recalculamos las variables ocultas

μ

$\mu$ y

σ

$\sigma$

iteraciones se repiten hasta que ocurre la convergencia. Podemos suponer que la convergencia ocurrió cuando los valores de las variables ocultas no cambian significativamente, por ejemplo, dentro de una constante dada. En el último ejemplo, que se considera en el artículo, no estableceremos constantes y, en consecuencia, no determinaremos cambios en los valores de las variables ocultas en cada iteración. Hagámoslo más simple: limitaremos nuestro algoritmo a un número fijo de pasos.

El paso M es lo suficientemente simple como para entenderlo, y lo veremos en el último ejemplo. Nos centraremos en el paso E del artículo.

El paso E se basa en el teorema bayesiano, en base al cual determinamos la probabilidad de que cada objeto sea asignado a uno u otro grupo. Recordemos de qué trata este teorema.

Recordemos la fórmula de Bayes

P_{x_{i}} (B_{j}) = \frac{P (B_{j}) \cdot P_{B_{j}} (x_{i})}{P (x_{i})}

$P_{x_i}(B_j) = \frac{P(B_j) \centerdot P_{B_j}(x_i)}{P(x_i)}$

dónde

P_{x_{i}} (B_{j})

$P_{x_i}(B_j)$ - la probabilidad de que con el evento

x_{i}

$x_i$ una hipótesis

B_{j}

$B_j$

P (B_{j})

$P(B_j)$ - probabilidad de hipótesis

B_{j}

$B_j$

P_{B_{j}} (x_{i})

$P_{B_j}(x_i)$ - probabilidad de ocurrencia de un evento

x_{i}

$x_i$ sujeto a la hipótesis

B_{j}

$B_j$

P (x_{i})

$P(x_i)$ - probabilidad de ocurrencia de un evento

x_{i}

$x_i$ siempre que cada hipótesis se cumple (calculado por la fórmula para la probabilidad total del evento)

A condición de que el evento

x_{i}

$x_i$ ya ha ocurrido, la probabilidad de que se cumpla la hipótesis

B_{j}

$B_j$ revalúa utilizando la fórmula anterior. Podemos decir eso

P (B_{j})

$P(B_j)$ es la probabilidad a priori de la hipótesis (antes de la prueba), y

P_{x_{i}} (B_{j})

$P_{x_i}(B_j)$ es la probabilidad posterior de la misma hipótesis evaluada después del evento

x_{i}

$x_i$ , es decir, teniendo en cuenta que el evento

x_{i}

$x_i$ fiable sucedió.

Ejemplo uno por el teorema de Bayes

Imagine que recibimos piezas producidas en dos máquinas diferentes. Se conocen las siguientes características de los productos recibidos en el almacén:

en total, piezas - 10,000 piezas (

N = 10000

$N = 10000$ ), de las cuales piezas producidas en la primera máquina - 6000 piezas (

N 1 = 6000

$N1 = 6000$ ), en el segundo - 4000 piezas (

N 2 = 4000

$N2 = 4000$ )

La proporción de productos estándar (es decir, no defectuosos) fabricados en la primera máquina es de 0.9, la proporción de productos estándar fabricados en la segunda máquina es de 0.8.

Extrajimos al azar una parte de las partes recibidas, y resultó ser estándar.

Necesitamos determinar la probabilidad de que:

1) la pieza se produzca en la primera máquina;

2) el artículo fue producido en la segunda máquina.

Decisión

1) Evento

x

$x$ en nuestro ejemplo, extrayendo un producto estándar.

2) Ahora, determinemos las hipótesis

B_{j}

$B_j$ . Tenemos dos partes, lo que significa dos hipótesis:

Hipótesis.

B_{1}

$B_1$ : un producto que fue retirado accidentalmente fue producido en la máquina No. 1.

Hipótesis

B_{2}

$B_2$ : un producto que fue retirado accidentalmente fue producido en la máquina No. 2.

En consecuencia, la probabilidad de extraer un producto fabricado en la primera máquina es

P (B_{1})

$P(B_1)$ constituye

N 1 / N = 6000 / 10000 = 0.6

$N1/N = 6000/10000 = 0.6$ .

La probabilidad de extraer un producto fabricado en una segunda máquina es

P (B_{2})

$P(B_2)$ constituye

N 2 / N = 0.4

$N2/N = 0.4$ . Además, no nos importa el producto estándar o el defectuoso, es importante de qué lote es.

3) La probabilidad de extraer un producto estándar de productos fabricados en la máquina No. 1 (es decir, según la primera hipótesis) corresponde a la proporción de productos estándar producidos en la máquina No. 1 y es 0.9,

P_{B_{1}} (x_{i}) = 0.9

$P_{B_1}(x_i) = 0.9$ .

La probabilidad de extraer el producto es un producto estándar sujeto a la hipótesis No. 2 -

P_{B_{2}} (x_{i}) = 0.8

$P_{B_2}(x_i) = 0.8$

4) Determinar la probabilidad de extraer un producto estándar de todo el conjunto de productos:

P (x_{i})

$P(x_i)$ . De acuerdo con la fórmula para la probabilidad total de un evento

P (x_{i}) = P (B_{1}) \cdot P_{B_{1}} (x_{i}) + P (B_{2}) \cdot P_{B_{2}} (x_{i}) = 0.6 \cdot 0.9 + 0.4 \cdot 0.8 = 0.54 + 0.32 = 0.86

$P(x_i) = P(B_1) \centerdot P_{B_1}(x_i) + P(B_2) \centerdot P_{B_2}(x_i) = 0.6 \centerdot 0.9 + 0.4 \centerdot 0.8 = 0.54 + 0.32 = 0.86$ . Aquí entendemos que entonces

0.14

$0.14$ - esta es la probabilidad de extraer la parte defectuosa de todas las partes que llegaron al almacén y en total resulta 1

(0.86 + 0.14 = 1)

$(0.86 + 0.14 = 1)$ .

5) Ahora tenemos todos los datos para resolver el problema.

5.1) Determine la probabilidad de que una pieza estándar extraída al azar se produzca en una máquina

№ 1

$№1$ :

P_{x_{i}} (B_{1}) = \frac{P (B_{1}) \cdot P_{B_{1}} (x_{i})}{P (x_{i})} = \frac{0.6 \cdot 0.9}{0.86} \approx 0.63

$P_{x_i}(B_1) = \frac{P(B_1) \centerdot P_{B_1}(x_i)}{P(x_i)} = \frac{0.6 \centerdot 0.9}{0.86} \approx 0.63$

5.2) Determinamos la probabilidad de que una pieza estándar extraída al azar se produzca en una máquina

№ 2

$№2$ :

P_{x_{i}} (B_{2}) = \frac{P (B_{2}) \cdot P_{B_{2}} (x_{i})}{P (x_{i})} = \frac{0.4 \cdot 0.8}{0.86} \approx 0.37

$P_{x_i}(B_2) = \frac{P(B_2) \centerdot P_{B_2}(x_i)}{P(x_i)} = \frac{0.4 \centerdot 0.8}{0.86} \approx 0.37$

Por lo tanto, reevaluamos las hipótesis

B_{1}

$B_1$ y

B_{2}

$B_2$ . Después de la reevaluación, las hipótesis también forman un grupo completo de eventos:

0.63 + 0.37 = 1

$0.63 + 0.37 = 1$ .

Bueno, ahora es el momento de pasar al segundo ejemplo.

Segundo ejemplo para el teorema de Bayes usando los parámetros de la distribución normal de una variable aleatoria: expectativa matemática $\mu$ y desviación estándar $\sigma$

Modifiquemos ligeramente las condiciones del ejemplo anterior. Suponemos que no conocemos la proporción de productos estándar fabricados en las máquinas No. 1 y No. 2, pero supongamos, en cambio, que conocemos las dimensiones promedio de los diámetros de los productos y la desviación estándar del diámetro de los productos en cada máquina.

La máquina No. 1 produce piezas de 64 mm de diámetro y una desviación estándar de 4 mm.

La máquina No. 2 produce piezas de 52 mm de diámetro y una desviación estándar de 2 mm.

Una condición importante es que todo el conjunto de productos se describe mediante la distribución normal o distribución gaussiana.

De lo contrario, las condiciones son las mismas, las escribimos:

N 1 = 6000, μ_{1} = 64, σ_{1} = 4

$N1 = 6000, \mu_1 = 64, \sigma_1 = 4$

N 2 = 4000, μ_{2} = 52, σ_{2} = 2

$N2 = 4000, \mu_2 = 52, \sigma_2 = 2$

Suponemos que en el proceso de aceptación del producto hubo un pequeño incidente, como resultado del cual todos los productos se mezclaron.

Nuestra tarea es ordenar los detalles y para cada uno determinar la probabilidad de que se haya producido en la máquina No. 1 o en la máquina No. 2. También asumiremos que la pieza se produjo en la máquina, cuya probabilidad es mayor.

Solución En

primer lugar, analizaremos el algoritmo de solución. Podemos calcular fácilmente la probabilidad de cada hipótesis, sin embargo, ya conocemos sus valores de un ejemplo pasado:

P (B_{1}) = 0.6

$P(B_1)=0.6$ ,

P (B_{2}) = 0.4

$P(B_2)=0.4$ . Pero, ¿cómo calculamos la probabilidad de confiscar un producto estándar individualmente para cada una de las hipótesis disponibles? ¡Todo es simple! No consideraremos la probabilidad, sino que determinaremos la densidad de probabilidad del elemento que se eliminará con su valor de diámetro para cada hipótesis por separado.

Para hacer esto, utilizamos la conocida fórmula:

ρ_{B_{j}} (x_{i}) = \frac{1}{σ_{j} \cdot \sqrt{2 π}} \cdot \exp (- \frac{(x_{i} - μ_{j})^{2}}{2 σ_{j}^{2}})

$\rho_{B_j}(x_i) = \frac{1}{\sigma_j \centerdot \sqrt{2 \pi}} \centerdot \exp(- \frac{(x_i-\mu_j)^2}{2 \sigma_j^2})$

dónde

x_{i}

$x_i$ - un valor aleatorio (en nuestro caso, el tamaño real del diámetro del producto),

μ_{j}

$\mu_j$ - expectativa matemática de variables aleatorias

j

$j$ hipótesis (en nuestro caso, el tamaño promedio del diámetro del producto producido en

j

$j$ herramienta de máquina)

σ_{j}

$\sigma_j$ - desviación estándar de variables aleatorias

j

$j$ hipótesis (en nuestro caso, la desviación del tamaño promedio del diámetro del producto producido en

j

$j$ herramienta de máquina).

Por lo tanto, hacemos una sustitución inteligente de la probabilidad de densidad de probabilidad en el numerador de la fórmula de Bayes, llevaremos a cabo reemplazos similares en el denominador.

Adaptamos las fórmulas específicamente para nuestro caso.

La probabilidad de que la pieza con dimensiones de diámetro

x_{i}

$x_i$ producido en la máquina No. 1 está determinado por la fórmula:

P_{x_{i}} (B_{1}) = \frac{w_{1} \cdot ρ_{B_{1}} (x_{i})}{w_{1} \cdot ρ_{B_{1}} (x_{i}) + w_{2} \cdot ρ_{B_{1}} (x_{i})}

$P_{x_i}(B_1) = \frac{w_1 \centerdot \rho_{B_1}(x_i)}{w_1 \centerdot \rho_{B_1}(x_i) + w_2 \centerdot \rho_{B_1}(x_i)}$

La probabilidad de que la pieza con dimensiones de diámetro

x_{i}

$x_i$ producido en la máquina No. 2:

P_{x_{i}} (B_{2}) = \frac{w_{2} \cdot ρ_{B_{2}} (x_{i})}{w_{1} \cdot ρ_{B_{1}} (x_{i}) + w_{2} \cdot ρ_{B_{1}} (x_{i})}

$P_{x_i}(B_2) = \frac{w_2 \centerdot \rho_{B_2}(x_i)}{w_1 \centerdot \rho_{B_1}(x_i) + w_2 \centerdot \rho_{B_1}(x_i)}$

Tenga en cuenta que hemos reemplazado la notación de la probabilidad de una hipótesis con

P (B_{j})

$P(B_j)$ sobre el

w_{j}

$w_j$ . Esto se debe a la designación correspondiente en el curso previamente designado. Además, usaremos tal notación.

Ahora estamos 100% completos y listos para resolver el problema.

Nos fijamos en el código

Importar bibliotecas

#  , 
import numpy as np
import random
import math
import matplotlib.pyplot as plt
%matplotlib inline
from sklearn.metrics import accuracy_score

Funciones utilizadas en el segundo ejemplo

#      ,    
#      : .,  . 
def gaus_func_01(mu,sigma,x):
    return math.e**(-(x-mu)**2/(2*sigma**2)) / (sigma*(2*math.pi)**0.5)

#        
def proba_x(X, w1, w2, mu_1, mu_2, sigma_1, sigma_2):
    for i in X:
        P1_x = gaus_func_01(mu_1,sigma_1,i)
        P2_x = gaus_func_01(mu_2,sigma_2,i)
        P_x = w1*P1_x + w2*P2_x
        P_x_1 = (w1*P1_x)/P_x
        P_x_2 = (w2*P2_x)/P_x
        proba_temp = []
        proba_temp.append(P_x_1)
        proba_temp.append(P_x_2)
        proba_X.append(proba_temp)
    return proba_X

#        
def pred_x(proba_X, limit_proba):
    pred_X = []
    for x in proba_X:
        if x[0] >= limit_proba:
            pred_X.append(1)
        else:
            pred_X.append(2)
    return np.array(pred_X)

#    
def graph_01(X, pred_X, mu_1, sigma_1, mu_2, sigma_2):
    true_pred = []
    false_pred_1 = []
    false_pred_2 = []
    for i in range(X.shape[0]):
        if pred_X[i] == y[i]:
            true_pred.append([X[i], -0.025])
        else:
            if y[i] == 1:
                false_pred_1.append([X[i], -0.0075])
            else:
                false_pred_2.append([X[i], -0.015])

    false_pred_1 = np.array(false_pred_1)            
    false_pred_2 = np.array(false_pred_2)
    true_pred = np.array(true_pred)

    x_theory = np.linspace(42, 85, 20000)
    y_theory_1 = []
    for x in x_theory:
        y_theory_1.append(gaus_func_01(mu_1,sigma_1,x))
    y_theory_2 = []
    for x in x_theory:
        y_theory_2.append(gaus_func_01(mu_2,sigma_2,x))

    plt.figure(figsize=(18, 8))    
    plt.plot(
        x_theory, y_theory_1, color = 'green', lw = 2, label = 'Theoretical probability density for machine 1')
    plt.plot(
        x_theory, y_theory_2, color = 'firebrick', lw = 2, label = 'Theoretical probability density for machine 2')
    plt.hist(
        X[:N1], bins = 'auto', color='#539caf', normed = True, alpha = 0.35, label = 'machine tool products 1')
    plt.hist(
        X[N1:N], bins = 'auto', color='sandybrown', normed = True, alpha = 0.75, label = 'machine tool products 2')
    plt.plot(mu_1, 0, 'o', markersize = 11, color = 'blue', label = 'Mu 1')
    plt.plot(mu_2, 0, 'o', markersize = 11, color = 'red', label = 'Mu 2')

    plt.plot([mu_1 - sigma_1, mu_1 - sigma_1], [0,0.85*np.max(y_theory_1)],
             ':', lw = 3, color = 'blue', alpha = 0.55, label = 'Mu1 - sigma1')
    plt.plot([mu_1 + sigma_1, mu_1 + sigma_1], [0,0.85*np.max(y_theory_1)],
             ':', lw = 3, color = 'blue', alpha = 0.55, label = 'Mu1 + sigma1')
    plt.plot([mu_2 - sigma_2, mu_2 - sigma_2], [0,0.85*np.max(y_theory_2)],
             ':', lw = 3, color = 'red', alpha = 0.55, label = 'Mu2 - sigma2')
    plt.plot([mu_2 + sigma_2, mu_2 + sigma_2], [0,0.85*np.max(y_theory_2)],
             ':', lw = 3, color = 'red', alpha = 0.55, label = 'Mu2 + sigma2')

    plt.plot([mu_1 - 2 * sigma_1, mu_1 - 2 * sigma_1], [0, 0.9*0.5 * np.max(y_theory_1)],
             ':', lw = 2.5, color = 'blue', alpha = 0.35, label = 'Mu1 - 2*sigma1')
    plt.plot([mu_1 + 2 * sigma_1, mu_1 + 2 * sigma_1], [0, 0.9*0.5 * np.max(y_theory_1)],
             ':', lw = 2.5, color = 'blue', alpha = 0.35, label = 'Mu1 + 2*sigma1')
    plt.plot([mu_2 - 2 * sigma_2, mu_2 - 2 * sigma_2], [0, 0.9*0.5 * np.max(y_theory_2)],
             ':', lw = 2.5, color = 'red', alpha = 0.35, label = 'Mu2 - 2*sigma2')
    plt.plot([mu_2 + 2 * sigma_2, mu_2 + 2 * sigma_2], [0, 0.9*0.5 * np.max(y_theory_2)],
             ':', lw = 2.5, color = 'red', alpha = 0.35, label = 'Mu2 + 2*sigma2')

    plt.plot(false_pred_1[:,0], false_pred_1[:,1], 'o', markersize = 2.5, color = 'blue', alpha = 0.2, label = 'errors1')
    plt.plot(false_pred_2[:,0], false_pred_2[:,1], 'o', markersize = 2.5, color = 'red', alpha = 0.3, label = 'errors2')
    plt.plot(true_pred[:,0], true_pred[:,1], 'o', markersize = 3, color = 'green', alpha = 0.2, label = 'right answers')

    plt.xlabel('Caliber')
    plt.ylabel('Probability density')
    plt.legend()
    plt.show()

Especificar los parámetros

#    
#      №1
N1 = 6000
#      №2
N2 = 4000
#      
N = N1+N2

#    №1
mu_1 = 64.
#        №1
sigma_1 = 3.5

#    №2
mu_2 = 52
#        №2
sigma_2 = 2.

Nosotros calcularemos

X = np.zeros((N))
np.random.seed(seed=42)
#    ,   №1
X[:N1] = np.random.normal(loc=mu_1, scale=sigma_1, size=N1)
#  ,   №2
X[N1:N] = np.random.normal(loc=mu_2, scale=sigma_2, size=N2)

#   
y = np.zeros((N))
y[:N1] = np.array((1))
y[N1:N] = np.array((2))

#     ,    №1
w1 = float(N1)/N
#     ,    №2
w2 = float(N2)/N

#           
proba_X = []
proba_X = proba_x(X, w1, w2, mu_1, mu_2, sigma_1, sigma_2)

#   ,   ,        
limit_proba = 0.5

#     
pred_X = []
pred_X = pred_x(proba_X, limit_proba)

#    
print '   :', round(accuracy_score(y, pred_X),3)
print

print ' №1'
graph_01(X, pred_X, mu_1, sigma_1, mu_2, sigma_2)

Cuadro # 1

¿Qué acabamos de hacer? Generamos de manera pseudoaleatoria 10,000 valores descritos por la distribución normal, de los cuales 6,000 valores con una media de 64 mm y una desviación estándar de 4 mm, 4,000 valores con una media de 52 mm y una desviación estándar de 2 mm. Además, de acuerdo con el algoritmo anterior, para cada parte, se determinó la probabilidad de su producción en la máquina No. 1 o No. 2. Después de eso, elegimos la hipótesis sobre la producción del producto en una u otra máquina, dependiendo de qué hipótesis es más probable. Finalmente, comparamos los resultados de nuestro algoritmo con las respuestas verdaderas.

En la salida, obtuvimos una parte de las respuestas correctas: 0.986. En general, esto es muy bueno, dado que llevamos a cabo la capacitación sin usar respuestas verdaderas.

Miremos la tabla. ¿A qué se recomienda prestar atención? Vea dónde se encuentran los productos que el algoritmo no identifica correctamente.

En primer lugar, vemos errores de algoritmos solo en el área donde los productos fabricados en diferentes máquinas tienen el mismo diámetro. Se ve bastante lógico.

En segundo lugar, vemos que el algoritmo se equivoca solo en aquellos objetos que están más distantes del verdadero valor promedio de los objetos y, al mismo tiempo, están bastante cerca del centro falso.

Pero lo más interesante es que los problemas comienzan principalmente después de cruzar valores aproximadamente iguales

μ + - 2 σ

$\mu +- 2\sigma$ , específicamente en nuestro caso en la intersección

μ_{1} - 2 σ_{1}

$\mu_1 - 2\sigma_1$ y

μ_{2} + 2 σ_{2}

$\mu_2 + 2\sigma_2$

Aquellos que lo deseen pueden "torcer" los parámetros y ver cómo cambia la calidad del algoritmo. Esto será útil para una mejor asimilación del material.

Bueno, pasamos al siguiente ejemplo

El tercer ejemplo. Caja bidimensional

Esta vez consideraremos el caso cuando tengamos datos no solo sobre el diámetro de los productos, sino también, por ejemplo, sobre el peso de cada producto. Deje que la máquina No. 1 produzca piezas que pesen 12 gy una desviación estándar de 1 g, la máquina No. 2 produzca productos que pesen 10 gy una desviación estándar de 0.8 mm.

Para el resto, usaremos los datos del ejemplo anterior, escríbalos.

N 1 = 6000, μ_{11} = 64, σ_{11} = 4, μ_{12} = 14, σ_{12} = 1

$N1 = 6000, \mu_{1 1} = 64, \sigma_{1 1} = 4, \mu_{1 2} = 14, \sigma_{1 2} = 1$

N 2 = 4000, μ_{21} = 52, σ_{21} = 2, μ_{22} = 10, σ_{22} = 0.9

$N2 = 4000, \mu_{2 1} = 52, \sigma_{2 1} = 2, \mu_{2 2} = 10, \sigma_{2 2} = 0.9$

La situación es la misma: todos los detalles se mezclaron en una gran pila y necesitamos clasificarlos y determinar para cada parte la probabilidad de su producción en la máquina No. 1 y la máquina No. 2.

Solución

En general, el enfoque de solución del ejemplo anterior no cambia, con la excepción de la fórmula para determinar la densidad de probabilidad de una variable aleatoria. Para el caso multidimensional, se utiliza la siguiente fórmula:

p_{B_{j}} (x_{i}) = \frac{1}{(2 π)^{k / 2} \cdot ∣ Σ_{j} ∣^{1 / 2}} \cdot \exp (- \frac{1}{2} \cdot (x_{i} - μ_{j})^{T} Σ_{j}^{- 1} (x_{i} - μ_{j}))

$p_{B_j}(x_i) = \frac{1}{(2\pi)^{k/2} \centerdot \mid \Sigma_j\mid ^{1/2}} \centerdot \exp(- \frac{1}{2} \centerdot (x_i-\mu_j)^T \Sigma_j^{-1} (x_i-\mu_j))$

dónde

Σ_{j}

$\Sigma_j$ Es una matriz de covarianzas de variables aleatorias.

j

$j$ hipótesis (en nuestro caso, la matriz de covarianzas de los parámetros de los productos producidos en

j

$j$ herramienta de máquina)

∣ Σ_{j} ∣

$\mid \Sigma_j\mid$ Es el determinante de la matriz de covarianza de variables aleatorias.

j

$j$ hipótesis

μ_{j}

$\mu_j$ - expectativa matemática de variables aleatorias

j

$j$ hipótesis (en nuestro caso, el valor promedio de los productos fabricados en

j

$j$ herramienta de máquina)

x_{i}

$x_i$ - variable aleatoria (en nuestro caso, los parámetros

i

$i$ -º producto)

Esta vez, para comprender el código, un lector respetado necesitará conocimientos sobre operaciones con matrices

Funciones utilizadas en el tercer ejemplo

#      ,    
#      : .,  . 
def gaus_func_02(k, m, x, mu, sigma):
    pj_xi = []
    for j in range(k):
        det_sigma_j = np.linalg.det(sigma[j])
        factor_1 = 1 / (((2 * math.pi)**(k/2)) * ((det_sigma_j)**0.5))
        factor_2 = []
        for i in x:
#    ,            
#          
#   .    ,      
#     
            factor_2.append(math.e**float(
            -0.5 * np.matrix(i - mu[j]) * np.matrix(np.linalg.inv(sigma[j])) * np.matrix(i - mu[j]).T))
        pj_xi.append(factor_1 * np.array(factor_2))
    return np.array(pj_xi)

#     ,       №1  №2
def proba_func_02(pjxi, w, k):
    #          
    P_X = []
    for j in range(k):
        P_X.append(w[j] * pjxi[j])
    P_X = np.sum(np.array(P_X), axis = 0)
    #    ,       №1  №2
    P_J_X = []
    for j in range(k):
        P_J_X.append(w[j] * pjxi[j] / P_X)
    return np.array(P_J_X)

#        
def pred_x_02(proba_X, limit_proba):
    pred_X = []
    for x in proba_X[0]:
        if x >= limit_proba:
            pred_X.append(1)
        else:
            pred_X.append(2)
    return np.array(pred_X)

#             
def graph_02_algorithm(pred_X, mu):
    #   
    pred_X = np.array(pred_X)

    #   ,         
    answers_1 = []
    answers_2 = []

    for i in range(pred_X.shape[0]):
        if pred_X[i] == 1:
            answers_1.append(X[i])
        else:
            answers_2.append(X[i])
    
    print ' "     "'
    plt.figure(figsize=(16, 6))  
    plt.plot(
        np.array(answers_1)[:,0], np.array(answers_1)[:,1], 'o', alpha = 0.7, color='sandybrown', label = 'Produced on machine #1')
    plt.plot(
        np.array(answers_2)[:,0], np.array(answers_2)[:,1], 'o', alpha = 0.45, color = 'darkblue', label = 'Produced on machine #2')
    plt.plot(mu[0][0], mu[0][1], 'o', markersize = 16, color = 'red', label = 'Mu 1')
    plt.plot(mu[1][0], mu[1][1], 'o',  markersize = 16, color = 'slateblue', label = 'Mu 2')
    plt.xlabel('Diameter')
    plt.ylabel('Weight')
    plt.legend()
    plt.show()
    
#          
def graph_02_true(X, mu):
    print ' "  "'
    plt.figure(figsize=(16, 6))  
    plt.plot(
        X[0:N1,0], X[0:N1,1], 'o', alpha = 0.7, color='sandybrown', label = 'Produced on machine #1')
    plt.plot(
        X[N1:N,0], X[N1:N,1], 'o', alpha = 0.45, color = 'darkblue', label = 'Produced on machine #2')
    plt.plot(mu[0][0], mu[0][1], 'o', markersize = 16, color = 'red', label = 'Mu 1')
    plt.plot(mu[1][0], mu[1][1], 'o',  markersize = 16, color = 'slateblue', label = 'Mu 2')
    plt.xlabel('Diameter')
    plt.ylabel('Weight')
    plt.legend()
    plt.show()

Especificar los parámetros

#  
k = 2

#      №1
N1 = 6000
#      №2
N2 = 4000
N = N1+N2

#   (  )
m = 2

#    №1
mu_1_1 = 64.
#    №1
mu_1_2 = 14.
#  
sigma_1_1 = 3.5
sigma_1_2 = 1.

#    №2
mu_2_1 = 52.
#    №2
mu_2_2 = 9.5
#  
sigma_2_1 = 2.
sigma_2_2 = 0.7

Nosotros calcularemos

X = np.zeros((N, m))
np.random.seed(seed=42)
#    ,   №1
X[:N1, 0] = np.random.normal(loc=mu_1_1, scale=sigma_1_1, size=N1)
X[:N1, 1] = np.random.normal(loc=mu_1_2, scale=sigma_1_2, size=N1)
#  ,   №2
X[N1:N, 0] = np.random.normal(loc=mu_2_1, scale=sigma_2_1, size=N2)
X[N1:N, 1] = np.random.normal(loc=mu_2_2, scale=sigma_2_2, size=N2)

#    (   ,    )
y = np.zeros((N))
y[:N1] = np.array((1))
y[N1:N] = np.array((2))

#           (  )
mu  = np.array(([mu_1_1, mu_1_2], [mu_2_1, mu_2_2]))

#        (  )
sigma = np.array(([sigma_1_1, 0.],[0., sigma_1_2], [sigma_2_1, 0.],[0., sigma_2_2]))
sigma = sigma.reshape(k, m, m)

#     ,    №1  №2
w = np.array([float(1./k), float(1./k)])

#   
pj_xi = gaus_func_02(k, m, X, mu, sigma)
proba_X = proba_func_02(pj_xi, w, m)

#   ,   ,        
limit_proba = 0.5

pred_X = pred_x_02(proba_X, limit_proba)
        
#    
print '   :', round(accuracy_score(y, pred_X),3)
print

graph_02_algorithm(pred_X, mu)

graph_02_true(X, mu)

Gráfico No. 2.1 “Distribución de productos de acuerdo con el algoritmo”

Gráfico No. 2.2 “Distribución verdadera de productos”

Por analogía con el ejemplo anterior, generamos 10,000 valores de acuerdo con los parámetros anteriores

μ

$\mu$ y

σ

$\sigma$ , escribió varias funciones para el funcionamiento de nuestro algoritmo y lo lanzó. La diferencia fundamental entre este código y el código del ejemplo anterior es que esta vez utilizamos expresiones matriciales para los cálculos.

Como resultado, nuestro algoritmo muestra la proporción de respuestas correctas igual a 0.998, que en realidad es bastante buena.

Los problemas, como vemos, son todos iguales: errores en los detalles que se produjeron en diferentes máquinas y al mismo tiempo tienen dimensiones y pesos similares.

El código está escrito para que el lector pueda sustituir cualquier valor de parámetro.

μ

$\mu$ y

σ

$\sigma$ y vea cómo funcionará el algoritmo: en qué casos la calidad se deteriorará, en qué mejorará. Lo principal es no olvidar estudiar el horario. Bueno, nos estamos moviendo a nuestro punto final al analizar el algoritmo EM, de hecho, el algoritmo EM mismo.

Conoce el algoritmo EM

Continuamos nuestro ejemplo con las piezas que llegaron al almacén. Pero esta vez solo sabremos que los productos se fabricaron en dos máquinas diferentes, hay 10,000 de ellos, cada parte tiene un diámetro y un tamaño, y no sabemos nada más. Pero la tarea no ha cambiado: nosotros, como antes, de toda la pila grande, elementos aleatoriamente mezclados, tendremos que determinar a qué máquina pertenece esta o aquella parte.

A primera vista, parece poco realista, pero de hecho tenemos en nuestras manos un poderoso conjunto de herramientas: la fórmula de Bayes y la fórmula de densidad de probabilidad de una variable aleatoria. Tomemos toda esta bondad.

Solución

¿Qué haremos? Como debería ser en el algoritmo EM, primero inicializamos los parámetros: la

probabilidad de la hipótesis es extraer la parte producida en la máquina No. 1 -

w_{1}

$w_1$ Determinamos la probabilidad de una hipótesis para extraer la pieza producida en la máquina No. 2 -

w_{2}

$w_2$ . Solo hay dos hipótesis, lo que significa que cada una de ellas en el primer paso será igual a 0.5.

Expectativa matemática de variables aleatorias

μ

$\mu$ definir de la siguiente manera. Mezclamos todos los productos usando la función aleatoria, dividimos el conjunto en dos partes iguales, para cada parte para cada parámetro (diámetro, peso) determinamos el valor promedio.

Desviación Estándar

σ

$\sigma$ toma lo que se llama desde el techo: ponlo igual a la unidad en todos los aspectos. Escribimos en el formato de la matriz de covarianza.

Estamos listos para hacer el primer paso E del algoritmo. Usando los parámetros inicializados de variables aleatorias, determinamos la probabilidad de que cada parte se asigne a la máquina No. 1 o a la máquina No. 2.

En realidad, de esta manera, tomamos el primer E-step.

Ahora depende del paso M. Todo es simple aquí. Después de determinar la probabilidad de que cada parte se produzca en una máquina en particular, podemos volver a calcular la probabilidad de cada hipótesis:

w_{1}

$w_1$ ,

w_{2}

$w_2$ y

μ

$\mu$ y

σ

$\sigma$ .

Haremos 15 de tales iteraciones, en dos pasos cada una.

Observamos el código.

Funciones utilizadas en el algoritmo EM

#   E-
def e_step(x, k, m, n, w, mu, sigma):
    #      i-     j-  
    pj_xi = []
    for j in range(k):
        det_sigma_j = np.linalg.det(sigma[j])
        factor_1 = 1 / (((2 * math.pi)**(k/2)) * ((det_sigma_j)**0.5))
        factor_2 = []
        for i in x:
            factor_2.append(math.e**float(
                -0.5 * np.matrix(i - mu[j]) * np.matrix(np.linalg.inv(sigma[j])) * np.matrix(i - mu[j]).T))
        pj_xi.append(factor_1 * np.array(factor_2))
    pj_xi = np.array(pj_xi)
    
    #     ,  i-    j- 
    pj_xi_w = []
    for j in range(k):
        pj_xi_w.append(pj_xi[j] * w[j])
    pj_xi_w = np.array(pj_xi_w)
    
    #     i-      
    sum_pj_xi_w = np.sum(pj_xi_w, axis = 0)
    
    #         
    proba_xi = []
    for j in range(k):
        proba_xi.append(pj_xi_w[j] / sum_pj_xi_w)
    
    return np.array(proba_xi)

#  ,    ,            ,
#       
def x_new(proba_xi):
    X1_new_ind = []
    X2_new_ind = []
    X_answers = []

    count = 0
    for x in proba_xi[0]:
        if x >= 0.5:
            X1_new_ind.append(count)
            X_answers.append(1)
        else:
            X2_new_ind.append(count)
            X_answers.append(2)
        count += 1
    #      ,    №1, №2   
    return X1_new_ind, X2_new_ind, X_answers


#   M-
def m_step(x, proba_xi,N):
    w_new = np.sum(proba_xi, axis = 1) / N
    
    #   
    mu_new = (np.array((np.matrix(proba_xi) * np.matrix(X))).T / np.sum(proba_xi, axis = 1)).T
    
    #  
    cov_new = []
    for mu in range(mu_new.shape[0]):
        X_cd = []
        X_cd_proba = []
        count = 0
        for x_i in x:
            cd = np.array(x_i - mu_new[mu])
            X_cd.append(cd)
            X_cd_proba.append(cd * proba_xi[mu][count])
            count += 1
        X_cd = np.array(X_cd)
        X_cd = X_cd.reshape(N, m)
        X_cd_proba = np.array(X_cd_proba)
        X_cd_proba = X_cd_proba.reshape(N, m)

        cov_new.append(np.matrix(X_cd.T) * np.matrix(X_cd_proba))
    cov_new = np.array((np.array(cov_new) / (np.sum(proba_xi, axis = 1)-1)))
    #             ,    
    #       ,        
    if cov_new[0][0][1] < 0:
        cov_new[0][0][1] = 0
    if cov_new[0][1][0] < 0:
        cov_new[0][1][0] = 0
    
    if cov_new[1][0][1] < 0:
        cov_new[1][0][1] = 0
    if cov_new[1][1][0] < 0:
        cov_new[1][1][0] = 0
    
    #   
    sigma_new = cov_new**0.5
    return w_new, mu_new, sigma_new

Especificar los parámetros

#   
#      №1 ( 1)
N1 = 6000
#      №2 ( 2)
N2 = 4000
#       
N = N1 + N2

#  
k = 2

#    №1
mu_samples_1_1 = 64.
#    №1
mu_samples_1_2 = 14.

#    №2
mu_samples_2_1 = 52.
#    №2
mu_samples_2_2 = 9.5

#      №1
sigma_samples_1_1 = 3.5
#      №1
sigma_samples_1_2 = 1.

#      №2
sigma_samples_2_1 = 2.
#      №2
sigma_samples_2_2 = 0.7

Nosotros calcularemos

X = np.zeros((N, 2))

np.random.seed(seed=42)
#    ,    №1
X[:N1, 0] = np.random.normal(loc=mu_samples_1_1, scale=sigma_samples_1_1, size=N1)
X[:N1, 1] = np.random.normal(loc=mu_samples_1_2, scale=sigma_samples_1_2, size=N1)
#    ,    №2
X[N1:N, 0] = np.random.normal(loc=mu_samples_2_1, scale=sigma_samples_2_1, size=N2)
X[N1:N, 1] = np.random.normal(loc=mu_samples_2_2, scale=sigma_samples_2_2, size=N2)

#   
m = X.shape[1]

#   
n = X.shape[0]

#        (   )
y = np.zeros((N))
y[:N1] = np.array((1))
y[N1:N] = np.array((2))

#     ,    №1  №2
w = np.array([float(1./k), float(1./k)])

np.random.seed(seed = None)
#        (   )
mu  = np.array(
    (np.mean(X[np.random.choice(n, n/k)], axis = 0), np.mean(X[np.random.choice(n, n/k)], axis = 0)))
# mu = np.array(([mu_samples_1_1, mu_samples_1_2],[mu_samples_2_1, mu_samples_2_2]))

#         (    )
sigma = np.array(([1., 0.],[0., 1.], [1., 0.],[0., 1.]))
# sigma = np.array(([sigma_samples_1_1, 0.],[0., sigma_samples_1_2], [sigma_samples_2_1, 0.],[0., sigma_samples_2_2]))
sigma = sigma.reshape(k, m, m)

#    EM-
steps = 15
#   EM-
for i in range(steps):
    proba_xi = e_step(X, k, m, n, w, mu, sigma)
    w, mu, sigma = m_step(X, proba_xi,N)
    X1_new_ind, X2_new_ind, X_answers = x_new(proba_xi)
    print ' №', i+1
    print
    print '   '
    print mu
    print
    print '   '
    print sigma
    print
    print '   '
    print round(accuracy_score(y, X_answers),3)
    
    plt.figure(figsize=(16, 6))  
    plt.plot(
        X[X1_new_ind,0], X[X1_new_ind,1], 'o', alpha = 0.7, color='sandybrown', label = 'Produced on machine #1')
    plt.plot(
        X[X2_new_ind,0], X[X2_new_ind,1], 'o', alpha = 0.45, color = 'darkblue', label = 'Produced on machine #2')
    plt.plot(mu[0][0], mu[0][1], 'o', markersize = 16, color = 'red', label = 'Mu 1')
    plt.plot(mu[1][0], mu[1][1], 'o',  markersize = 16, color = 'purple', label = 'Mu 2')
    plt.xlabel('Diameter')
    plt.ylabel('Weight')
    plt.legend()
    plt.show()

Gráfico 3.1 “Distribución de partes en cada iteración del algoritmo EM”

Todas las iteraciones

¿Que pasó al final? Y resultó lo siguiente, al no tener información preliminar sobre la probabilidad de hipótesis para la extracción de piezas producidas en la máquina No. 1:

w_{1}

$w_1$ o máquina número 2

w_{2}

$w_2$ , expectativa matemática y desviación estándar de variables aleatorias -

μ

$\mu$ y

σ

$\sigma$ , logramos una calidad comparable al ejemplo anterior, en el que se conocían todos los parámetros anteriores.

Compare los parámetros que establecemos al generar el conjunto de productos, los llamados parámetros verdaderos y los parámetros que el algoritmo EM determinó para nosotros:

parámetros verdaderos

μ_{11} = 64, μ_{12} = 14, σ_{11} = 3.5, σ_{12} = 1.

$\mu_{11} = 64, \mu_{12} = 14, \sigma_{11} = 3.5, \sigma_{12} = 1.$

μ_{21} = 52, μ_{12} = 9.5, σ_{21} = 2.0, σ_{22} = 0.7

$\mu_{21} = 52, \mu_{12} = 9.5, \sigma_{21} = 2.0, \sigma_{22} = 0.7$

Parámetros específicos

μ_{11} \approx 64.03, μ_{12} \approx 13.99, σ_{11} \approx 3.44, σ_{12} \approx 1.23

$\mu_{11} \approx 64.03, \mu_{12} \approx 13.99, \sigma_{11} \approx 3.44, \sigma_{12} \approx 1.23$

μ_{21} \approx 52.06, μ_{12} \approx 9.54, σ_{21} \approx 1.65, σ_{22} \approx 0.77

$\mu_{21} \approx 52.06, \mu_{12} \approx 9.54, \sigma_{21} \approx 1.65, \sigma_{22} \approx 0.77$

En este artículo llega a su fin. En conclusión, agregaré que para el análisis de datos, no tiene sentido reinventar la rueda y escribir el algoritmo EM por su cuenta, simplemente puede usar la función de biblioteca sklearn proporcionada.

Veamos cómo la biblioteca sklearn GaussianMixture funciona con nuestro ejemplo.

Módulo de importación

from sklearn.mixture import GaussianMixture as GMM

Especificar los parámetros

#   
#      №1 ( 1)
N1 = 6000
#      №2 ( 2)
N2 = 4000
#       
N = N1 + N2

#    №1
mu_samples_1_1 = 64.
#    №1
mu_samples_1_2 = 14.

#    №2
mu_samples_2_1 = 52.
#    №2
mu_samples_2_2 = 9.5

#      №1
sigma_samples_1_1 = 3.5
#      №1
sigma_samples_1_2 = 1.

#      №2
sigma_samples_2_1 = 2.
#      №2
sigma_samples_2_2 = 0.7


#  
k = 2

X = np.zeros((N, 2))

np.random.seed(seed=42)
#    ,    №1
X[:N1, 0] = np.random.normal(loc=mu_samples_1_1, scale=sigma_samples_1_1, size=N1)
X[:N1, 1] = np.random.normal(loc=mu_samples_1_2, scale=sigma_samples_1_2, size=N1)
#    ,    №2
X[N1:N, 0] = np.random.normal(loc=mu_samples_2_1, scale=sigma_samples_2_1, size=N2)
X[N1:N, 1] = np.random.normal(loc=mu_samples_2_2, scale=sigma_samples_2_2, size=N2)

#        (   )
y = np.zeros((N))
y[:N1] = np.array((1))
y[N1:N] = np.array((2))

Ejecute el algoritmo EM

np.random.seed(1)
model = GMM(n_components=k, covariance_type='full')
model.fit(X)


temp_predict_X = model.predict(X)
X_answers = []
for i in range(X.shape[0]):
    if temp_predict_X[i] == 0:
        X_answers.append(1)
    else:
        X_answers.append(2)
        

print '   '
print round(accuracy_score(y, X_answers),3)

Trabajado con una explosión! ¡Sin errores! Como era de esperar, la biblioteca GaussianMixture de skussian es más rápida y mejor que nuestro algoritmo. Sin embargo, el propósito del artículo, como recordamos, era comprender cómo funcionan los algoritmos EM. Espero que se haya logrado el objetivo :)

Fuentes de información

Literatura

1) Estadística aplicada: clasificaciones y reducción de dimensionalidad, S.A. Ayvazyan, V.M. Buchstaber, I.S. Enyukov, L.D. Meshalkin, Moscú, Finanzas y Estadísticas, 1989

Recursos en Internet

1) Artículo 'Un tutorial apacible del algoritmo EM y su aplicación a la estimación de parámetros para la mezcla gaussiana y los modelos ocultos de Markov', Instituto Internacional de Informática, Jeff A. Bilmes, 1998

2) Presentación 'Modelos de mezcla y algoritmo EM', Universidad de Nueva York, David Sontag

3) Curso 'Búsqueda de estructuras de datos'

4) Fórmula para probabilidad completa y fórmulas de Bayes

5) Función 'sklearn.mixture.GaussianMixture'

6) Distribución normal, Wikipedia

7) Multidimensional distribución normal, Wikipedia

8)github.com/diefimov/MTH594_MachineLearning/blob/master/ipython/Lecture10.ipynb

9) Un artículo sobre el algoritmo EM

10) Otro artículo sobre el algoritmo EM

11) Github github.com/AlexanderPetrenko83/Articles

Desmontamos el algoritmo EM en pequeños ladrillos

Algoritmo EM en general

Recordemos la fórmula de Bayes

Ejemplo uno por el teorema de Bayes

Segundo ejemplo para el teorema de Bayes usando los parámetros de la distribución normal de una variable aleatoria: expectativa matemática μ\muy desviación estándar σ\sigma

El tercer ejemplo. Caja bidimensional

Conoce el algoritmo EM

Fuentes de información

More articles:

Segundo ejemplo para el teorema de Bayes usando los parámetros de la distribución normal de una variable aleatoria: expectativa matemática $\mu$ y desviación estándar $\sigma$