Una red que reconoce dígitos, escrita entera por nosotras

Trece capítulos diciendo dónde la red pierde. Este es donde gana 🎉

Vamos a escribir una red completa, de principio a fin, sin importar ni una librería de deep learning. Y va a reconocer números escritos a mano.

El encargo

1.797 imágenes de dígitos de 8 por 8, las mismas del capítulo 9. Hay que decir qué número es cada una, de 0 a 9.

Y ojo con eso último, porque es nuevo: hasta ahora todo el libro fue sí o no. Aquí hay diez respuestas posibles, y eso cambia dos piezas.

import numpy as np
from sklearn.datasets import load_digits
from sklearn.linear_model import LogisticRegression
from sklearn.metrics import accuracy_score, confusion_matrix
from sklearn.model_selection import train_test_split

digitos = load_digits()
X = digitos.images.reshape(len(digitos.images), -1) / 16.0
y = digitos.target

i_tr, i_te = train_test_split(np.arange(len(X)), test_size=0.25,
                              random_state=42, stratify=y)
X_tr, X_te = X[i_tr], X[i_te]
y_tr, y_te = y[i_tr], y[i_te]

Y = np.zeros((len(y_tr), 10))              # el objetivo, una columna por dígito
Y[np.arange(len(y_tr)), y_tr] = 1

print('entreno con', len(X_tr), 'imágenes de', X.shape[1], 'píxeles')
print('el objetivo pasa de', y_tr.shape, 'a', Y.shape)
print('la primera imagen es un', y_tr[0], 'y su fila objetivo es', Y[0].astype(int))

entreno con 1347 imágenes de 64 píxeles
el objetivo pasa de (1347,) a (1347, 10)
la primera imagen es un 2 y su fila objetivo es [0 0 1 0 0 0 0 0 0 0]

Ese Y con una columna por dígito es el one-hot del capítulo 10, ahora aplicado a la respuesta en vez de a las entradas. La red va a sacar diez números y compararlos con esa fila 🔟

Las dos piezas nuevas

Para diez clases cambia la activación de salida y cambia la pérdida:

def softmax(z):
    z = z - z.max(axis=-1, keepdims=True)
    e = np.exp(z)
    return e / e.sum(axis=-1, keepdims=True)

def perdida(P, Y):
    return float(-np.mean(np.sum(Y * np.log(np.clip(P, 1e-12, 1)), axis=1)))

ejemplo = softmax(np.array([[2.0, 0.5, -1.0]]))
print('tres puntajes se vuelven:', np.round(ejemplo, 4), ' suman', round(float(ejemplo.sum()), 4))

tres puntajes se vuelven: [[0.7856 0.1753 0.0391]]  suman 1.0

La softmax del capítulo 3 en vez de la sigmoide, para que las diez salidas sumen 1 y se puedan leer como probabilidades que compiten.

Y la entropía cruzada categórica, que es la del capítulo 5 generalizada: en vez de mirar una probabilidad, mira la de la clase correcta.

Y aquí está lo bonito: aunque las dos fórmulas cambien, P - Y sigue siendo el gradiente, exactamente igual que con una sola salida. La pareja softmax y entropía cruzada se simplifica igual de bien que sigmoide y entropía cruzada. Por eso van siempre juntas 🎯

La red entera

def entrena(ocultas=64, vueltas=1000, paso=0.5, semilla=0, decay=0.0):
    r = np.random.default_rng(semilla)
    entradas = X_tr.shape[1]
    # arranque con la escala de He, que es la que va con ReLU
    W1 = r.normal(0, np.sqrt(2 / entradas), (entradas, ocultas))
    b1 = np.zeros(ocultas)
    W2 = r.normal(0, np.sqrt(2 / ocultas), (ocultas, 10))
    b2 = np.zeros(10)

    for _ in range(vueltas):
        h = np.maximum(0, X_tr @ W1 + b1)          # ReLU, capítulo 3
        P = softmax(h @ W2 + b2)
        d2 = (P - Y) / len(Y)                      # capítulo 5
        d1 = (d2 @ W2.T) * (h > 0)                 # capítulo 6
        W2 -= paso * (h.T @ d2 + decay * W2)       # capítulo 8
        b2 -= paso * d2.sum(axis=0)
        W1 -= paso * (X_tr.T @ d1 + decay * W1)
        b1 -= paso * d1.sum(axis=0)
    return W1, b1, W2, b2

def predice(Z, W1, b1, W2, b2):
    return softmax(np.maximum(0, Z @ W1 + b1) @ W2 + b2)

