Skip to article frontmatterSkip to article content
Site not loading correctly?

This may be due to an incorrect BASE_URL configuration. See the MyST Documentation for reference.

Un auto-encodeur est un réseau de neurones entraîné à reproduire sa propre entrée. Cela paraît vain tant qu’on n’ajoute pas une contrainte : le réseau doit faire passer les données par un goulot d’étranglement étroit avant de les reconstruire. Pour y parvenir, il lui faut apprendre une représentation compacte qui conserve la structure des données et rejette le reste.

Auto-encoder

L’architecture comporte trois parties et elle est habituellement symétrique :

  • l’encodeur comprime l’entrée en un petit jeu de caractéristiques (couches linéaires, couches convolutives, etc.),
  • le goulot d’étranglement (ou espace latent) est la couche la plus petite, la représentation en basse dimension des données,
  • le décodeur prend les caractéristiques latentes et reconstruit les données d’origine.

Parce que la cible d’entraînement est l’entrée elle-même, aucune étiquette n’est nécessaire. C’est la forme la plus simple d’apprentissage auto-supervisé, et c’est pourquoi les auto-encodeurs comptent : après entraînement, l’encodeur est un extracteur de caractéristiques appris sur des données non étiquetées. On emploie les auto-encodeurs pour la compression, le débruitage (denoising), l’extraction de caractéristiques et la détection d’anomalies.

Dans ce carnet, nous construisons trois auto-encodeurs sur des spectrogrammes de sismogrammes synthétiques : un auto-encodeur dense, un auto-encodeur convolutif et un auto-encodeur débruiteur qui transforme des spectrogrammes bruités en spectrogrammes propres. Nous terminons par une démonstration d’auto-encodeur masqué — exécutée sur deux domaines, spectrogrammes sismiques et champs climatiques maillés — et par des expériences de transfert, dont une sur de vraies formes d’onde miniPNW, qui montrent pourquoi l’encodeur entraîné est la partie qui vaut d’être gardée, et jusqu’où ses caractéristiques portent.

Un bon panorama des variantes d’auto-encodeurs : le billet de blog de Lilian Weng.

🖥️ Diapositives du cours — Séance 26 (mer. 2 déc., contexte d’enrichissement)

import os
import numpy as np
import matplotlib.pyplot as plt
import scipy.signal
import torch
import torch.nn as nn
from torch.utils.data import DataLoader, TensorDataset
from torchinfo import summary

from mlgeo_synth import seismogram_dataset, synthetic_seismogram

device = torch.device("cuda" if torch.cuda.is_available()
                      else "mps" if torch.backends.mps.is_available()
                      else "cpu")
print("device:", device)

torch.manual_seed(0)
rng = np.random.default_rng(0)
device: mps

1. Un jeu de spectrogrammes issu de sismogrammes synthétiques

Nous générons 600 sismogrammes d’événements et 600 sismogrammes de bruit avec mlgeo_synth (30 s à 100 Hz, soit 3 000 échantillons chacun). Plutôt que d’alimenter le réseau avec des formes d’onde brutes, nous travaillons sur des log-spectrogrammes : des images temps-fréquence calculées par transformée de Fourier à court terme. Un spectrogramme transforme un signal 1-D en image 2-D, ce qui nous permet de réutiliser tout ce que nous savons des réseaux convolutifs. C’est aussi ce que font les débruiteurs de recherche comme DeepDenoiser (nous y revenons à la section 4).

Les paramètres du spectrogramme sont choisis pour que toutes les images aient la même taille : nperseg=128 donne 65 bandes de fréquence (nous supprimons la bande continue pour en obtenir 64) et noverlap=83 donne exactement 64 trames temporelles pour une trace de 3 000 échantillons. Nous prenons le log10 de la puissance avec un petit plancher pour éviter log(0), puis nous normalisons avec un unique minimum et un unique maximum globaux, de sorte que toutes les images vivent dans [0, 1] à la même échelle. Cette échelle partagée comptera plus tard, quand l’entrée et la cible du débruiteur devront être comparables.

# Generate the waveform dataset
X_wave, y_lab, metas = seismogram_dataset(n_events=600, n_noise=600, fs=100.0,
                                          duration_s=30.0, seed=0)
print("waveforms:", X_wave.shape, "| events:", int(y_lab.sum()), "| noise:", int((1 - y_lab).sum()))

FS = 100.0
NPERSEG, NOVERLAP = 128, 83
EPS = 1e-12

def log_spectrogram(trace):
    """Return a 64x64 log-power spectrogram of a 3000-sample trace."""
    f, t, Sxx = scipy.signal.spectrogram(trace, fs=FS, nperseg=NPERSEG, noverlap=NOVERLAP)
    return np.log10(Sxx[1:, :] + EPS), f[1:], t   # drop the DC bin -> 64 x 64

# Compute all spectrograms
S0, freqs, times = log_spectrogram(X_wave[0])
X_spec = np.zeros((len(X_wave),) + S0.shape, dtype=np.float32)
for i, tr in enumerate(X_wave):
    X_spec[i] = log_spectrogram(tr)[0]

# Global normalization to [0, 1] -- store the constants for reuse
S_MIN, S_MAX = X_spec.min(), X_spec.max()
def normalize(spec):
    return np.clip((spec - S_MIN) / (S_MAX - S_MIN), 0.0, 1.0).astype(np.float32)

X_img = normalize(X_spec)
print("images:", X_img.shape, "| range:", X_img.min(), "-", X_img.max())
waveforms: (1200, 3000) | events: 600 | noise: 600
images: (1200, 64, 64) | range: 0.0 - 1.0
# Show a few spectrograms with their labels
fig, axs = plt.subplots(2, 4, figsize=(10, 4.5), sharex=True, sharey=True)
idx_ev = np.where(y_lab == 1)[0][:4]
idx_no = np.where(y_lab == 0)[0][:4]
for k in range(4):
    axs[0, k].pcolormesh(times, freqs, X_img[idx_ev[k]], cmap="viridis", vmin=0, vmax=1)
    axs[0, k].set_title(f"event (snr={metas[idx_ev[k]]['snr']:.1f})", fontsize=9)
    axs[1, k].pcolormesh(times, freqs, X_img[idx_no[k]], cmap="viridis", vmin=0, vmax=1)
    axs[1, k].set_title("noise", fontsize=9)
for ax in axs[1, :]:
    ax.set_xlabel("time (s)")
for ax in axs[:, 0]:
    ax.set_ylabel("frequency (Hz)")
plt.tight_layout()
plt.show()
<Figure size 1000x450 with 8 Axes>

Les événements présentent la signature classique : un début impulsif, une énergie concentrée à la fréquence coin de la source, et une coda qui décroît avec le temps. Le bruit remplit toute l’image plus uniformément. Notez que les événements à faible rapport signal sur bruit sont difficiles à repérer à l’œil : c’est tout l’intérêt de construire un débruiteur plus loin.

1.1 Découpage entraînement/validation et chargeurs de données

Nous découpons 80/20 et enveloppons les images dans des DataLoader PyTorch. Chaque chargeur produit des paires (input, target). Pour un auto-encodeur simple, la cible est l’entrée elle-même.

