La clasificación binaria asigna cada muestra a una de dos clases. Nuestra tarea en este capítulo es estándar en sismología: dada una ventana corta de datos de sismograma, decidir si contiene un evento sísmico (etiqueta 1) o solo ruido (etiqueta 0). En lugar de formas de onda crudas, trabajamos con cuatro características físicas calculadas a partir de cada ventana:
sta_lta: la razón entre una amplitud promedio de corto plazo y una amplitud promedio de largo plazo. Una llegada impulsiva eleva el promedio de corto plazo y empuja la razón muy por encima de 1. La distribución es de cola pesada, y abarca aproximadamente de 0.4 a 40.kurtosis: qué tan impulsiva es la distribución de amplitudes; las llegadas en forma de pico la elevan.spectral_centroid_hz: la frecuencia media de la ventana ponderada por amplitud, aproximadamente de 0.2 a 24 Hz. Los eventos locales llevan más energía de alta frecuencia que el ruido microsísmico generado por el océano.dominant_freq_hz: la frecuencia del pico espectral.
Para los gráficos bidimensionales de este cuaderno usamos sta_lta y spectral_centroid_hz: una ventana con evento tiende a tener a la vez una razón STA/LTA alta y un centroide espectral alto, mientras que el ruido se ubica en valores bajos de ambos.
1. Datos de evento contra ruido¶
Cargamos un conjunto balanceado de 2000 ventanas, mitad eventos y mitad ruido, del paquete de datos sintéticos del curso.
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import mlgeo_synth
df = mlgeo_synth.detector_features(n=2000, event_fraction=0.5, seed=42)
df.head()fig, ax = plt.subplots()
for label, name, color in [(0, "noise", "tab:gray"), (1, "event", "tab:red")]:
sub = df[df["label"] == label]
ax.scatter(sub["sta_lta"], sub["spectral_centroid_hz"],
s=14, alpha=0.5, color=color, edgecolors="k", linewidths=0.3, label=name)
plt.xscale('log')
ax.set_xlabel("STA/LTA ratio")
ax.set_ylabel("Spectral centroid (Hz)")
ax.legend()
ax.grid(True)
plt.show()
El eje x en escala logarítmica ya nos dice cómo construir la matriz de características: tomar el logaritmo de una razón de cola pesada como STA/LTA es ingeniería de características estándar, y lo hacemos antes de entregar las características a cualquier modelo.
X = np.column_stack([np.log10(df["sta_lta"]), df["spectral_centroid_hz"]])
y = df["label"].valuesAntes de cualquier clasificador, el modelo de referencia (baseline). Primero dividimos los datos: todo lo que se ajusta — un escalador, un modelo de referencia trivial, un clasificador — ve solo el conjunto de entrenamiento, y el conjunto de prueba queda intacto hasta la puntuación. El predictor trivial asigna cada ventana a la clase más frecuente; en este conjunto balanceado eso da una exactitud cercana a 0.5, y todo clasificador de aquí en adelante debe superarla.
from sklearn.dummy import DummyClassifier
from sklearn.model_selection import train_test_split
# split first: the baseline, like any model, is fit on the training split only
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.4, random_state=42)
baseline = DummyClassifier(strategy="most_frequent").fit(X_train, y_train)
print(f"Majority-class baseline accuracy: {baseline.score(X_test, y_test):.2f}")Majority-class baseline accuracy: 0.49
Comenzaremos con el fundamental LDA (análisis discriminante lineal).
from sklearn.discriminant_analysis import LinearDiscriminantAnalysis
from sklearn.preprocessing import StandardScaler
# define ML
clf = LinearDiscriminantAnalysis()
# normalize data: fit the scaler on the training split only, then apply it to
# both splits. Fitting on the full dataset would leak the test set's statistics
# into training — the exact flaw the audit exercise in 3.10 teaches you to catch.
scaler = StandardScaler().fit(X_train)
X_train = scaler.transform(X_train)
X_test = scaler.transform(X_test)
X_scaled = scaler.transform(X) # transform-only, for plotting the full set
# Fit the model.
clf.fit(X_train, y_train)
# calculate the mean accuracy on the given test data and labels.
score = clf.score(X_test, y_test)
print("The mean accuracy on the given test and labels is %f" % score)The mean accuracy on the given test and labels is 0.952500
from sklearn.inspection import DecisionBoundaryDisplay
ax = plt.subplot()
# plot the decision boundary as a background
DecisionBoundaryDisplay.from_estimator(clf, X_scaled, cmap='PiYG', alpha=0.8, ax=ax, eps=0.5)
ax.scatter(X_scaled[:, 0], X_scaled[:, 1], c=y, cmap='PiYG', alpha=0.6, edgecolors="k")
ax.set_xlabel("log10(STA/LTA), scaled")
ax.set_ylabel("Spectral centroid, scaled")
plt.show()
El LDA traza una sola línea recta a través del espacio de características. Aquí funciona bien porque las dos clases son aproximadamente separables por una línea.
Probemos un clasificador distinto: KNN (k vecinos más cercanos).
from sklearn.neighbors import KNeighborsClassifier
# define ML
K = 5
clf = KNeighborsClassifier(K)
# same split, same train-fit scaling as above
# Fit the model.
clf.fit(X_train, y_train)
# calculate the mean accuracy on the given test data and labels.
score = clf.score(X_test, y_test)
print("The mean accuracy on the given test and labels is %f" % score)The mean accuracy on the given test and labels is 0.960000
# plot the decision boundary as a background
ax = plt.subplot()
DecisionBoundaryDisplay.from_estimator(clf, X_scaled, cmap='PiYG', alpha=0.8, ax=ax, eps=0.5)
ax.scatter(X_scaled[:, 0], X_scaled[:, 1], c=y, cmap='PiYG', alpha=0.6, edgecolors="k")
ax.set_xlabel("log10(STA/LTA), scaled")
ax.set_ylabel("Spectral centroid, scaled")
plt.show()
Características con unidades físicas distintas¶
Ahora probamos qué ocurre cuando no normalizamos las características. Nuestras dos características tienen significados físicos distintos: el log-STA/LTA es adimensional, el centroide espectral es una frecuencia. Suponga que un colega nos entrega el centroide en milihercios en lugar de hercios. Nada físico ha cambiado, pero un eje es ahora numéricamente mil veces más grande.
# rebuild the raw features, spectral centroid now in millihertz
X_raw = np.column_stack([np.log10(df["sta_lta"]), df["spectral_centroid_hz"] * 1000.0])
# same KNN, but no scaling this time
clf_raw = KNeighborsClassifier(K)
Xr_train, Xr_test, yr_train, yr_test = train_test_split(X_raw, y, test_size=0.4, random_state=42)
clf_raw.fit(Xr_train, yr_train)
score = clf_raw.score(Xr_test, yr_test)
print("The mean accuracy without scaling is %f" % score)
# plot the decision boundary as a background
ax = plt.subplot()
DecisionBoundaryDisplay.from_estimator(clf_raw, X_raw, cmap='PiYG', alpha=0.8, ax=ax, eps=0.5)
ax.scatter(X_raw[:, 0], X_raw[:, 1], c=y, cmap='PiYG', alpha=0.6, edgecolors="k")
ax.set_xlabel("log10(STA/LTA)")
ax.set_ylabel("Spectral centroid (mHz)")
plt.show()The mean accuracy without scaling is 0.865000