W1, b1, W2, b2 = entrena()
p_tr = predice(X_tr, W1, b1, W2, b2).argmax(axis=1)
p_te = predice(X_te, W1, b1, W2, b2).argmax(axis=1)
print('acierto entrenando:', round(accuracy_score(y_tr, p_tr), 4))
print('acierto en prueba :', round(accuracy_score(y_te, p_te), 4))

acierto entrenando: 1.0
acierto en prueba : 0.9756

Ahí está, y con eso ya terminó el libro en lo esencial 💛

Una red de una capa oculta, escrita en veinte líneas, reconoce números escritos a mano con casi 98% de acierto. Sin TensorFlow, sin PyTorch, sin nada.

Fíjate en el arranque de los pesos, que es lo único que no habíamos visto: sqrt(2 / entradas). Se llama inicialización de He y está pensada para ReLU. Si arrancas con la escala 0,1 de los capítulos anteriores, con 64 entradas por neurona la suma sale demasiado chica y la red tarda mucho más en arrancar.

Y ahora sí: contra la logística

logistica = LogisticRegression(max_iter=3000).fit(X_tr, y_tr)
print('logística:', round(accuracy_score(y_te, logistica.predict(X_te)), 4))
print('nuestra red:', round(accuracy_score(y_te, p_te), 4))

logística: 0.9622
nuestra red: 0.9756

0,9756 la red contra 0,9622 la logística 🏆

Después de trece capítulos perdiendo, gana. Y gana en el sitio donde el capítulo 1 dijo que ganaría: cuando el dato tiene estructura.

Es un punto y tres décimas de diferencia, que parece poco y no lo es: sobre 450 imágenes son seis errores menos. En un sistema que procesa cheques o formularios, eso es la diferencia entre revisar seis cosas a mano o no.

Y ahora fíjate en algo que en el capítulo 7 era la señal de alarma: el acierto en entrenamiento es 1,0000, o sea perfecto. Allá eso venía con una prueba de 0,5870 y significaba memorización pura.

Aquí viene con 0,9756 en prueba. La red se aprendió las 1.347 imágenes y además generaliza.

Eso es lo que cambia cuando el dato tiene estructura: no es que memorizar deje de ocurrir, es que lo que aprende de memoria también sirve para las imágenes nuevas. Por eso aquí no hicieron falta ni parada temprana ni frenos 🧠

En qué se equivoca

M = confusion_matrix(y_te, p_te)
print('errores:', int((p_te != y_te).sum()), 'de', len(y_te))
print()
for i in range(10):
    for j in range(10):
        if i != j and M[i, j] >= 2:
            print(f'   confundió el {i} con el {j}: {M[i, j]} veces')

errores: 11 de 450

   confundió el 0 con el 4: 2 veces
   confundió el 8 con el 1: 3 veces

Once errores en 450 imágenes. Y los que repite son el 8 con el 1 (tres veces) y el 0 con el 4 (dos veces).

El 8 con el 1 tiene sentido: en 8 por 8 píxeles, un 8 mal escrito es un trazo vertical con dos bultitos, y un 1 es un trazo vertical. El 0 con el 4 es más raro, y es lo que yo iría a mirar imagen por imagen antes de entregar esto 🔍

¿Sabe cuándo duda?

La pregunta del capítulo 13 del libro de machine learning, ahora aquí:

probabilidades = predice(X_te, W1, b1, W2, b2)
confianza = probabilidades.max(axis=1)
fallo = p_te != y_te

print('confianza media cuando acierta:', round(float(confianza[~fallo].mean()), 4))
print('confianza media cuando falla  :', round(float(confianza[fallo].mean()), 4))
print()
for corte in [0.5, 0.9, 0.99]:
    seguros = confianza >= corte
    print(f'si solo contesto con confianza >= {corte}: '
          f'contesto {seguros.mean():.4f} de las veces y acierto '
          f'{accuracy_score(y_te[seguros], p_te[seguros]):.4f}')

confianza media cuando acierta: 0.9742
confianza media cuando falla  : 0.7159

si solo contesto con confianza >= 0.5: contesto 0.9844 de las veces y acierto 0.9819
si solo contesto con confianza >= 0.9: contesto 0.9178 de las veces y acierto 0.9952
si solo contesto con confianza >= 0.99: contesto 0.7867 de las veces y acierto 0.9972

Sí sabe 🎯

Cuando acierta, su confianza media es 0,9742. Cuando falla, 0,7159. Son veintiséis puntos de diferencia, y eso convierte la confianza en una señal utilizable.