n = len(X_img)
perm = rng.permutation(n)
n_train = int(0.8 * n)
itrain, ival = perm[:n_train], perm[n_train:]

# Tensors with a channel dimension: (N, 1, 64, 64)
T_train = torch.from_numpy(X_img[itrain]).unsqueeze(1)
T_val = torch.from_numpy(X_img[ival]).unsqueeze(1)
y_train = torch.from_numpy(y_lab[itrain]).long()
y_val = torch.from_numpy(y_lab[ival]).long()

train_loader = DataLoader(TensorDataset(T_train, T_train), batch_size=64, shuffle=True)
val_loader = DataLoader(TensorDataset(T_val, T_val), batch_size=64, shuffle=False)

# validation examples used for display: the highest-SNR events, easiest to see
snr_val = np.array([metas[j]["snr"] if y_lab[j] == 1 else 0.0 for j in ival])
disp = torch.from_numpy(np.argsort(snr_val)[::-1][:6].copy())
print("train:", T_train.shape, "| val:", T_val.shape)
train: torch.Size([960, 1, 64, 64]) | val: torch.Size([240, 1, 64, 64])

2. Auto-encodeur dense

Le premier auto-encodeur n’emploie que des couches entièrement connectées. L’image 64x64 est aplatie en un vecteur de dimension 4096, comprimée en un vecteur latent de 32 nombres (une compression d’un facteur 128), puis redéployée. La Sigmoid finale maintient la sortie dans [0, 1], conformément aux images normalisées.

LATENT = 32

class DenseEncoder(nn.Module):
    def __init__(self, latent=LATENT):
        super().__init__()
        self.net = nn.Sequential(
            nn.Flatten(),
            nn.Linear(64 * 64, 256), nn.SELU(),
            nn.Linear(256, latent), nn.SELU(),   # bottleneck
        )
    def forward(self, x):
        return self.net(x)

class DenseDecoder(nn.Module):
    def __init__(self, latent=LATENT):
        super().__init__()
        self.net = nn.Sequential(
            nn.Linear(latent, 256), nn.SELU(),
            nn.Linear(256, 64 * 64), nn.Sigmoid(),
        )
    def forward(self, x):
        return self.net(x).view(-1, 1, 64, 64)

class AutoEncoder(nn.Module):
    """Generic wrapper: any encoder followed by any decoder."""
    def __init__(self, encoder, decoder):
        super().__init__()
        self.encoder = encoder
        self.decoder = decoder
    def forward(self, x):
        return self.decoder(self.encoder(x))

dense_ae = AutoEncoder(DenseEncoder(), DenseDecoder()).to(device)
summary(dense_ae, input_size=(64, 1, 64, 64), device=device)
========================================================================================== Layer (type:depth-idx) Output Shape Param # ========================================================================================== AutoEncoder [64, 1, 64, 64] -- ├─DenseEncoder: 1-1 [64, 32] -- │ └─Sequential: 2-1 [64, 32] -- │ │ └─Flatten: 3-1 [64, 4096] -- │ │ └─Linear: 3-2 [64, 256] 1,048,832 │ │ └─SELU: 3-3 [64, 256] -- │ │ └─Linear: 3-4 [64, 32] 8,224 │ │ └─SELU: 3-5 [64, 32] -- ├─DenseDecoder: 1-2 [64, 1, 64, 64] -- │ └─Sequential: 2-2 [64, 4096] -- │ │ └─Linear: 3-6 [64, 256] 8,448 │ │ └─SELU: 3-7 [64, 256] -- │ │ └─Linear: 3-8 [64, 4096] 1,052,672 │ │ └─Sigmoid: 3-9 [64, 4096] -- ========================================================================================== Total params: 2,118,176 Trainable params: 2,118,176 Non-trainable params: 0 Total mult-adds (Units.MEGABYTES): 135.56 ========================================================================================== Input size (MB): 1.05 Forward/backward pass size (MB): 2.38 Params size (MB): 8.47 Estimated Total Size (MB): 11.90 ==========================================================================================

2.1 Fonction d’entraînement

Une seule fonction entraîne tous les modèles de ce carnet. Elle lit des paires (input, target) dans le chargeur : le même code traite donc la reconstruction simple (cible = entrée), le débruitage (cible = image propre) et le masquage (cible = image complète). La perte est l’erreur quadratique moyenne entre la reconstruction et la cible.

def train_ae(model, train_loader, val_loader=None, n_epochs=8, lr=1e-3, print_every=1):
    criterion = nn.MSELoss()
    optimizer = torch.optim.Adam(model.parameters(), lr=lr)
    loss_train = np.zeros(n_epochs)
    loss_val = np.zeros(n_epochs)
    for epoch in range(n_epochs):
        model.train()
        running = 0.0
        for xb, tb in train_loader:
            xb, tb = xb.to(device), tb.to(device)
            optimizer.zero_grad()
            loss = criterion(model(xb), tb)
            loss.backward()
            optimizer.step()
            running += loss.item()
        loss_train[epoch] = running / len(train_loader)
        if val_loader is not None:
            model.eval()
            running = 0.0
            with torch.no_grad():
                for xb, tb in val_loader:
                    xb, tb = xb.to(device), tb.to(device)
                    running += criterion(model(xb), tb).item()
            loss_val[epoch] = running / len(val_loader)
            if (epoch + 1) % print_every == 0:
                print(f"[epoch {epoch+1:2d}] train loss: {loss_train[epoch]:.4f}"
                      f" - val loss: {loss_val[epoch]:.4f}")
        elif (epoch + 1) % print_every == 0:
            print(f"[epoch {epoch+1:2d}] train loss: {loss_train[epoch]:.4f}")
    return loss_train, loss_val
loss_d, loss_dv = train_ae(dense_ae, train_loader, val_loader, n_epochs=8, lr=1e-3)

plt.figure(figsize=(5, 3))
plt.plot(np.arange(1, len(loss_d) + 1), loss_d, label="training loss")
plt.plot(np.arange(1, len(loss_dv) + 1), loss_dv, label="validation loss")
plt.xlabel("epoch"); plt.ylabel("MSE loss"); plt.legend(); plt.title("Dense autoencoder")
plt.tight_layout(); plt.show()
[epoch  1] train loss: 0.0293 - val loss: 0.0054
[epoch  2] train loss: 0.0049 - val loss: 0.0040
[epoch  3] train loss: 0.0039 - val loss: 0.0036
[epoch  4] train loss: 0.0035 - val loss: 0.0034
[epoch  5] train loss: 0.0032 - val loss: 0.0033
[epoch  6] train loss: 0.0032 - val loss: 0.0031
[epoch  7] train loss: 0.0032 - val loss: 0.0032
[epoch  8] train loss: 0.0048 - val loss: 0.0039
<Figure size 500x300 with 1 Axes>

Une remarque sur le calcul : 8 époques sur 960 petites images suffisent à observer le comportement. Sur votre propre machine, augmentez le nombre d’époques ou la taille du jeu de données pour des reconstructions plus nettes.

2.2 Entrée et reconstruction

