Redes convolucionales

import numpy as np
import matplotlib.pyplot as plt
from matplotlib.colors import LinearSegmentedColormap
import torch
from torch import nn
from torchvision import datasets

torch.manual_seed(0)
device = "cuda" if torch.cuda.is_available() else "cpu"

# La primera vez descarga unos 30 MB; después queda guardado en «datos»
entrenamiento = datasets.FashionMNIST(root="datos", train=True,
                                      download=True)
prueba = datasets.FashionMNIST(root="datos", train=False, download=True)
nombres = ["camiseta", "pantalón", "suéter", "vestido", "abrigo",
           "sandalia", "camisa", "zapatilla", "bolso", "botín"]

print(entrenamiento.data.shape, entrenamiento.data.dtype)
prenda = entrenamiento.data[0].numpy() / 255   # números entre 0 y 1
print("Etiqueta:", nombres[entrenamiento.targets[0]])
# Filas 6 a 10, columnas 9 a 16: un pedazo del borde de la prenda
print(np.round(prenda[6:11, 9:17], 1))

def convolucion2d(imagen, filtro, padding=0, stride=1):
    """Desliza `filtro` sobre `imagen` y suma los productos."""
    if padding > 0:
        # Un marco de ceros alrededor de la imagen
        imagen = np.pad(imagen, padding)
    k = filtro.shape[0]
    alto = (imagen.shape[0] - k) // stride + 1
    ancho = (imagen.shape[1] - k) // stride + 1
    salida = np.zeros((alto, ancho))
    for i in range(alto):
        for j in range(ancho):
            fila, col = i * stride, j * stride
            ventana = imagen[fila:fila + k, col:col + k]
            salida[i, j] = np.sum(ventana * filtro)
    return salida

# Filtro Sobel horizontal: resta la izquierda de la derecha
sobel_x = np.array([[-1, 0, 1],
                    [-2, 0, 2],
                    [-1, 0, 1]])
# Una imagen chiquita: oscura a la izquierda, clara a la derecha
mini = np.array([[0, 0, 0, 9, 9],
                 [0, 0, 0, 9, 9],
                 [0, 0, 0, 9, 9],
                 [0, 0, 0, 9, 9],
                 [0, 0, 0, 9, 9]])
print(convolucion2d(mini, sobel_x))

from sklearn.datasets import load_sample_image

foto = load_sample_image("china.jpg")
gris = foto @ np.array([0.299, 0.587, 0.114])
# Uno de cada dos píxeles: menos ventanas que recorrer
gris = gris[::2, ::2] / 255
sobel_y = sobel_x.T

bordes_x = convolucion2d(gris, sobel_x)
bordes_y = convolucion2d(gris, sobel_y)
# Los dos bordes juntos, como la longitud de un vector
magnitud = np.sqrt(bordes_x ** 2 + bordes_y ** 2)
print("Foto en gris:", gris.shape, "  Mapa de bordes:", bordes_x.shape)
print(f"Valores de bordes_x: de {bordes_x.min():.2f} "
      f"a {bordes_x.max():.2f}")

# Blanco para el cero y tinta para los valores más altos
blanco_tinta = LinearSegmentedColormap.from_list("bt", ["white", "C0"])
fig, ejes = plt.subplots(2, 2, figsize=(4.8, 4.8))
ejes[0, 0].imshow(foto)
ejes[0, 0].set_title("Original")
# Valor absoluto: un borde es borde, suba o baje el brillo
ejes[0, 1].imshow(np.abs(bordes_x), cmap=blanco_tinta, vmax=2)
ejes[0, 1].set_title("Bordes verticales")
ejes[1, 0].imshow(np.abs(bordes_y), cmap=blanco_tinta, vmax=2)
ejes[1, 0].set_title("Bordes horizontales")
ejes[1, 1].imshow(magnitud, cmap=blanco_tinta, vmax=2)
ejes[1, 1].set_title("Los dos juntos")
for eje in ejes.flat:
    eje.axis("off")
plt.show()

for p, s in [(0, 1), (1, 1), (1, 2), (0, 2)]:
    salida = convolucion2d(prenda, sobel_x, padding=p, stride=s)
    formula = (28 + 2 * p - 3) // s + 1
    print(f"padding {p}, stride {s}: salida {salida.shape}, "
          f"fórmula {formula}")

def pooling2d(mapa, tam=2, modo="max"):
    """Resume cada bloque de tam × tam con su máximo o su promedio."""
    alto, ancho = mapa.shape[0] // tam, mapa.shape[1] // tam
    recorte = mapa[:alto * tam, :ancho * tam]
    # Ejes: (fila de bloque, fila dentro, columna de bloque, col. dentro)
    bloques = recorte.reshape(alto, tam, ancho, tam)
    if modo == "max":
        return bloques.max(axis=(1, 3))
    return bloques.mean(axis=(1, 3))

# Un punto brillante y el mismo punto movido un píxel
punto = np.zeros((4, 4))
punto[0, 0] = 1
movido = np.zeros((4, 4))
movido[1, 1] = 1
print("Iguales antes del pooling:", np.array_equal(punto, movido))
print("Iguales después del max-pooling:",
      np.array_equal(pooling2d(punto), pooling2d(movido)))
