En este cuaderno clasificamos registros sísmicos en cuatro tipos de fuente — sismo, explosión, evento superficial y ruido — a partir de 61 características físicas de la forma de onda (forma espectral, estadísticas de la envolvente, curtosis, energías por banda). El conjunto de datos es una colección curada de eventos sísmicos del noroeste del Pacífico de Estados Unidos, 1000 por clase, archivada en Zenodo: DOI 10.5281/zenodo.14025693.
Comparamos tres clasificadores clásicos: la máquina de vectores de soporte (SVM), los k vecinos más cercanos (KNN) y el bosque aleatorio, y los evaluamos por clase con reportes de clasificación, matrices de confusión y curvas ROC uno-contra-el-resto (one-vs-rest).
Este conjunto de datos regresa en el cuaderno 3.9, y ancla la tabla de clasificación del curso (leaderboard) definida al final de este cuaderno.
1. Cargar los datos¶
El cargador siguiente descarga y guarda en caché los cuatro archivos de clase, los concatena en una sola tabla y elimina la única columna de características con valores faltantes.
import numpy as np
import pandas as pd
import matplotlib.pyplot as pltimport pooch
SEISMIC_FILES = {
"1000_earthquakes_physical_features.csv": "md5:28129c8dd1b3e14f655d489577b841b5",
"1000_explosion_physical_features.csv": "md5:af1342d32e163e961e043364136359b0",
"1000_noise_physical_features.csv": "md5:16cdb992fed6cf6273d5624f5df905da",
"1000_surface_physical_features.csv": "md5:9a2c2643030cf058704d68e130654e9d",
}
frames = []
for fname, checksum in SEISMIC_FILES.items():
path = pooch.retrieve(
url=f"https://zenodo.org/api/records/14025693/files/{fname}/content",
known_hash=checksum,
fname=fname,
path=pooch.os_cache("mlgeo"),
)
frames.append(pd.read_csv(path, index_col=0))
seismic = pd.concat(frames, ignore_index=True)
seismic = seismic.dropna(axis=1) # drops the one feature column with missing values2. Explorar¶
Verifique primero el balance de clases. Las cuatro clases tienen 1000 eventos cada una, así que este conjunto curado está balanceado. Tenga presente que los catálogos reales no lo están.
seismic["source"].value_counts().plot(kind="bar")
plt.ylabel("Number of events")
plt.title("Class balance")
plt.tight_layout()
Construya la matriz de características X (61 características numéricas de la forma de onda) y el vector de etiquetas y (cuatro cadenas de tipo de fuente).
X = seismic.drop(columns=["source", "serial_no"])
y = seismic["source"]
print("X:", X.shape, " y:", y.shape)
print("A few feature names:", list(X.columns[:8]))X: (4000, 61) y: (4000,)
A few feature names: ['Window_Length', 'RappMaxMean', 'RappMaxMedian', 'AsDec', 'KurtoSig', 'KurtoEnv', 'SkewSig', 'SkewEnv']
3. División entrenamiento/prueba¶
La celda siguiente es la división canónica de este capítulo: esta línea exacta define la división de la tabla de clasificación al final del cuaderno. Todos entrenan sobre el mismo X_train y predicen sobre el mismo X_test.
Dos detalles importan:
- Una versión anterior de esta lección usaba
shuffle=False. Con la tabla ordenada por clase, eso colocaba clases enteras en el conjunto de prueba y ninguna en el de entrenamiento — el clasificador nunca veía las clases sobre las que era calificado. stratify=yconserva las proporciones de clase en ambas mitades. Eso importa sobre todo cuando las clases están desbalanceadas, y los catálogos sísmicos reales lo están: las ventanas de ruido superan en número a los sismos por órdenes de magnitud.
from sklearn.model_selection import train_test_split
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.25, random_state=2026, stratify=y)4. Escalar después de dividir¶
Ajustamos el escalador solo sobre el conjunto de entrenamiento y luego transformamos ambos conjuntos. Ajustar el escalador sobre todos los datos dejaría que las estadísticas del conjunto de prueba se fugaran hacia el entrenamiento.
from sklearn.preprocessing import StandardScaler
scaler = StandardScaler()
X_train_scaled = scaler.fit_transform(X_train)
X_test_scaled = scaler.transform(X_test)5. Modelo de referencia¶
Antes de cualquier modelo, califique un modelo de referencia (baseline) trivial. DummyClassifier(strategy="most_frequent") siempre predice la clase mayoritaria; con cuatro clases balanceadas obtiene 0.25. Todo modelo de aquí en adelante debe superar este número.
from sklearn.dummy import DummyClassifier
from sklearn import metrics
dummy = DummyClassifier(strategy="most_frequent")
dummy.fit(X_train, y_train)
baseline_acc = metrics.accuracy_score(y_test, dummy.predict(X_test))
print("Baseline accuracy:", baseline_acc)Baseline accuracy: 0.25
6. Tres clasificadores¶
La SVM y los k vecinos más cercanos dependen de distancias o de márgenes en el espacio de características, así que usan las características escaladas.
from sklearn.svm import SVC
from sklearn.neighbors import KNeighborsClassifier
from sklearn.ensemble import RandomForestClassifier
# Support Vector Machine classifier
clf = SVC(gamma='scale') # model design
clf.fit(X_train_scaled, y_train) # learn
svc_prediction = clf.predict(X_test_scaled) # predict on test
print("SVC test accuracy:", metrics.accuracy_score(y_true=y_test, y_pred=svc_prediction))
# K-nearest Neighbors
knn_clf = KNeighborsClassifier() # model design
knn_clf.fit(X_train_scaled, y_train) # learn
knn_prediction = knn_clf.predict(X_test_scaled) # predict on test
print("K-nearest Neighbors test accuracy:", metrics.accuracy_score(y_true=y_test, y_pred=knn_prediction))SVC test accuracy: 0.877
K-nearest Neighbors test accuracy: 0.83
El bosque aleatorio trabaja sobre las características sin escalar: los árboles dividen por umbrales, y un reescalado monótono no cambia de qué lado de un umbral cae un punto.
# Random Forest, on the unscaled features
rf_clf = RandomForestClassifier(random_state=42) # model design
rf_clf.fit(X_train, y_train) # learn
rf_prediction = rf_clf.predict(X_test) # predict on test
print("Random Forest test accuracy:", metrics.accuracy_score(y_true=y_test, y_pred=rf_prediction))Random Forest test accuracy: 0.884
7. Evaluación por clase¶
La exactitud es un solo número. El reporte de clasificación da la precisión, la exhaustividad (recall) y el F1 por clase, y la matriz de confusión muestra qué clases se confunden entre sí.
from sklearn.metrics import ConfusionMatrixDisplay
print("Support Vector Machine")
print(f"Classification report for classifier {clf}:\n"
f"{metrics.classification_report(y_test, svc_prediction)}\n")
disp = ConfusionMatrixDisplay.from_estimator(clf, X_test_scaled, y_test, xticks_rotation=45)
disp.figure_.suptitle("Confusion Matrix: SVC")
plt.tight_layout()
plt.show()Support Vector Machine
Classification report for classifier SVC():
precision recall f1-score support
earthquake 0.82 0.87 0.84 250
explosion 0.88 0.79 0.83 250
noise 0.90 0.93 0.92 250
surface event 0.92 0.92 0.92 250
accuracy 0.88 1000
macro avg 0.88 0.88 0.88 1000
weighted avg 0.88 0.88 0.88 1000
print("K-nearest neighbors")
print(f"Classification report for classifier {knn_clf}:\n"
f"{metrics.classification_report(y_test, knn_prediction)}\n")
disp = ConfusionMatrixDisplay.from_estimator(knn_clf, X_test_scaled, y_test, xticks_rotation=45)
disp.figure_.suptitle("Confusion Matrix: KNN")
plt.tight_layout()
plt.show()K-nearest neighbors
Classification report for classifier KNeighborsClassifier():
precision recall f1-score support
earthquake 0.78 0.84 0.81 250
explosion 0.79 0.68 0.73 250
noise 0.92 0.91 0.92 250
surface event 0.83 0.89 0.86 250
accuracy 0.83 1000
macro avg 0.83 0.83 0.83 1000
weighted avg 0.83 0.83 0.83 1000
print("Random Forest")
print(f"Classification report for classifier {rf_clf}:\n"
f"{metrics.classification_report(y_test, rf_prediction)}\n")
disp = ConfusionMatrixDisplay.from_estimator(rf_clf, X_test, y_test, xticks_rotation=45)
disp.figure_.suptitle("Confusion Matrix: Random Forest")
plt.tight_layout()
plt.show()Random Forest
Classification report for classifier RandomForestClassifier(random_state=42):
precision recall f1-score support
earthquake 0.83 0.87 0.85 250
explosion 0.88 0.80 0.83 250
noise 0.93 0.94 0.93 250
surface event 0.90 0.94 0.92 250
accuracy 0.88 1000
macro avg 0.88 0.88 0.88 1000
weighted avg 0.88 0.88 0.88 1000
¿Qué pares de clases se confunden más? Compare los términos fuera de la diagonal en las tres matrices de confusión. Una tabulación cruzada de las etiquetas verdaderas contra las predicciones del bosque aleatorio hace que los conteos sean fáciles de leer.
pd.crosstab(y_test, rf_prediction, rownames=["true"], colnames=["predicted"])La confusión dominante es entre explosiones y sismos. Eso es físicamente plausible: las voladuras de cantera y los sismos someros excitan un contenido de frecuencias similar a distancias regionales, así que sus características de forma de onda se traslapan.
8. Curvas ROC uno-contra-el-resto¶
Las curvas ROC están definidas para problemas binarios. Para un problema multiclase usamos la estrategia uno-contra-el-resto (one-vs-rest): binarizar las etiquetas (una columna por clase) y ajustar un clasificador binario por clase.
Para mantener la división idéntica a la canónica, binarizamos los y_train y y_test obtenidos arriba en lugar de volver a dividir los datos.
from sklearn.multiclass import OneVsRestClassifier
from sklearn.preprocessing import label_binarize
from sklearn import svm
from sklearn.metrics import roc_curve, auc
classes = ['earthquake', 'explosion', 'noise', 'surface event']
y_train_bin = label_binarize(y_train, classes=classes)
y_test_bin = label_binarize(y_test, classes=classes)
ovr_classifier = OneVsRestClassifier(svm.SVC(kernel='linear'))
y_score = ovr_classifier.fit(X_train_scaled, y_train_bin).decision_function(X_test_scaled)
plt.figure(figsize=(7, 6))
plt.plot([0, 1], [0, 1], 'k--', label='chance')
for i, name in enumerate(classes):
fpr, tpr, _ = roc_curve(y_test_bin[:, i], y_score[:, i])
plt.plot(fpr, tpr, label=f'{name} (AUC = {auc(fpr, tpr):.2f})')
plt.xlim([0.0, 1.0])
plt.ylim([0.0, 1.05])
plt.grid(True)
plt.xlabel('False Positive Rate')
plt.ylabel('True Positive Rate')
plt.title('One-vs-rest ROC curves, linear SVC')
plt.legend(loc="lower right")
Tabla de clasificación del curso¶
Entrene el clasificador que quiera sobre X_train. La ingeniería de características está permitida. No espíe y_test mientras desarrolla: seleccione su modelo y sus hiperparámetros con validación cruzada solo sobre el conjunto de entrenamiento — calificando los candidatos sobre subconjuntos de validación rotativos tallados dentro de los datos de entrenamiento (lección 3.8).
Cuando termine:
- Prediga sobre el
X_testcanónico (lo define la celda de división de la sección 3). - Guarde sus predicciones en
results/predictions_<uwnetid>.csvcon el formato de abajo. - Envíe el archivo mediante un pull request al repositorio del curso. La integración continua (CI) califica el F1 macro — los scores F1 por clase promediados con igual peso para cada clase — contra las etiquetas retenidas y publica una tabla de clasificación.
La columna row_id es la posición de la fila en la tabla concatenada canónica; la división la conserva, así que identifica cada muestra de prueba. La demostración de abajo escribe un archivo de envío a partir de las predicciones del bosque aleatorio.
import os
os.makedirs("results", exist_ok=True)
pred_df = pd.DataFrame({"row_id": X_test.index, "prediction": rf_prediction})
pred_df.to_csv("results/predictions_example.csv", index=False)
pred_df.head()- Kharita, A. (2024). Physical features for small sample of data (1000 events per class). Zenodo. 10.5281/ZENODO.14025693