def show_reconstruction(model, inputs, targets=None, n_images=5, row_titles=None):
    """Top row: inputs. Middle (optional): targets. Bottom: model reconstructions."""
    model.eval()
    with torch.no_grad():
        recon = model(inputs[:n_images].to(device)).cpu().numpy().squeeze(1)
    rows = [inputs[:n_images].numpy().squeeze(1), recon]
    if targets is not None:
        rows.append(targets[:n_images].numpy().squeeze(1))
    if row_titles is None:
        row_titles = ["input", "reconstruction", "target"][: len(rows)]
    fig, axs = plt.subplots(len(rows), n_images, figsize=(2 * n_images, 2 * len(rows)))
    for r in range(len(rows)):
        for c in range(n_images):
            axs[r, c].imshow(rows[r][c], origin="lower", cmap="viridis", vmin=0, vmax=1)
            axs[r, c].set_xticks([]); axs[r, c].set_yticks([])
        axs[r, 0].set_ylabel(row_titles[r])
    plt.tight_layout(); plt.show()

show_reconstruction(dense_ae, T_val[disp], n_images=6)
<Figure size 1200x400 with 12 Axes>

L’auto-encodeur dense reproduit la partie lisse de chaque image, le niveau de fond et la bande brillante à basse fréquence, mais les arrivées franches ont disparu. Aplatir l’image jette la structure de voisinage 2-D, et 32 nombres latents ne peuvent pas mémoriser où se situe une fine bande verticale.

3. Auto-encodeur convolutif

Les couches convolutives respectent la structure 2-D du spectrogramme. L’encodeur divise trois fois la taille de l’image par deux au moyen de convolutions à pas (stride) 2 (64 -> 32 -> 16 -> 8), et le décodeur le reflète avec des couches ConvTranspose2d. Nous ajoutons une BatchNorm2d après chaque convolution, pratique standard qui accélère et stabilise l’entraînement.

Le goulot d’étranglement est maintenant une carte de caractéristiques 64x8x8. Comptez les nombres : 4096, autant que l’entrée. La compression est ici spatiale, pas en nombre brut : chacune des 8x8 positions doit résumer une imagette 8x8 en 64 caractéristiques, si bien que le réseau ne peut toujours pas recopier les pixels d’un bout à l’autre. Des formes de goulot différentes imposent des contraintes différentes, et nous y revenons ci-dessous.

class ConvEncoder(nn.Module):
    def __init__(self):
        super().__init__()
        self.net = nn.Sequential(
            nn.Conv2d(1, 16, 3, stride=2, padding=1),   # -> 16 x 32 x 32
            nn.BatchNorm2d(16), nn.ReLU(),
            nn.Conv2d(16, 32, 3, stride=2, padding=1),  # -> 32 x 16 x 16
            nn.BatchNorm2d(32), nn.ReLU(),
            nn.Conv2d(32, 64, 3, stride=2, padding=1),  # -> 64 x 8 x 8
            nn.BatchNorm2d(64), nn.ReLU(),
        )
    def forward(self, x):
        return self.net(x)

class ConvDecoder(nn.Module):
    def __init__(self):
        super().__init__()
        self.net = nn.Sequential(
            nn.ConvTranspose2d(64, 32, 3, stride=2, padding=1, output_padding=1),
            nn.BatchNorm2d(32), nn.ReLU(),
            nn.ConvTranspose2d(32, 16, 3, stride=2, padding=1, output_padding=1),
            nn.BatchNorm2d(16), nn.ReLU(),
            nn.ConvTranspose2d(16, 1, 3, stride=2, padding=1, output_padding=1),
            nn.Sigmoid(),
        )
    def forward(self, x):
        return self.net(x)

conv_ae = AutoEncoder(ConvEncoder(), ConvDecoder()).to(device)
summary(conv_ae, input_size=(64, 1, 64, 64), device=device)
========================================================================================== Layer (type:depth-idx) Output Shape Param # ========================================================================================== AutoEncoder [64, 1, 64, 64] -- ├─ConvEncoder: 1-1 [64, 64, 8, 8] -- │ └─Sequential: 2-1 [64, 64, 8, 8] -- │ │ └─Conv2d: 3-1 [64, 16, 32, 32] 160 │ │ └─BatchNorm2d: 3-2 [64, 16, 32, 32] 32 │ │ └─ReLU: 3-3 [64, 16, 32, 32] -- │ │ └─Conv2d: 3-4 [64, 32, 16, 16] 4,640 │ │ └─BatchNorm2d: 3-5 [64, 32, 16, 16] 64 │ │ └─ReLU: 3-6 [64, 32, 16, 16] -- │ │ └─Conv2d: 3-7 [64, 64, 8, 8] 18,496 │ │ └─BatchNorm2d: 3-8 [64, 64, 8, 8] 128 │ │ └─ReLU: 3-9 [64, 64, 8, 8] -- ├─ConvDecoder: 1-2 [64, 1, 64, 64] -- │ └─Sequential: 2-2 [64, 1, 64, 64] -- │ │ └─ConvTranspose2d: 3-10 [64, 32, 16, 16] 18,464 │ │ └─BatchNorm2d: 3-11 [64, 32, 16, 16] 64 │ │ └─ReLU: 3-12 [64, 32, 16, 16] -- │ │ └─ConvTranspose2d: 3-13 [64, 16, 32, 32] 4,624 │ │ └─BatchNorm2d: 3-14 [64, 16, 32, 32] 32 │ │ └─ReLU: 3-15 [64, 16, 32, 32] -- │ │ └─ConvTranspose2d: 3-16 [64, 1, 64, 64] 145 │ │ └─Sigmoid: 3-17 [64, 1, 64, 64] -- ========================================================================================== Total params: 46,849 Trainable params: 46,849 Non-trainable params: 0 Total mult-adds (Units.MEGABYTES): 805.85 ========================================================================================== Input size (MB): 1.05 Forward/backward pass size (MB): 56.62 Params size (MB): 0.19 Estimated Total Size (MB): 57.86 ==========================================================================================
loss_c, loss_cv = train_ae(conv_ae, train_loader, val_loader, n_epochs=15, lr=1e-3)

plt.figure(figsize=(5, 3))
plt.plot(np.arange(1, len(loss_c) + 1), loss_c, label="training loss")
plt.plot(np.arange(1, len(loss_cv) + 1), loss_cv, label="validation loss")
plt.xlabel("epoch"); plt.ylabel("MSE loss"); plt.legend(); plt.title("Convolutional autoencoder")
plt.tight_layout(); plt.show()
[epoch  1] train loss: 0.0644 - val loss: 0.0683
[epoch  2] train loss: 0.0327 - val loss: 0.0327
[epoch  3] train loss: 0.0214 - val loss: 0.0190
[epoch  4] train loss: 0.0145 - val loss: 0.0115
[epoch  5] train loss: 0.0100 - val loss: 0.0078
[epoch  6] train loss: 0.0071 - val loss: 0.0059
[epoch  7] train loss: 0.0056 - val loss: 0.0049
[epoch  8] train loss: 0.0046 - val loss: 0.0042
[epoch  9] train loss: 0.0040 - val loss: 0.0037
[epoch 10] train loss: 0.0035 - val loss: 0.0030
[epoch 11] train loss: 0.0030 - val loss: 0.0027
[epoch 12] train loss: 0.0027 - val loss: 0.0026
[epoch 13] train loss: 0.0025 - val loss: 0.0025
[epoch 14] train loss: 0.0024 - val loss: 0.0024
[epoch 15] train loss: 0.0024 - val loss: 0.0024
<Figure size 500x300 with 1 Axes>
show_reconstruction(conv_ae, T_val[disp], n_images=6)
print(f"validation MSE  dense: {loss_dv[-1]:.4f}   conv: {loss_cv[-1]:.4f}")
<Figure size 1200x400 with 12 Axes>
validation MSE  dense: 0.0039   conv: 0.0024