print(pooling2d(movido))

bordes_prenda = np.abs(convolucion2d(prenda, sobel_x))
maximo = pooling2d(bordes_prenda, modo="max")
promedio = pooling2d(bordes_prenda, modo="promedio")
print(prenda.shape, bordes_prenda.shape, maximo.shape, promedio.shape)

# Los tres mapas de bordes con la misma escala, para poder compararlos
tope = bordes_prenda.max()
fig, ejes = plt.subplots(2, 2, figsize=(4.8, 4.8))
paneles = [(prenda, "Prenda 28 × 28", 1),
           (bordes_prenda, "Bordes 26 × 26", tope),
           (maximo, "Máximo 13 × 13", tope),
           (promedio, "Promedio 13 × 13", tope)]
for eje, (mapa, titulo, vmax) in zip(ejes.flat, paneles):
    eje.imshow(mapa, cmap=blanco_tinta, vmin=0, vmax=vmax)
    eje.set_title(titulo)
    eje.axis("off")
plt.show()

import torch.nn.functional as F

imagen_t = torch.tensor(prenda, dtype=torch.float32).reshape(1, 1, 28, 28)
filtro_t = torch.tensor(sobel_x, dtype=torch.float32).reshape(1, 1, 3, 3)
salida_t = F.conv2d(imagen_t, filtro_t)
print(salida_t.shape)
print("¿Igual a la nuestra?",
      np.allclose(salida_t[0, 0].numpy(), convolucion2d(prenda, sobel_x),
                  atol=1e-5))

capa = nn.Conv2d(in_channels=1, out_channels=16, kernel_size=3,
                 padding=1)
lote = torch.zeros(64, 1, 28, 28)            # (N, C, H, W)
print("Pesos de la capa:", tuple(capa.weight.shape))
print("Después de la convolución:", tuple(capa(lote).shape))
print("Después de MaxPool2d(2):  ",
      tuple(nn.MaxPool2d(2)(capa(lote)).shape))
print("Después de AvgPool2d(2):  ",
      tuple(nn.AvgPool2d(2)(capa(lote)).shape))
print("Después de AdaptiveAvgPool2d(1):",
      tuple(nn.AdaptiveAvgPool2d(1)(capa(lote)).shape))

def contar(modelo):
    """Cantidad total de números entrenables de un modelo."""
    return sum(p.numel() for p in modelo.parameters())

# Una capa densa que diera la misma salida: 16 × 28 × 28 números
densa = nn.Linear(28 * 28, 16 * 28 * 28)
print(f"Parámetros de la convolución: {contar(capa)}")
print(f"Parámetros de la capa densa:  {contar(densa)}")

from torch.utils.data import DataLoader, TensorDataset

g = torch.Generator().manual_seed(0)
idx_ent = torch.randperm(len(entrenamiento), generator=g)[:10000]
idx_val = torch.randperm(len(prueba), generator=g)[:2000]

x_ent = entrenamiento.data[idx_ent].float() / 255
x_val = prueba.data[idx_val].float() / 255
y_ent = entrenamiento.targets[idx_ent]
y_val = prueba.targets[idx_val]

# Media y desviación del entrenamiento, para las dos partes
media, desv = x_ent.mean().item(), x_ent.std().item()
print(f"Media de los píxeles: {media:.3f}   desviación: {desv:.3f}")
# unsqueeze(1) agrega el eje del canal: (N, 28, 28) -> (N, 1, 28, 28)
x_ent = ((x_ent - media) / desv).unsqueeze(1)
x_val = ((x_val - media) / desv).unsqueeze(1)
print("Entrenamiento:", tuple(x_ent.shape), " Validación:",
      tuple(x_val.shape))

cargador_ent = DataLoader(TensorDataset(x_ent, y_ent), batch_size=64,
                          shuffle=True,
                          generator=torch.Generator().manual_seed(0))
cargador_val = DataLoader(TensorDataset(x_val, y_val), batch_size=500)

import torch
from torch import nn


def evaluar(modelo, cargador, funcion_perdida, dispositivo="cpu"):
    """Pérdida media y exactitud de `modelo` sobre `cargador`."""
    modelo.eval()
    suma, aciertos, total = 0.0, 0, 0
    with torch.no_grad():
        for x, y in cargador:
            x, y = x.to(dispositivo), y.to(dispositivo)
            salida = modelo(x)
            suma += funcion_perdida(salida, y).item() * len(y)
            aciertos += (salida.argmax(dim=1) == y).sum().item()
            total += len(y)
    return suma / total, aciertos / total