La exactitud cae y la frontera de decisión es ahora un conjunto de bandas horizontales: el clasificador está usando, en la práctica, solo el centroide. Los clasificadores basados en distancias ven solo números, no unidades. La distancia euclidiana en el plano (log-STA/LTA, milihercios) está dominada por el eje de los milihercios, así que la información del STA/LTA se ignora. Escalar con StandardScaler pone ambas características en pie de igualdad, y por eso normalizamos antes de ajustar.
2. Métricas de desempeño del clasificador¶
En un clasificador binario, designamos una de las dos clases como positiva y la otra como negativa. Consideremos N muestras de datos.
| Clase verdadera \ Clase predicha | Negativa | Positiva | Total |
|---|---|---|---|
| Negativa | Verdadero negativo | Falso positivo | n |
| Positiva | Falso negativo | Verdadero positivo | p |
| Total | n’ | p’ | N |
En este ejemplo, había originalmente un total de etiquetas positivas y etiquetas negativas. Terminamos con muestras predichas como positivas y predichas como negativas.
Verdadero positivo TP (true positive): el número de datos predichos como positivos que eran originalmente positivos.
Verdadero negativo TN (true negative): el número de datos predichos como negativos que eran originalmente negativos.
Falso positivo FP (false positive): el número de datos predichos como positivos pero que eran originalmente negativos.
Falso negativo FN (false negative): el número de datos predichos como negativos pero que eran originalmente positivos.
Matriz de confusión:
Cuenta las instancias en que un elemento de la clase A se clasifica en la clase B: la entrada es el número de muestras de clase verdadera predichas como clase . scikit-learn ordena las clases por etiqueta, así que con ruido = 0 (negativa) y evento = 1 (positiva) la matriz que imprime confusion_matrix es:
La primera fila son las ventanas de ruido verdadero, la segunda fila los eventos verdaderos; las predicciones correctas quedan en la diagonal. La matriz de confusión puede extenderse a una clasificación multiclase, y la matriz es entonces de KxK en lugar de 2x2. La mejor matriz de confusión es la que se acerca a una matriz diagonal, con términos fuera de la diagonal pequeños.
Otras métricas de desempeño del modelo El desempeño del modelo puede evaluarse con lo siguiente:
Error: la fracción de los datos que fue mal clasificada
-> 0
Exactitud (accuracy): la fracción de los datos que fue clasificada correctamente:
--> 1
Tasa de TP: la fracción de las muestras verdaderamente positivas que se clasifican correctamente:
--> 1
Esta razón, , es también la exhaustividad (recall) o sensibilidad. La cantidad es la misma en los tres casos y la fórmula la fija; la palabra depende del campo. La recuperación de información y el aprendizaje automático dicen «exhaustividad»; la medicina, la epidemiología y la literatura de detección dicen «sensibilidad»; en la práctica hispanohablante de ML también se dice recall. Este libro escribe «exhaustividad» y glosa recall al usarla por primera vez en cada lección.
Tasa de TN: la fracción de las muestras verdaderamente negativas que se clasifican correctamente:
--> 1
Esta razón es también la especificidad.
Precisión (precision): la razón entre las muestras predichas en la clase positiva que eran en efecto positivas y el número total de muestras predichas como positivas.
--> 1
Score F1:
--> 1.
La media armónica da más peso al menor de los dos, así que el score F1 es alto solo si tanto la exhaustividad como la precisión son altas.
¿Cómo covarían la precisión y la exhaustividad?
def returnPrecisionAndRecall(TP, FP, TN, FN):
precision = TP / (TP + FP) if (TP + FP) else 0
recall = TP / (TP + FN) if (TP + FN) else 0
return {'precision': precision, 'recall': recall}
# Okay, set how many true and false values are in the original dataset
actualTrueValues = 10
actualFalseValues = 10
sumValues = actualTrueValues + actualFalseValues
precisionRecallCollector = []
# Now, run a set of simulations
for n in range(5000):
TP = 0
FP = 0
TN = 0
FN = 0
# Begin by randomly setting the total number of true values returned
totalTrue = np.random.randint(0, high=sumValues)
# The total number of false values is sumValues - totalTrue
totalFalse = sumValues - totalTrue
# Partition totalTrue and totalFalse into TP, FP and TN, FN, respectively
if totalTrue != 0:
TP = np.random.randint(0, high=totalTrue)
FP = totalTrue - TP
if totalFalse != 0:
TN = np.random.randint(0, high=totalFalse)
FN = totalFalse - TN
thisPandR = returnPrecisionAndRecall(TP, FP, TN, FN)
precisionRecallCollector.append([thisPandR['precision'], thisPandR['recall']])pAndRArray = np.asarray(precisionRecallCollector)
plt.scatter(pAndRArray[:,0], pAndRArray[:,1], color='k', alpha=0.25)
plt.xlabel('Precision')
plt.ylabel('Recall')
De lo anterior vemos que cabe esperar un rango amplio de comportamientos entre la precisión y la exhaustividad. Que un clasificador prediga bien una clase no significa que prediga bien la otra. Esto demuestra la importancia de reportar ambas métricas.
Imprimamos estas medidas para nuestro clasificador KNN con características escaladas, evaluado sobre el conjunto de prueba, usando scikit-learn.
from sklearn.metrics import confusion_matrix,precision_score,recall_score,f1_score
# Fit the model.
y_test_pred=clf.predict(X_test)
print("confusion matrix")
print(confusion_matrix(y_test,y_test_pred))
print("precison, recall")
print(precision_score(y_test,y_test_pred),recall_score(y_test,y_test_pred))
print("F1 score")
print(f1_score(y_test,y_test_pred))confusion matrix
[[384 12]
[ 20 384]]
precison, recall
0.9696969696969697 0.9504950495049505
F1 score
0.96
Un reporte completo y bien formateado del desempeño puede obtenerse con la función classification_report:
from sklearn.metrics import classification_report
print(f"Classification report for classifier {clf}:\n"
f"{classification_report(y_test, y_test_pred)}\n")Classification report for classifier KNeighborsClassifier():
precision recall f1-score support
0 0.95 0.97 0.96 396
1 0.97 0.95 0.96 404
accuracy 0.96 800
macro avg 0.96 0.96 0.96 800
weighted avg 0.96 0.96 0.96 800
La precisión y la exhaustividad se intercambian: aumentar la precisión reduce la exhaustividad.
El clasificador usa un valor de umbral para decidir si un dato pertenece a una clase. Aumentar el umbral da scores de precisión más altos; disminuir el umbral da scores de exhaustividad más altos. Veamos los distintos valores de los scores.
Característica Operativa del Receptor (ROC)
Grafica la tasa de verdaderos positivos contra la tasa de falsos positivos. La curva ROC es visual, pero podemos cuantificar el desempeño del clasificador con el área bajo la curva (AUC). Idealmente, el AUC es 1.