L’auto-encodeur convolutif atteint une perte de validation plus basse que le dense, et les reconstructions gardent maintenant une trace des arrivées : une faible bande verticale à l’onde S et la tache brillante au début du signal. Cela reste flou, et c’est inhérent aux architectures à goulot d’étranglement : le détail fin est éliminé par construction.

Le goulot d’étranglement est un bouton de réglage, pas un choix figé. Le modèle dense a comprimé jusqu’à 32 nombres, une compression brutale, et a perdu les arrivées. Le modèle convolutif conserve 4 096 nombres disposés en une carte spatiale grossière, contrainte bien plus douce, ce qui explique une sortie plus nette. Rétrécissez le goulot et les reconstructions deviennent plus floues mais la représentation plus abstraite ; élargissez-le et la reconstruction s’améliore jusqu’à ce que, à l’extrême, le réseau puisse recopier l’entrée et n’apprenne plus rien d’utile. La bonne taille latente dépend de ce que vous voulez faire de la représentation, idée sur laquelle nous revenons à la section 6. Comme exercice à la maison, ajoutez une quatrième convolution à pas 2 (goulot 128x4x4) ou projetez la carte de caractéristiques sur un vecteur de 32 nombres, réentraînez et observez comment change la qualité de reconstruction.

4. Auto-encodeur débruiteur : le gain sismologique

Jusqu’ici le réseau reproduit son entrée. Un petit changement lui fait faire quelque chose de véritablement utile : donnez-lui un spectrogramme bruité en entrée et le spectrogramme propre du même événement comme cible. Le réseau ne peut plus apprendre l’identité ; il lui faut apprendre à quoi ressemble un signal de séisme et à quoi ressemble le bruit, puis ne garder que le premier.

Pour l’entraîner, il nous faut des paires appariées propre/bruité, ce à quoi sert précisément un générateur synthétique. synthetic_seismogram appelé deux fois avec les mêmes paramètres de source et la même graine produit exactement le même événement sous-jacent ; seule l’amplitude du bruit change avec snr. Nous employons snr=1e6 (un bruit un million de fois plus faible que le signal, donc pratiquement propre) pour la cible et un snr aléatoire entre 5 et 30 pour l’entrée : bruité, mais avec l’événement encore présent dans les données. La cellule ci-dessous vérifie l’appariement : le résidu entre les deux traces a l’amplitude attendue d’après le rapport signal sur bruit et n’est pas corrélé au signal.

n_pairs = 400
noisy_w = np.zeros((n_pairs, 3000))
clean_w = np.zeros((n_pairs, 3000))
prng = np.random.default_rng(42)
for i in range(n_pairs):
    mag = prng.uniform(1.5, 3.5)
    dist = prng.uniform(5, 60)
    snr_low = 10 ** prng.uniform(np.log10(5.0), np.log10(30.0))   # noisy input
    kw = dict(duration_s=30.0, fs=100.0, magnitude=mag, distance_km=dist, seed=1000 + i)
    _, clean_w[i], mc = synthetic_seismogram(snr=1e6, **kw)
    t, noisy_w[i], mn = synthetic_seismogram(snr=snr_low, **kw)

# Verify: same signal, different noise
resid = noisy_w[0] - clean_w[0]
print(f"identical source? peak amplitudes: {mc['peak_amplitude']:.4f} vs {mn['peak_amplitude']:.4f}")
print(f"residual std: {resid.std():.4f} (expected noise std: {mn['peak_amplitude']/mn['snr']:.4f})")
print(f"correlation of residual with clean signal: {np.corrcoef(resid, clean_w[0])[0, 1]:+.3f}")

fig, axs = plt.subplots(2, 1, figsize=(8, 3.5), sharex=True)
axs[0].plot(t, noisy_w[0], lw=0.5, color="gray"); axs[0].set_ylabel("noisy")
axs[1].plot(t, clean_w[0], lw=0.5, color="C0"); axs[1].set_ylabel("clean")
axs[1].set_xlabel("time (s)")
plt.tight_layout(); plt.show()
identical source? peak amplitudes: 0.4918 vs 0.4918
residual std: 0.0244 (expected noise std: 0.0264)
correlation of residual with clean signal: +0.004
<Figure size 800x350 with 2 Axes>
# Spectrograms of both, normalized with the SAME global constants as before
S_noisy = np.stack([normalize(log_spectrogram(w)[0]) for w in noisy_w])
S_clean = np.stack([normalize(log_spectrogram(w)[0]) for w in clean_w])

n_tr = int(0.8 * n_pairs)
Tn_train = torch.from_numpy(S_noisy[:n_tr]).unsqueeze(1)
Tc_train = torch.from_numpy(S_clean[:n_tr]).unsqueeze(1)
Tn_val = torch.from_numpy(S_noisy[n_tr:]).unsqueeze(1)
Tc_val = torch.from_numpy(S_clean[n_tr:]).unsqueeze(1)

den_train = DataLoader(TensorDataset(Tn_train, Tc_train), batch_size=64, shuffle=True)
den_val = DataLoader(TensorDataset(Tn_val, Tc_val), batch_size=64, shuffle=False)
print("denoiser training pairs:", len(Tn_train), "| validation pairs:", len(Tn_val))
denoiser training pairs: 320 | validation pairs: 80

L’architecture est inchangée : une instance neuve de l’auto-encodeur convolutif de la section 3. Seules les données ont changé. La tâche est plus difficile que la reconstruction simple, car l’entrée et la cible se ressemblent désormais très peu, et il existe un raccourci facile : produire partout le fond sombre moyen. Ce raccourci est un minimum local puissant (sans normalisation par lots, le réseau s’y enlise), aussi entraînons-nous plus longtemps. Le modèle est minuscule et chaque époque prend une fraction de seconde.

denoiser = AutoEncoder(ConvEncoder(), ConvDecoder()).to(device)
loss_n, loss_nv = train_ae(denoiser, den_train, den_val, n_epochs=80, lr=1e-3,
                           print_every=10)
[epoch 10] train loss: 0.0933 - val loss: 0.0891
[epoch 20] train loss: 0.0355 - val loss: 0.0354
[epoch 30] train loss: 0.0208 - val loss: 0.0208
[epoch 40] train loss: 0.0151 - val loss: 0.0157
[epoch 50] train loss: 0.0122 - val loss: 0.0129
[epoch 60] train loss: 0.0102 - val loss: 0.0117
[epoch 70] train loss: 0.0090 - val loss: 0.0109
[epoch 80] train loss: 0.0079 - val loss: 0.0103
# Noisy input / denoised output / clean target triplets
show_reconstruction(denoiser, Tn_val, targets=Tc_val, n_images=5,
                    row_titles=["noisy input", "denoised output", "clean target"])