Y la tabla de abajo es el producto de verdad: pidiéndole confianza de al menos 0,9 contesta el 91,78% de las veces y en esas acierta el 99,52%. El 8% restante se manda a una persona.

Ese es exactamente el diseño de un sistema de lectura de formularios que funciona en producción. No uno que acierta siempre: uno que sabe cuáles pasar a revisión 🤝

Y esto para qué le sirve a una distribuidora

Porque el libro empezó con las ventas de una distribuidora peruana y termina con dígitos, y conviene cerrar el círculo 🔄

En el capítulo 1 medimos que sobre esa tabla de 3.000 ventas la red sacaba 0,5834 y la regresión logística 0,7214. Aquí sobre imágenes la red saca 0,9756 y la logística 0,9622. El mismo tipo de modelo, resultados invertidos, y la diferencia no está en el modelo: está en si el dato tiene forma.

Así que la pregunta útil no es "¿usamos deep learning?", es "¿qué tenemos que no entre en una tabla?". En una distribuidora eso suele existir y suele estar tirado:

• 📝 Pedidos anotados a mano por el vendedor, que alguien transcribe.
• 🧾 Facturas y guías escaneadas de las que hay que sacar números.
• 📸 Fotos de anaquel para ver si el producto está donde debería.
• 💬 Comentarios y reclamos de clientes en texto libre.

Los tres primeros son este capítulo con imágenes más grandes. El cuarto son los capítulos 10, 11 y 13.

Y para lo que sí entra en una tabla, o sea a quién visitar y quién va a comprar, la respuesta sigue siendo el libro anterior: una regresión logística, que se entrena en dos décimas de segundo y se explica leyendo sus coeficientes 💛

Lo que se entrega

Como en el proyecto final del libro de machine learning, esto no se entrega en un cuaderno:

• 💾 Los cuatro arrays (W1, b1, W2, b2) guardados con np.savez. Pesan poco y no dependen de ninguna librería.
• 📄 Una página con qué predice, con qué datos, cuánto acierta (0,9756), contra qué se compara (0,9622) y en qué se equivoca (8 con 1, 0 con 4).
• 🚦 El umbral de confianza elegido y qué porcentaje manda a revisión.
• 📐 La lista de tamaños de entrada, porque el capítulo 9 nos enseñó que ahí es donde revientan las cosas.

Ejercicios

import os
import tempfile

ruta = os.path.join(tempfile.gettempdir(), 'red-digitos.npz')
np.savez(ruta, W1=W1, b1=b1, W2=W2, b2=b2)
print('pesa', round(os.path.getsize(ruta) / 1024, 1), 'KB')

cargado = np.load(ruta)
p2 = predice(X_te, cargado['W1'], cargado['b1'],
             cargado['W2'], cargado['b2']).argmax(axis=1)
print('predice igual:', bool((p2 == p_te).all()))

pesa 38.5 KB
predice igual: True

for ocultas in [16, 64, 128]:
    ws = entrena(ocultas=ocultas)
    pr = predice(X_te, *ws).argmax(axis=1)
    print(f'{ocultas:3d} neuronas: {accuracy_score(y_te, pr):.4f}')
print()
for vueltas in [200, 1000, 3000]:
    ws = entrena(vueltas=vueltas)
    pr = predice(X_te, *ws).argmax(axis=1)
    print(f'{vueltas:5d} vueltas: {accuracy_score(y_te, pr):.4f}')

 16 neuronas: 0.9556
 64 neuronas: 0.9756
128 neuronas: 0.9711

  200 vueltas: 0.9622
 1000 vueltas: 0.9756
 3000 vueltas: 0.9711

def entrena_escala(escala, vueltas=1000, ocultas=64):
    r = np.random.default_rng(0)
    W1 = r.normal(0, escala, (X_tr.shape[1], ocultas)); b1 = np.zeros(ocultas)
    W2 = r.normal(0, escala, (ocultas, 10)); b2 = np.zeros(10)
    for _ in range(vueltas):
        h = np.maximum(0, X_tr @ W1 + b1)
        P = softmax(h @ W2 + b2)
        d2 = (P - Y) / len(Y)
        d1 = (d2 @ W2.T) * (h > 0)
        W2 -= 0.5 * (h.T @ d2); b2 -= 0.5 * d2.sum(axis=0)
        W1 -= 0.5 * (X_tr.T @ d1); b1 -= 0.5 * d1.sum(axis=0)
    return accuracy_score(y_te, softmax(np.maximum(0, X_te @ W1 + b1) @ W2 + b2).argmax(axis=1))