def entrenar(modelo, cargador_ent, cargador_val, epocas, optimizador,
             funcion_perdida=None, planificador=None, dispositivo="cpu",
             mostrar=True):
    """Entrena `modelo` y devuelve el historial de cada época."""
    if funcion_perdida is None:
        funcion_perdida = nn.CrossEntropyLoss()
    modelo.to(dispositivo)
    historial = {"perdida_ent": [], "exactitud_ent": [],
                 "perdida_val": [], "exactitud_val": []}
    for epoca in range(1, epocas + 1):
        modelo.train()
        suma, aciertos, total = 0.0, 0, 0
        for x, y in cargador_ent:
            x, y = x.to(dispositivo), y.to(dispositivo)
            optimizador.zero_grad()
            salida = modelo(x)
            perdida = funcion_perdida(salida, y)
            perdida.backward()
            optimizador.step()
            suma += perdida.item() * len(y)
            aciertos += (salida.argmax(dim=1) == y).sum().item()
            total += len(y)
        if planificador is not None:
            planificador.step()
        perdida_val, exactitud_val = evaluar(
            modelo, cargador_val, funcion_perdida, dispositivo)
        historial["perdida_ent"].append(suma / total)
        historial["exactitud_ent"].append(aciertos / total)
        historial["perdida_val"].append(perdida_val)
        historial["exactitud_val"].append(exactitud_val)
        if mostrar:
            print(f"Época {epoca}: pérdida {suma / total:.3f}, "
                  f"validación {perdida_val:.3f}, "
                  f"exactitud de validación {exactitud_val:.3f}")
    return historial

torch.manual_seed(0)
cnn = nn.Sequential(
    nn.Conv2d(1, 16, kernel_size=3, padding=1),   # (16, 28, 28)
    nn.ReLU(),
    nn.MaxPool2d(2),                              # (16, 14, 14)
    nn.Conv2d(16, 32, kernel_size=3, padding=1),  # (32, 14, 14)
    nn.ReLU(),
    nn.MaxPool2d(2),                              # (32, 7, 7)
    nn.Flatten(),                                 # 1568
    nn.Linear(32 * 7 * 7, 10),
)
torch.manual_seed(0)
mlp = nn.Sequential(
    nn.Flatten(),                                 # 784
    nn.Linear(28 * 28, 26),
    nn.ReLU(),
    nn.Linear(26, 10),
)
print(f"Parámetros de la CNN: {contar(cnn)}")
print(f"Parámetros del MLP:   {contar(mlp)}")

torch.manual_seed(0)
hist_cnn = entrenar(cnn, cargador_ent, cargador_val, epocas=8,
                    optimizador=torch.optim.Adam(cnn.parameters(),
                                                 lr=1e-3),
                    dispositivo=device)

torch.manual_seed(0)
hist_mlp = entrenar(mlp, cargador_ent, cargador_val, epocas=8,
                    optimizador=torch.optim.Adam(mlp.parameters(),
                                                 lr=1e-3),
                    dispositivo=device)

epocas = range(1, 9)
plt.plot(epocas, hist_mlp["exactitud_val"], "o-", color="C0",
         label="MLP")
plt.plot(epocas, hist_cnn["exactitud_val"], "o-", color="C1",
         label="CNN")
plt.xlabel("Época")
plt.ylabel("Exactitud de validación")
plt.title("Mismo tamaño, distinta forma")
plt.legend()
plt.show()

# cpu(): por si la red se entrenó en una GPU
filtros = cnn[0].weight.detach().cpu()        # (16, 1, 3, 3)
print(tuple(filtros.shape))

fig, ejes = plt.subplots(4, 4, figsize=(4.8, 4.8))
for k, eje in enumerate(ejes.flat):
    eje.imshow(filtros[k, 0], cmap=blanco_tinta)
    eje.set_title(f"Filtro {k}", fontsize=11)
    eje.axis("off")
plt.suptitle("Filtros de la primera capa")
plt.show()

# Una prenda de validación: la que está en la posición 1
imagen = x_val[1:2].to(device)                # (1, 1, 28, 28)
print("Etiqueta:", nombres[y_val[1]])
with torch.no_grad():
    mapas = cnn[1](cnn[0](imagen))[0].cpu()   # (16, 28, 28)
print(tuple(mapas.shape))

fig, ejes = plt.subplots(4, 4, figsize=(4.8, 4.8))
for k, eje in enumerate(ejes.flat):
    eje.imshow(mapas[k], cmap=blanco_tinta)
    eje.set_title(f"Mapa {k}", fontsize=11)
    eje.axis("off")
plt.suptitle("Mapas de activación de la primera capa")
plt.show()

from torchvision import transforms

def ruido(x):
    """Suma a cada píxel un número al azar de desviación 0,1."""
    return x + 0.1 * torch.randn_like(x)

aumentar = transforms.Compose([
    transforms.RandomHorizontalFlip(p=0.5),
    transforms.RandomCrop(28, padding=3),
    transforms.Lambda(ruido),
])

torch.manual_seed(0)
# Una imagen de entrenamiento, sin normalizar: (1, 28, 28)
original = entrenamiento.data[0].float().unsqueeze(0) / 255
fig, ejes = plt.subplots(3, 3, figsize=(4.8, 4.8))
for k, eje in enumerate(ejes.flat):
    version = original if k == 0 else aumentar(original)
    eje.imshow(version[0], cmap=blanco_tinta, vmin=0, vmax=1)
    eje.set_title("Original" if k == 0 else f"Versión {k}", fontsize=11)
    eje.axis("off")
plt.show()
print("Forma de cada versión:", tuple(aumentar(original).shape))