<Figure size 1000x600 with 15 Axes>

Le débruiteur supprime le plancher de bruit à large bande et retrouve l’énergie de l’événement : l’arrivée S, la coda, la concentration à basse fréquence. La sortie n’est pas parfaite ; les arrivées sont étalées et les exemples à faible rapport signal sur bruit perdent du détail. Mais souvenez-vous de ce qu’est ce réseau : six couches convolutives entraînées 80 époques sur 320 paires d’images.

Cette miniature a des équivalents directs à l’échelle de la recherche. DeepDenoiser (Zhu et al., 2019) travaille sur la transformée de Fourier à court terme des sismogrammes, exactement notre représentation d’entrée, et prédit des masques temps-fréquence qui séparent le signal du bruit ; il est employé dans des chaînes opérationnelles de surveillance sismique. WaveDecompNet (Yin et al., 2022) pousse l’idée plus loin avec un décodeur à deux branches qui décompose un enregistrement en une composante séisme et une composante bruit, en gardant les deux, parce que le « bruit » (le champ d’ondes ambiant) est lui-même scientifiquement utile. Tous deux sont au fond des réseaux encodeur-décodeur ; ils diffèrent de notre jouet par la profondeur, par la taille de l’ensemble d’entraînement et par les connexions de saut évoquées à la section 7.

5. Auto-encodeur masqué : dix lignes de changement

Le débruitage est une façon de corrompre l’entrée ; le masquage en est une autre. Mettez à zéro des imagettes carrées choisies au hasard et demandez au réseau de reconstruire l’image entière. Pour combler une imagette manquante, le réseau ne peut pas recopier des pixels ; il lui faut comprendre le contexte, par exemple que la coda d’un événement se poursuit régulièrement dans le temps et que le bruit a une forme spectrale constante. Le changement par rapport à la section 3 tient en une dizaine de lignes : une fonction de masquage et une nouvelle paire de chargeurs.

def mask_patches(imgs, patch=8, n_patches=12, seed=0):
    """Zero out n_patches random patch x patch squares in each image."""
    g = np.random.default_rng(seed)
    masked = imgs.clone()
    N, _, H, W = imgs.shape
    for i in range(N):
        for _ in range(n_patches):
            r, c = g.integers(0, H - patch), g.integers(0, W - patch)
            masked[i, 0, r:r + patch, c:c + patch] = 0.0
    return masked

M_train = mask_patches(T_train, seed=1)
M_val = mask_patches(T_val, seed=2)
mask_train = DataLoader(TensorDataset(M_train, T_train), batch_size=64, shuffle=True)
mask_val = DataLoader(TensorDataset(M_val, T_val), batch_size=64, shuffle=False)

masked_ae = AutoEncoder(ConvEncoder(), ConvDecoder()).to(device)
loss_m, loss_mv = train_ae(masked_ae, mask_train, mask_val, n_epochs=8, lr=1e-3)
[epoch  1] train loss: 0.0614 - val loss: 0.0478
[epoch  2] train loss: 0.0298 - val loss: 0.0262
[epoch  3] train loss: 0.0167 - val loss: 0.0131
[epoch  4] train loss: 0.0097 - val loss: 0.0071
[epoch  5] train loss: 0.0062 - val loss: 0.0045
[epoch  6] train loss: 0.0044 - val loss: 0.0035
[epoch  7] train loss: 0.0035 - val loss: 0.0031
[epoch  8] train loss: 0.0031 - val loss: 0.0029
show_reconstruction(masked_ae, M_val[disp], targets=T_val[disp], n_images=5,
                    row_titles=["masked input", "reconstruction", "original"])
<Figure size 1000x600 with 15 Axes>

Avec 12 imagettes 8x8 cachées au hasard, soit au plus 18,75 % de chaque image 64x64 (moins là où les imagettes se recouvrent), le réseau remplit les trous avec un contenu temps-fréquence plausible inféré du voisinage.

Cet objectif — corrompre l’entrée, prédire ce qui a été retiré — est le cœur du pré-entraînement auto-supervisé moderne. La modélisation du langage masquée (cacher des mots, les prédire) est la façon dont on pré-entraîne les modèles de langage de la famille BERT, et les auto-encodeurs masqués pour les images (He et al., 2022) ont montré que cacher 75 % des imagettes puis les reconstruire produit des caractéristiques qui se transfèrent bien à la classification et à la détection. Aucun humain n’a jamais rien étiqueté : les données se supervisent elles-mêmes, si bien que le pré-entraînement peut consommer des archives non étiquetées arbitrairement grandes.

Ce dernier point explique pourquoi cela compte en géosciences. Les réseaux sismiques enregistrent en continu et accumulent des pétaoctets de formes d’onde non étiquetées, alors que les catalogues pointés par des analystes n’en couvrent qu’une mince fraction. Un modèle pré-entraîné par masquage ou par débruitage sur l’archive brute apprend à quoi ressemblent les signaux sismiques avant même de voir une seule étiquette. Les travaux de modèles de fondation en sismologie et en télédétection suivent exactement cette recette : pré-entraînement auto-supervisé sur l’archive, puis un petit affinage supervisé pour chaque tâche aval.

5.1 Les dix mêmes lignes sur un second domaine : les champs climatiques

Rien dans l’objectif de masquage n’est sismologique. Pour le prouver, l’architecture identique et les dix lignes identiques s’exécutent sur un second domaine : les champs maillés d’anomalies climatiques du chapitre 4.3. Nous générons deux champs avec mlgeo_synth.climate_field (grille latitude-longitude 40 x 80, 30 ans d’anomalies mensuelles) et traitons la carte de chaque mois comme une image : 720 images, normalisées dans [0, 1] avec des constantes globales, exactement comme nous avons traité les spectrogrammes. L’auto-encodeur convolutif est entièrement convolutif, si bien que les cartes 40 x 80 le traversent sans aucune modification du code — le goulot d’étranglement vaut simplement 64 x 5 x 10 au lieu de 64 x 8 x 8.

from mlgeo_synth import climate_field

C_maps = np.concatenate([
    climate_field(n_lat=40, n_lon=80, n_months=360, trend_c_per_decade=0.3, seed=sd)[0]
    for sd in (0, 1)]).astype(np.float32)          # (720, 40, 80) monthly anomaly maps
C_imgs = ((C_maps - C_maps.min()) / (C_maps.max() - C_maps.min())).astype(np.float32)

perm_c = np.random.default_rng(7).permutation(len(C_imgs))
n_ct = int(0.8 * len(C_imgs))
Ct_train = torch.from_numpy(C_imgs[perm_c[:n_ct]]).unsqueeze(1)
Ct_val = torch.from_numpy(C_imgs[perm_c[n_ct:]]).unsqueeze(1)

Cm_train, Cm_val = mask_patches(Ct_train, seed=3), mask_patches(Ct_val, seed=4)
cmask_train = DataLoader(TensorDataset(Cm_train, Ct_train), batch_size=64, shuffle=True)
cmask_val = DataLoader(TensorDataset(Cm_val, Ct_val), batch_size=64, shuffle=False)