for escala in [0.01, 0.1, float(np.sqrt(2 / 64))]:
    print(f'escala {escala:.4f}: {entrena_escala(escala):.4f}')

escala 0.0100: 0.9644
escala 0.1000: 0.9689
escala 0.1768: 0.9756

fallos = np.where(p_te != y_te)[0]
for k in fallos[:2]:
    imagen = (X_te[k].reshape(8, 8) * 16).astype(int)
    print(f'era un {y_te[k]} y dijo {p_te[k]} con confianza '
          f'{confianza[k]:.3f}')
    for fila in imagen:
        print('   ' + ''.join('#' if v > 8 else ('+' if v > 3 else '.') for v in fila))
    print()

era un 8 y dijo 1 con confianza 0.880
   ..+##+..
   ..+###..
   ..+##...
   ..+#+...
   ..+#....
   ..+#+...
   ..###...
   ..+##...

era un 0 y dijo 4 con confianza 0.839
   ...##...
   ...##+..
   ..##.#..
   ..#+.#..
   ..##.#..
   ..#+##..
   ..###+..
   ...##...

print(f'{"corte":>6} {"contesta":>10} {"acierta":>9} {"revisan":>9} {"errores":>9}')
for corte in [0.0, 0.8, 0.95, 0.99, 0.999]:
    seguros = confianza >= corte
    if seguros.sum() == 0:
        continue
    acierto = accuracy_score(y_te[seguros], p_te[seguros])
    print(f'{corte:6.3f} {seguros.mean():10.4f} {acierto:9.4f} '
          f'{(~seguros).sum():9d} {int((p_te[seguros] != y_te[seguros]).sum()):9d}')

 corte   contesta   acierta   revisan   errores
 0.000     1.0000    0.9756         0        11
 0.800     0.9511    0.9860        22         6
 0.950     0.8756    0.9949        56         2
 0.990     0.7867    0.9972        96         1
 0.999     0.5644    1.0000       196         0

PIEZAS = [
    ('la neurona: multiplicar, sumar, aplastar', 2),
    ('ReLU en la capa oculta', 3),
    ('la capa oculta que dobla', 4),
    ('la pérdida y el descenso de gradiente', 5),
    ('la retropropagación con (P - Y)', 6),
    ('mirar entrenamiento y prueba a la vez', 7),
    ('el weight decay, aquí sin necesitarlo', 8),
    ('los píxeles como dato con forma', 9),
    ('el one-hot del objetivo', 10),
    ('la softmax estable', 12),
]
for pieza, capitulo in PIEZAS:
    print(f'capítulo {capitulo:2d}: {pieza}')

capítulo  2: la neurona: multiplicar, sumar, aplastar
capítulo  3: ReLU en la capa oculta
capítulo  4: la capa oculta que dobla
capítulo  5: la pérdida y el descenso de gradiente
capítulo  6: la retropropagación con (P - Y)
capítulo  7: mirar entrenamiento y prueba a la vez
capítulo  8: el weight decay, aquí sin necesitarlo
capítulo  9: los píxeles como dato con forma
capítulo 10: el one-hot del objetivo
capítulo 12: la softmax estable

P = predice(X_te, W1, b1, W2, b2)
P - y_te

ValueError: operands could not be broadcast together with shapes (450,10) (450,) 

Lo que te llevas

• 🏆 Una red escrita entera en numpy clasifica dígitos al 97,56% contra el 96,22% de una regresión logística. Por fin gana, y gana donde el capítulo 1 dijo.
• 🔟 Para diez clases cambian dos piezas: softmax en la salida y entropía cruzada categórica. Y el gradiente sigue siendo P - Y.
• 🧠 El acierto en entrenamiento es 1,0000 y aquí eso no es alarma: con estructura de verdad, lo que aprende de memoria también sirve para lo nuevo.
• 🎬 La inicialización de He vale un punto entero de acierto, y es de lo que menos se habla.
• 🎯 Su confianza media es 0,9742 cuando acierta y 0,7159 cuando falla: el modelo sabe cuándo duda.
• 🤝 Con confianza mínima de 0,9 contesta el 91,78% de las veces y acierta el 99,52%. Y con 0,999 no pasa ni un error, revisando 196 de 450.
• 💾 Se entrega en cuatro arrays de numpy, sin dependencias que se rompan.

En el capítulo 15 te dejo dónde seguir, qué librerías tocan ahora y qué cosas este libro no cubre.

Que tengas lindo día! 🌸