[fuente: https://
Calculamos la curva sobre el conjunto de prueba: ediciones anteriores de este libro graficaban la ROC del conjunto de entrenamiento, y eso favorece artificialmente al modelo.
y_scores = clf.predict_proba(X_test)
print(y_scores[:5])[[1. 0. ]
[1. 0. ]
[0. 1. ]
[0.8 0.2]
[0. 1. ]]
from sklearn.metrics import roc_curve
fpr,tpr,thresholds=roc_curve(y_test,y_scores[:,1])
plt.plot(fpr,tpr,linewidth=2);plt.grid(True)
plt.xlabel('False Positive Rate')
plt.ylabel('True Positive Rate')
plt.plot([0,1],[0,1],'k--')
Ahora exploramos los distintos clasificadores empaquetados en scikit-learn. Podemos probar sistemáticamente su desempeño y guardar los scores de precisión, exhaustividad y F1.
Una advertencia sobre una de las entradas: el entrenamiento de GaussianProcessClassifier escala como con el número de muestras. Aquí, con , no hay problema, pero es una trampa en los catálogos de detección reales, que rutinariamente contienen de 105 a 106 ventanas.
from sklearn.svm import SVC
from sklearn.gaussian_process import GaussianProcessClassifier
from sklearn.gaussian_process.kernels import RBF
from sklearn.tree import DecisionTreeClassifier
from sklearn.ensemble import RandomForestClassifier, AdaBoostClassifier
from sklearn.naive_bayes import GaussianNB
from sklearn.discriminant_analysis import QuadraticDiscriminantAnalysis
# define models
names = [
"Nearest Neighbors",
"Linear SVM",
"RBF SVM",
"Gaussian Process",
"Decision Tree",
"Random Forest",
"AdaBoost",
"Naive Bayes",
"QDA",
]
classifiers = [
KNeighborsClassifier(7),
SVC(kernel="linear", C=0.025),
SVC(gamma=2, C=1),
GaussianProcessClassifier(1.0 * RBF(1.0)),
DecisionTreeClassifier(max_depth=5),
RandomForestClassifier(max_depth=5, n_estimators=10, max_features=1),
AdaBoostClassifier(),
GaussianNB(),
QuadraticDiscriminantAnalysis(),
]3. Exploración de modelos¶
Explore cómo se desempeña cada uno de estos modelos sobre los datos sintéticos.
Guarde en un arreglo los valores de precisión, exhaustividad y score F1.
Encuentre el modelo con mejor desempeño
# reuse the train/test split and train-fit scaling from section 1
pre=np.zeros(len(classifiers))
rec=np.zeros(len(classifiers))
f1=np.zeros(len(classifiers))
for ii,iclass in enumerate(classifiers):
iclass.fit(X_train, y_train)
y_test_pred=iclass.predict(X_test)
pre[ii] =precision_score(y_test,y_test_pred)
rec[ii] =recall_score(y_test,y_test_pred)
f1[ii] =f1_score(y_test,y_test_pred)
df_scores=pd.DataFrame({'CLF name':names,'precision':pre,'recall':rec,'f1_score':f1})
print(df_scores) CLF name precision recall f1_score
0 Nearest Neighbors 0.969620 0.948020 0.958698
1 Linear SVM 0.969309 0.938119 0.953459
2 RBF SVM 0.967005 0.943069 0.954887
3 Gaussian Process 0.971939 0.943069 0.957286
4 Decision Tree 0.966837 0.938119 0.952261
5 Random Forest 0.962500 0.952970 0.957711
6 AdaBoost 0.959900 0.948020 0.953923
7 Naive Bayes 0.966921 0.940594 0.953576
8 QDA 0.966921 0.940594 0.953576
Detección desbalanceada: cuando la exactitud deja de significar algo¶
El conjunto balanceado de arriba era una ficción de aula conveniente. En una estación real, las ventanas con evento son raras: la mayor parte del día es ruido. Ahora extraemos 5000 ventanas con una fracción de eventos del 2 por ciento, unos 100 eventos entre 4900 ventanas de ruido — un desbalance de 1:50. Las mismas características, la misma transformación logarítmica, el mismo escalado. La división es estratificada para que la clase rara conserve la misma proporción en entrenamiento y prueba.
df_imb = mlgeo_synth.detector_features(n=5000, event_fraction=0.02, seed=42)
print(df_imb["label"].value_counts())
X_imb = np.column_stack([np.log10(df_imb["sta_lta"]), df_imb["spectral_centroid_hz"]])
y_imb = df_imb["label"].values
Xi_train, Xi_test, yi_train, yi_test = train_test_split(
X_imb, y_imb, stratify=y_imb, test_size=0.4, random_state=42
)
# same discipline as before: the scaler is fit on the training split only
scaler_imb = StandardScaler().fit(Xi_train)
Xi_train = scaler_imb.transform(Xi_train)
Xi_test = scaler_imb.transform(Xi_test)
baseline = DummyClassifier(strategy="most_frequent").fit(Xi_train, yi_train)
print(f"Majority-class baseline accuracy: {baseline.score(Xi_test, yi_test):.3f}")label
0 4900
1 100
Name: count, dtype: int64
Majority-class baseline accuracy: 0.980
El predictor trivial que nunca declara un evento ya alcanza 0.98. La exactitud es ahora inútil como métrica. Entrenemos dos clasificadores reales y miremos las métricas que importan para la clase de los eventos.
from sklearn.linear_model import LogisticRegression
from sklearn.metrics import accuracy_score
logreg = LogisticRegression().fit(Xi_train, yi_train)
forest = RandomForestClassifier(n_estimators=100, random_state=42).fit(Xi_train, yi_train)
for name, model in [("Logistic regression", logreg), ("Random forest", forest)]:
yp = model.predict(Xi_test)
print(name)
print(f" accuracy : {accuracy_score(yi_test, yp):.3f}")
print(f" precision (event): {precision_score(yi_test, yp):.3f}")
print(f" recall (event) : {recall_score(yi_test, yp):.3f}")
print(f" F1 (event) : {f1_score(yi_test, yp):.3f}")Logistic regression
accuracy : 0.992
precision (event): 0.929
recall (event) : 0.650
F1 (event) : 0.765
Random forest
accuracy : 0.992
precision (event): 0.853
recall (event) : 0.725
F1 (event) : 0.784
Ambos modelos superan la exactitud del modelo de referencia por menos de un punto porcentual, y sin embargo pierden aproximadamente un tercio de los eventos. La exhaustividad sobre la clase de los eventos es el número que le importa a la analista: cuenta los sismos que el detector realmente atrapa. El intercambio es operativo. Un detector afinado para alta exhaustividad inunda a la analista con disparos falsos por revisar; un detector afinado para alta precisión mantiene limpia la lista de disparos pero pierde los eventos pequeños. Dónde ubicarse sobre esa curva es una decisión sobre el tiempo de la analista y el costo científico, no un detalle de modelado.
Podemos ver el intercambio completo de una sola vez graficando la curva ROC y la curva de precisión–exhaustividad lado a lado.
from sklearn.metrics import PrecisionRecallDisplay, RocCurveDisplay
fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(11, 4.5))
for name, model in [("Logistic regression", logreg), ("Random forest", forest)]:
RocCurveDisplay.from_estimator(model, Xi_test, yi_test, ax=ax1, name=name)
PrecisionRecallDisplay.from_estimator(model, Xi_test, yi_test, ax=ax2, name=name)
ax1.plot([0, 1], [0, 1], 'k--', linewidth=1)
ax1.set_title("ROC curve")
ax2.set_title("Precision-recall curve")
ax1.grid(True)
ax2.grid(True)
plt.tight_layout()
plt.show()
La curva ROC se ve excelente, y ese es exactamente el problema. Con un desbalance de 1:50, la tasa de falsos positivos divide entre las aproximadamente 1960 ventanas de ruido del conjunto de prueba, así que se mantiene diminuta incluso cuando los falsos positivos superan en número a las detecciones verdaderas. La curva de precisión–exhaustividad divide, en cambio, entre el número de positivos predichos, así que sigue la carga de trabajo real de la analista: cae en cuanto la lista de disparos se llena de ruido. Cuando la clase positiva es rara, lea la curva PR, no la ROC.
Tratar el desbalance¶
Diagnosticar el desbalance es solo la mitad del trabajo. Dos perillas estándar mueven el punto de operación sin recolectar datos nuevos:
class_weight="balanced"repondera la pérdida de entrenamiento de modo que cada error sobre un evento cuente tanto como unos 49 errores sobre ruido, empujando la frontera ajustada hacia la clase rara.- Desplazamiento del umbral:
predictcorta la probabilidad predicha en 0.5, pero nada obliga a esa elección. Bajar el umbral intercambia precisión por exhaustividad a lo largo de la curva que acabamos de graficar.
logreg_bal = LogisticRegression(class_weight="balanced").fit(Xi_train, yi_train)
# threshold moving: same model as before, decision cut lowered from 0.5 to 0.25
proba = logreg.predict_proba(Xi_test)[:, 1]
yp_thresh = (proba >= 0.25).astype(int)
rows = []
for name, yp in [
("logreg, threshold 0.5", logreg.predict(Xi_test)),
("logreg, class_weight='balanced'", logreg_bal.predict(Xi_test)),
("logreg, threshold 0.25", yp_thresh),
]:
rows.append({"model": name,
"precision (event)": precision_score(yi_test, yp),
"recall (event)": recall_score(yi_test, yp),
"F1 (event)": f1_score(yi_test, yp)})
pd.DataFrame(rows).round(3)Ambas perillas elevan la exhaustividad sobre los eventos pagando con precisión: la lista de disparos se alarga, pero se pierden menos eventos. Ninguna es «mejor» — ambas se deslizan a lo largo de la misma curva de precisión–exhaustividad, y la elección del punto de operación pertenece a la analista, no al modelo.
Cierre¶
La elección de la métrica es una decisión de diseño, no una ocurrencia tardía. Codifica el costo operativo de cada tipo de error: un evento perdido cuesta ciencia, un disparo falso cuesta tiempo de analista. Decida cuánto cuesta cada error en su aplicación, elija la métrica que le ponga precio, y solo entonces compare modelos.