masked_ae_clim = AutoEncoder(ConvEncoder(), ConvDecoder()).to(device)
loss_mc, loss_mcv = train_ae(masked_ae_clim, cmask_train, cmask_val, n_epochs=8, lr=1e-3)
[epoch  1] train loss: 0.0517 - val loss: 0.0481
[epoch  2] train loss: 0.0297 - val loss: 0.0468
[epoch  3] train loss: 0.0199 - val loss: 0.0456
[epoch  4] train loss: 0.0119 - val loss: 0.0320
[epoch  5] train loss: 0.0068 - val loss: 0.0098
[epoch  6] train loss: 0.0040 - val loss: 0.0045
[epoch  7] train loss: 0.0027 - val loss: 0.0023
[epoch  8] train loss: 0.0020 - val loss: 0.0021
show_reconstruction(masked_ae_clim, Cm_val[:5], targets=Ct_val[:5], n_images=5,
                    row_titles=["masked input", "reconstruction", "original"])
<Figure size 1000x600 with 15 Axes>

Le réseau remplit les carrés masqués d’une structure d’anomalie lisse et cohérente en latitude, et finit près d’une MSE de validation de 0,002 sur ce domaine : les bandes saisonnières et le motif zonal lui offrent un contexte solide pour inférer. Une architecture, un objectif, deux domaines des sciences de la Terre — c’est cette portabilité, et non une reconstruction particulière, qui plaide pour le pré-entraînement auto-supervisé sur des archives maillées, et c’est exactement la recette derrière les modèles de fondation de télédétection et de météorologie.

6. C’est l’encodeur que l’on réutilise

Après entraînement, on jette habituellement le décodeur. L’objet de valeur est l’encodeur : un extracteur de caractéristiques appris sans étiquettes. Voici le test. L’encodeur de l’auto-encodeur convolutif de la section 3 n’a été entraîné que sur la reconstruction ; il n’a jamais vu d’étiquette événement/bruit. Nous examinons d’abord son espace latent, puis nous l’utilisons pour classifier avec très peu d’étiquettes.

# PCA of the frozen encoder's latent space, colored by (held-out) labels
from sklearn.decomposition import PCA

conv_ae.eval()
with torch.no_grad():
    Z_val = conv_ae.encoder(T_val.to(device)).flatten(1).cpu().numpy()
Z2 = PCA(n_components=2).fit_transform(Z_val)

plt.figure(figsize=(5, 4))
for lab, name, col in [(0, "noise", "C0"), (1, "event", "C1")]:
    m = y_val.numpy() == lab
    plt.scatter(Z2[m, 0], Z2[m, 1], s=8, c=col, label=name, alpha=0.6)
plt.xlabel("PC 1"); plt.ylabel("PC 2"); plt.legend()
plt.title("Latent space of the conv encoder (PCA)")
plt.tight_layout(); plt.show()
<Figure size 500x400 with 1 Axes>

Les deux classes se séparent déjà, au moins en partie, alors qu’aucune étiquette n’a servi à l’entraînement. Voici maintenant l’expérience : supposons qu’un analyste n’ait étiqueté que 10 % de l’ensemble d’entraînement (96 spectrogrammes). Nous comparons deux classificateurs d’architecture identique, encodeur + petite tête linéaire, entraînés sur ces mêmes 96 exemples étiquetés :

  1. Sonde linéaire (linear probe) : l’encodeur pré-entraîné, gelé ; seule la tête s’entraîne. (Il s’agit d’une sonde, pas d’un affinage : l’affinage mettrait aussi à jour les poids de l’encodeur, alors qu’une sonde les maintient fixes et mesure ce que les seules caractéristiques pré-entraînées portent.)
  2. À partir de zéro : un encodeur et une tête initialisés au hasard, entraînés en entier.

La tête est une unique couche linéaire précédée d’une BatchNorm1d qui standardise les 4 096 caractéristiques de l’encodeur ; sans cette standardisation, une sonde linéaire s’entraîne mal sur des activations ReLU brutes. Une subtilité : un encodeur gelé qui contient des couches de normalisation par lots doit rester en mode eval() pendant l’entraînement, faute de quoi ses statistiques glissantes continuent de se mettre à jour bien que ses poids soient gelés. Nous redéfinissons train() pour l’imposer. Les deux classificateurs reçoivent la tête identique et le budget d’entraînement identique : la seule différence tient donc à la provenance des poids de l’encodeur.

import copy

# 10% of the training labels
n_lab = int(0.10 * len(T_train))
lab_idx = torch.randperm(len(T_train), generator=torch.Generator().manual_seed(3))[:n_lab]
X_lab, y_lab_small = T_train[lab_idx], y_train[lab_idx]
print(f"labeled examples: {n_lab} ({int(y_lab_small.sum())} events, {n_lab - int(y_lab_small.sum())} noise)")
lab_loader = DataLoader(TensorDataset(X_lab, y_lab_small), batch_size=32, shuffle=True)

class EncoderClassifier(nn.Module):
    def __init__(self, encoder, freeze=False):
        super().__init__()
        self.encoder = encoder
        self.frozen = freeze
        if freeze:
            for p in self.encoder.parameters():
                p.requires_grad = False
        self.head = nn.Sequential(nn.Flatten(),
                                  nn.BatchNorm1d(64 * 8 * 8),
                                  nn.Linear(64 * 8 * 8, 2))
    def train(self, mode=True):
        super().train(mode)
        if self.frozen:
            self.encoder.eval()   # frozen batch-norm stats stay fixed
        return self
    def forward(self, x):
        return self.head(self.encoder(x))

def train_classifier(model, loader, n_epochs=30, lr=1e-3, X_eval=None, y_eval=None):
    if X_eval is None:
        X_eval, y_eval = T_val, y_val
    criterion = nn.CrossEntropyLoss()
    optimizer = torch.optim.Adam([p for p in model.parameters() if p.requires_grad], lr=lr)
    acc_hist = np.zeros(n_epochs)
    for epoch in range(n_epochs):
        model.train()
        for xb, yb in loader:
            xb, yb = xb.to(device), yb.to(device)
            optimizer.zero_grad()
            loss = criterion(model(xb), yb)
            loss.backward()
            optimizer.step()
        model.eval()
        with torch.no_grad():
            pred = model(X_eval.to(device)).argmax(1).cpu()
        acc_hist[epoch] = (pred == y_eval).float().mean().item()
    return acc_hist

# 1) linear probe: frozen pretrained encoder + trainable head
torch.manual_seed(5)
probe = EncoderClassifier(copy.deepcopy(conv_ae.encoder), freeze=True).to(device)
acc_probe = train_classifier(probe, lab_loader)

# 2) same architecture from scratch on the same 96 examples
torch.manual_seed(5)
scratch = EncoderClassifier(ConvEncoder(), freeze=False).to(device)
acc_sc = train_classifier(scratch, lab_loader)

print(f"final validation accuracy  linear probe: {acc_probe[-1]:.3f}   from scratch: {acc_sc[-1]:.3f}")
labeled examples: 96 (45 events, 51 noise)
final validation accuracy  linear probe: 0.725   from scratch: 0.721
plt.figure(figsize=(5.5, 3.5))
ep = np.arange(1, len(acc_probe) + 1)
plt.plot(ep, acc_probe, "o-", label="pretrained encoder (frozen) + head")
plt.plot(ep, acc_sc, "s-", label="same architecture from scratch")
plt.xlabel("epoch"); plt.ylabel("validation accuracy")
plt.title("Event vs noise with 10% of the labels")
plt.legend(loc="lower right"); plt.ylim(0.4, 1.02)
plt.tight_layout(); plt.show()
<Figure size 550x350 with 1 Axes>

Lisez les courbes en partant de la gauche. L’encodeur pré-entraîné est utile dès la première époque : ses caractéristiques organisent déjà les données, si bien que la tête atteint son plateau presque immédiatement. Le modèle entraîné à partir de zéro passe l’essentiel de son budget au niveau du hasard pendant qu’il apprend des filtres convolutifs à partir de 96 exemples seulement, puis il grimpe et, sur cette petite tâche, finit par s’en approcher. C’est cette trajectoire qui compte : avec le pré-entraînement, vous payez l’apprentissage des caractéristiques une seule fois, sur des données non étiquetées, au lieu de le repayer à chaque tâche pauvre en étiquettes. Les chiffres exacts varient d’une exécution à l’autre ; l’écart des premières époques se creuse à mesure que les étiquettes se raréfient ou que les modèles grossissent. Essayez 5 % ou 2 % des étiquettes sur votre propre machine.

C’est là l’argument pratique de l’apprentissage auto-supervisé. Les objectifs de reconstruction, de débruitage et de masquage extraient de la structure de données non étiquetées, et les données non étiquetées sont ce dont les géosciences disposent en abondance. La ressource coûteuse, les étiquettes d’expert, n’est alors dépensée que sur un petit ensemble étiqueté, soit par une sonde linéaire comme ici, soit par un affinage complet qui met aussi à jour l’encodeur. Quand vous lisez que des modèles de fondation ont été pré-entraînés sur des archives sismiques continues ou sur des piles d’images satellitaires, c’est de cette recette en deux temps — pré-entraînement auto-supervisé puis adaptation supervisée légère — qu’il s’agit.

6.1 L’encodeur survit-il aux données réelles ? Une sonde de transfert miniPNW

La sonde ci-dessus va du synthétique au synthétique : un encodeur pré-entraîné sur des spectrogrammes synthétiques, sondé avec des étiquettes synthétiques, noté sur des données de validation synthétiques. Le chapitre 4.3 a montré ce qu’il advient quand un classificateur entraîné sur du synthétique rencontre de vraies formes d’onde : il s’effondre au niveau du hasard. La question auto-supervisée est plus subtile et plus encourageante : même si le classificateur meurt, les caractéristiques pré-entraînées traversent-elles l’écart ?

Le test : construire de vrais spectrogrammes à partir des formes d’onde miniPNW étiquetées mises en cache par le chapitre 2.11 (Ni et al., 2023) — 300 fenêtres de vrais séismes (pointé P à 7 s dans une fenêtre de 30 s) et 300 fenêtres de bruit pré-événement, chacune normalisée à un pic unitaire comme les formes d’onde synthétiques, puis passées dans la chaîne de spectrogrammes identique avec la normalisation globale identique. Puis entraîner la sonde linéaire identique — l’encodeur gelé pré-entraîné sur du synthétique plus une tête neuve — sur 96 exemples étiquetés réels, le même budget de 10 % d’étiquettes que précédemment, et évaluer sur un ensemble de test réservé de fenêtres réelles. L’exactitude synthétique-vers-synthétique figure à côté pour comparaison directe.

import urllib.request
import h5py
import pandas as pd

# Same files and loader pattern as notebook 2.11, which caches them in its data/ folder
pnw_dir = os.path.join('..', 'Chapter2-DataManipulation', 'data')
os.makedirs(pnw_dir, exist_ok=True)
metadata_path = os.path.join(pnw_dir, 'miniPNW_metadata.csv')
waveform_path = os.path.join(pnw_dir, 'miniPNW_waveforms.hdf5')
base_url = 'https://dasway.ess.washington.edu/shared/niyiyu/PNW-ML'

have_pnw = True
try:
    if not os.path.exists(metadata_path):
        urllib.request.urlretrieve(f'{base_url}/miniPNW_metadata.csv', metadata_path)
    if not os.path.exists(waveform_path):
        print('Downloading miniPNW waveforms (about 670 MB, one-time; shared with notebook 2.11)...')
        urllib.request.urlretrieve(f'{base_url}/miniPNW_waveforms.hdf5', waveform_path)
except Exception as err:
    have_pnw = False
    print('miniPNW cache is missing and the download failed, so the transfer probe below is skipped.\n'
          'Run notebook 2.11 first (it downloads and caches the files), then rerun this section.\n'
          f'Reason: {err!r}')
have_pnw = have_pnw and os.path.exists(waveform_path)
print('miniPNW available:', have_pnw)
miniPNW available: True
if have_pnw:
    meta = pd.read_csv(metadata_path)
    eq = meta[(meta['source_type'] == 'earthquake') & meta['trace_P_arrival_sample'].notna()].copy()
    n_win, pre = 3000, 700  # 30 s at 100 Hz, P pick 7 s into the window
    p_samp = eq['trace_P_arrival_sample'].astype(int)
    eq = eq[(p_samp - pre >= n_win) & (p_samp - pre + n_win <= 15001)].head(300)

    def read_z(f, trace_name):
        """Vertical component of one miniPNW trace (same reader as notebook 2.11)."""
        bucket, narray = trace_name.split('$')
        x, _, z = (int(v) for v in narray.split(',:'))
        return f['/data/' + bucket][x, 2, :z]  # channel order N, E, Z

    R_ev = np.zeros((len(eq), n_win))
    R_no = np.zeros((len(eq), n_win))
    with h5py.File(waveform_path, 'r') as f:
        for i, (_, row) in enumerate(eq.iterrows()):
            tr = read_z(f, row['trace_name']).astype(np.float64)
            pk = int(row['trace_P_arrival_sample'])
            R_ev[i] = tr[pk - pre : pk - pre + n_win]
            R_no[i] = tr[:n_win]  # the trace starts 50 s before the pick: pre-event noise

    R_wave = np.concatenate([R_ev, R_no])
    yR = np.concatenate([np.ones(len(R_ev)), np.zeros(len(R_no))]).astype(int)
    alive = R_wave.std(axis=1) > 0
    R_wave, yR = R_wave[alive], yR[alive]
    R_wave = R_wave - R_wave.mean(axis=1, keepdims=True)
    R_wave = R_wave / np.abs(R_wave).max(axis=1, keepdims=True)  # unit peak, like the synthetics

    # identical spectrogram pipeline, identical global normalization constants
    R_spec = np.stack([normalize(log_spectrogram(w)[0]) for w in R_wave])
    perm_r = np.random.default_rng(11).permutation(len(R_spec))
    n_rt = int(0.8 * len(R_spec))
    TR_train = torch.from_numpy(R_spec[perm_r[:n_rt]]).unsqueeze(1)
    TR_val = torch.from_numpy(R_spec[perm_r[n_rt:]]).unsqueeze(1)
    yR_train = torch.from_numpy(yR[perm_r[:n_rt]]).long()
    yR_val = torch.from_numpy(yR[perm_r[n_rt:]]).long()
    print(f"real spectrograms: {len(R_spec)} ({int(yR.sum())} earthquake, {int((1 - yR).sum())} noise)"
          f" | train {len(TR_train)}, val {len(TR_val)}")
real spectrograms: 600 (300 earthquake, 300 noise) | train 480, val 120
if have_pnw:
    # same 10% label budget as the synthetic probe: 96 real labeled examples
    lab_r = torch.randperm(len(TR_train), generator=torch.Generator().manual_seed(3))[:n_lab]
    real_lab_loader = DataLoader(TensorDataset(TR_train[lab_r], yR_train[lab_r]),
                                 batch_size=32, shuffle=True)

    torch.manual_seed(5)
    probe_real = EncoderClassifier(copy.deepcopy(conv_ae.encoder), freeze=True).to(device)
    acc_real = train_classifier(probe_real, real_lab_loader, X_eval=TR_val, y_eval=yR_val)

    # zero-shot for reference: the probe trained on synthetic labels, applied to real data
    probe.eval()
    with torch.no_grad():
        zeroshot = (probe(TR_val.to(device)).argmax(1).cpu() == yR_val).float().mean().item()

    print(f"linear probe on synthetic val (synthetic labels):  {acc_probe[-1]:.3f}")
    print(f"linear probe on real miniPNW val (real labels):    {acc_real[-1]:.3f}")
    print(f"synthetic-to-real transfer gap:                    {acc_probe[-1] - acc_real[-1]:.3f}")
    print(f"zero-shot (synthetic-label probe on real val):     {zeroshot:.3f}")

    plt.figure(figsize=(5.5, 3.5))
    ep = np.arange(1, len(acc_probe) + 1)
    plt.plot(ep, acc_probe, "o-", label="probe on synthetic data")
    plt.plot(ep, acc_real, "s-", label="same encoder, probe on real miniPNW")
    plt.xlabel("epoch"); plt.ylabel("validation accuracy")
    plt.title("Linear probes on the synthetic-pretrained encoder")
    plt.legend(loc="lower right"); plt.ylim(0.4, 1.02)
    plt.grid(alpha=0.3)
    plt.tight_layout(); plt.show()
linear probe on synthetic val (synthetic labels):  0.725
linear probe on real miniPNW val (real labels):    0.908
synthetic-to-real transfer gap:                    -0.183
zero-shot (synthetic-label probe on real val):     0.633
<Figure size 550x350 with 1 Axes>

Trois nombres à lire ensemble, et un signe qui pourrait surprendre. La sonde synthétique-vers-synthétique se situe à 72,5 %. Le même encodeur gelé, sondé avec 96 étiquettes réelles, atteint 90,8 % sur de vraies fenêtres miniPNW — l’« écart » de transfert est négatif, -18 points. Cela ne veut pas dire que le pré-entraînement synthétique bat les données réelles ; cela veut dire que les deux tâches ne sont pas également difficiles. L’ensemble de validation synthétique inclut délibérément des événements descendant jusqu’à un rapport signal sur bruit de 0,5, dont beaucoup sont irréductiblement indétectables, alors que les séismes miniPNW sont pointés par des analystes et le plus souvent nets, et que le vrai bruit pré-événement a une signature spectrale distinctive. Les conclusions honnêtes sont les conclusions relatives. Premièrement, le zero-shot échoue : la tête entraînée sur des étiquettes synthétiques ne fait que 63,3 % sur des fenêtres réelles, en écho à l’effondrement du classificateur entraîné sur du synthétique en 4.3 — les frontières de décision ne se transfèrent pas. Deuxièmement, les caractéristiques sous-jacentes, elles, se transfèrent : sans qu’un seul poids d’encodeur ait été mis à jour sur des données réelles, 96 étiquettes réelles suffisent à atteindre 91 %. Cette asymétrie est l’argument pratique du pré-entraînement auto-supervisé — les représentations survivent à un changement de domaine auquel les classificateurs ne survivent pas, et le prix de la traversée est d’une centaine d’étiquettes, pas d’un réseau réentraîné.

7. Au-delà du goulot d’étranglement : connexions de saut et U-Nets

Nos reconstructions sont floues parce que tout doit passer par un goulot d’étranglement en basse dimension, qui élimine le détail fin par construction. Le correctif standard est la connexion de saut (skip connection) : alimenter directement chaque étage du décodeur avec la sortie de l’étage correspondant de l’encodeur, de sorte que le détail à haute résolution contourne le goulot pendant que le chemin profond porte le contexte. Un encodeur-décodeur muni de connexions de saut à tous les niveaux est un U-Net (Ronneberger et al., 2015), le cheval de trait de la segmentation d’images et de la plupart des modèles d’apprentissage profond en sismologie.

Unet

DeepDenoiser et WaveDecompNet, vus à la section 4, emploient tous deux des connexions de saut, et les pointeurs de phases comme PhaseNet sont des U-Nets appliqués aux formes d’onde. Les encodeurs-décodeurs multitâches poussent l’idée plus loin : l’Earthquake Transformer (Mousavi et al., 2020) décode les probabilités de détection, de pointé P et de pointé S à partir d’un encodeur unique. Nous n’implémentons pas de U-Net ici ; ce qu’il faut retenir, c’est qu’un U-Net est un auto-encodeur plus des raccourcis.

Résumé

  • Un auto-encodeur apprend une représentation latente compacte en reconstruisant sa propre entrée ; aucune étiquette n’est nécessaire.
  • Les auto-encodeurs convolutifs battent les denses sur des données de type image comme les spectrogrammes.
  • Changer la cible transforme la reconstruction en quelque chose d’utile : des spectrogrammes propres donnent un débruiteur (la miniature de DeepDenoiser et de WaveDecompNet) ; reconstruire des imagettes masquées donne l’objectif qui sous-tend le pré-entraînement masqué moderne.
  • L’encodeur entraîné est le produit réutilisable : gelé derrière une sonde linéaire entraînée sur 10 % des étiquettes, il a atteint son plateau d’exactitude en quelques époques, tandis que la même architecture entraînée à partir de zéro a passé l’essentiel de son budget à réapprendre des caractéristiques.
  • L’objectif de masquage s’est porté sur un second domaine — des cartes mensuelles d’anomalies climatiques — sans autre changement qu’un nouveau chargeur de données ; c’est la portabilité entre domaines qui fait tout l’intérêt de l’auto-supervision.
  • Sondé sur de vrais spectrogrammes miniPNW, l’encodeur gelé pré-entraîné sur du synthétique soutient une exactitude de 90,8 % à partir de 96 étiquettes réelles, alors que la tête à étiquettes synthétiques ne transfère qu’à 63,3 % : les caractéristiques traversent l’écart synthétique-réel, les frontières de décision non.

Ensuite : les réseaux de neurones informés par la physique, où la fonction de perte encode elle-même les équations qui gouvernent le problème.

References
  1. Ni, Y., Hutko, A., Skene, F., Denolle, M., Malone, S., Bodin, P., Hartog, R., & Wright, A. (2023). Curated Pacific Northwest AI-ready Seismic Dataset. Seismica, 2(1). 10.26443/seismica.v2i1.368