2. Modelo final con validación segura

import os
# joblib no puede contar los núcleos físicos en este Windows y lanza un aviso; se le indica cuántos usar
os.environ["LOKY_MAX_CPU_COUNT"] = str(max(1, (os.cpu_count() or 2) // 2))

import json
import time
import warnings
from pathlib import Path

import joblib
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import matplotlib.ticker as mtick
from matplotlib.colors import LinearSegmentedColormap
from scipy import stats
from IPython.display import display

import sklearn
from sklearn.compose import ColumnTransformer
from sklearn.impute import SimpleImputer
from sklearn.preprocessing import MinMaxScaler, OneHotEncoder
from sklearn.pipeline import Pipeline
from sklearn.model_selection import (train_test_split, GridSearchCV, StratifiedKFold,
                                     cross_validate, cross_val_predict)
from sklearn.metrics import (roc_auc_score, accuracy_score, precision_score, recall_score, f1_score,
                             confusion_matrix, classification_report, roc_curve)
from sklearn.inspection import permutation_importance
from sklearn.svm import SVC
from sklearn.linear_model import LogisticRegression
from sklearn.neighbors import KNeighborsClassifier
from sklearn.ensemble import RandomForestClassifier, GradientBoostingClassifier

warnings.filterwarnings("ignore")
pd.set_option("display.max_columns", 60)
pd.set_option("display.float_format", "{:.4f}".format)

SEED = 42
print(f"pandas {pd.__version__} | numpy {np.__version__} | scikit-learn {sklearn.__version__} | joblib {joblib.__version__}")

COLOR = {
    "principal":  "#2a78d6",
    "secundario": "#eb6834",
    "texto":      "#0b0b0b",
    "texto_2":    "#52514e",
    "tenue":      "#898781",
    "rejilla":    "#e1e0d9",
    "eje":        "#c3c2b7",
}

plt.rcParams.update({
    "font.family": "sans-serif",
    "font.sans-serif": ["Segoe UI", "Calibri", "Arial", "DejaVu Sans"],
    "font.size": 10,
    "axes.titlesize": 11,
    "axes.titleweight": "semibold",
    "axes.titlelocation": "left",
    "axes.titlepad": 10,
    "axes.labelsize": 10,
    "axes.labelcolor": COLOR["texto_2"],
    "axes.edgecolor": COLOR["eje"],
    "axes.linewidth": 0.8,
    "axes.spines.top": False,
    "axes.spines.right": False,
    "axes.grid": True,
    "axes.axisbelow": True,
    "grid.color": COLOR["rejilla"],
    "grid.linewidth": 0.6,
    "xtick.color": COLOR["texto_2"],
    "ytick.color": COLOR["texto_2"],
    "xtick.labelsize": 9,
    "ytick.labelsize": 9,
    "text.color": COLOR["texto"],
    "figure.facecolor": "white",
    "axes.facecolor": "white",
    "figure.dpi": 110,
    "savefig.dpi": 150,
    "savefig.bbox": "tight",
    "legend.frameon": False,
    "legend.fontsize": 9,
})

CLASE_COLOR = {0: COLOR["principal"], 1: COLOR["secundario"]}
CLASE_NOMBRE = {0: "Sin enfermedad", 1: "Con enfermedad"}
CMAP_SEC = LinearSegmentedColormap.from_list("azul", ["#f7f9fc", "#9ec5f4", "#2a78d6", "#0d366b"])


def titulo_figura(fig, texto):
    fig.suptitle(texto, x=0.01, ha="left", fontsize=13, fontweight="semibold")

ROOT = Path.cwd().parent if Path.cwd().name == "notebooks" else Path.cwd()
DATA_DIR = ROOT / "data"
APP_DIR = ROOT / "app"
FIG_DIR = ROOT / "figuras"
FIG_DIR.mkdir(exist_ok=True)

df = pd.read_csv(DATA_DIR / "heart.csv")

TARGET = "HeartDisease"
X = df.drop(columns=TARGET)
y = df[TARGET]
print(f"{len(df)} pacientes | {X.shape[1]} variables | {y.mean():.1%} con enfermedad")
display(X.head(3))

X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, stratify=y, random_state=SEED)

display(pd.DataFrame({
    "n": [len(y_train), len(y_test)],
    "pct_con_enfermedad": [y_train.mean() * 100, y_test.mean() * 100],
    "ceros_cholesterol": [(X_train["Cholesterol"] == 0).sum(), (X_test["Cholesterol"] == 0).sum()],
}, index=["train", "test"]).style.format({"pct_con_enfermedad": "{:.2f} %"}))

pd.concat([X_train, y_train], axis=1).to_csv(DATA_DIR / "train.csv", index=False)
pd.concat([X_test, y_test], axis=1).to_csv(DATA_DIR / "test.csv", index=False)
print(f"Guardados: {DATA_DIR / 'train.csv'} y {DATA_DIR / 'test.csv'}")

COLS_NUM = ["Age", "MaxHR", "Oldpeak"]
comparacion = pd.DataFrame({
    "min_train": X_train[COLS_NUM].min(), "min_todo": X[COLS_NUM].min(),
    "max_train": X_train[COLS_NUM].max(), "max_todo": X[COLS_NUM].max(),
})
comparacion["¿igual?"] = np.where((comparacion["min_train"] == comparacion["min_todo"])
                                  & (comparacion["max_train"] == comparacion["max_todo"]), "Sí", "No")
display(comparacion)

# Consecuencia: en test pueden aparecer valores fuera del rango [0, 1] después de escalar, y eso es correcto
escalador = MinMaxScaler().fit(X_train[COLS_NUM])
X_test_esc = pd.DataFrame(escalador.transform(X_test[COLS_NUM]), columns=COLS_NUM)
print("Rango de test después de escalar con los parámetros de train:")
display(X_test_esc.agg(["min", "max"]))

COLS_ESCALAR = ["Age", "MaxHR", "Oldpeak"]
COLS_CEROS = ["RestingBP", "Cholesterol"]
COLS_ONEHOT = ["Sex", "ChestPainType", "RestingECG", "ExerciseAngina", "ST_Slope"]
COLS_BINARIAS = ["FastingBS"]
COLUMNAS_ENTRADA = list(X.columns)          # orden en que el modelo espera las variables


def crear_preprocesador():
    return ColumnTransformer([
        ("num", MinMaxScaler(), COLS_ESCALAR),
        ("num_ceros", Pipeline([
            ("imputar", SimpleImputer(missing_values=0, strategy="median", add_indicator=True)),
            ("escalar", MinMaxScaler()),
        ]), COLS_CEROS),
        ("cat", OneHotEncoder(handle_unknown="ignore", sparse_output=False), COLS_ONEHOT),
        ("bin", "passthrough", COLS_BINARIAS),
    ], verbose_feature_names_out=False)


prep = crear_preprocesador().fit(X_train)
print(f"Variables de salida ({len(prep.get_feature_names_out())}):")
print(list(prep.get_feature_names_out()))
print("Medianas aprendidas con train:", dict(zip(COLS_CEROS, prep.named_transformers_["num_ceros"]["imputar"].statistics_)))

pipe = Pipeline([("prep", crear_preprocesador()), ("clf", LogisticRegression())])

param_grid = [
    {"clf": [LogisticRegression(max_iter=5000, solver="liblinear", random_state=SEED)],
     "clf__C": [0.01, 0.1, 1, 10], "clf__penalty": ["l1", "l2"]},
    {"clf": [SVC(probability=True, random_state=SEED)],
     "clf__C": [0.1, 1, 10, 100], "clf__gamma": ["scale", 0.01, 0.1]},
    {"clf": [KNeighborsClassifier()],
     "clf__n_neighbors": [5, 11, 21, 31, 51, 71, 101], "clf__weights": ["uniform", "distance"], "clf__p": [1, 2]},
    {"clf": [RandomForestClassifier(n_estimators=300, random_state=SEED)],
     "clf__max_depth": [None, 5, 10], "clf__min_samples_leaf": [1, 3, 5], "clf__max_features": ["sqrt", 0.5]},
    {"clf": [GradientBoostingClassifier(random_state=SEED)],
     "clf__n_estimators": [100, 300], "clf__learning_rate": [0.01, 0.03, 0.1],
     "clf__max_depth": [2, 3], "clf__subsample": [0.8, 1.0]},
]

CV = StratifiedKFold(n_splits=5, shuffle=True, random_state=SEED)
grid = GridSearchCV(pipe, param_grid, cv=CV, scoring={"AUC": "roc_auc", "Accuracy": "accuracy"},
                    refit="AUC", n_jobs=-1, return_train_score=True)

t0 = time.time()
grid.fit(X_train, y_train)
print(f"{len(grid.cv_results_['params'])} combinaciones x 5 folds = {len(grid.cv_results_['params']) * 5} ajustes "
      f"en {time.time() - t0:.1f} s")
print(f"Mejor modelo: {type(grid.best_params_['clf']).__name__}")
print(f"Mejores hiperparámetros: { {k.replace('clf__', ''): v for k, v in grid.best_params_.items() if k != 'clf'} }")
print(f"AUC CV del mejor: {grid.best_score_:.4f}")

resultados = pd.DataFrame(grid.cv_results_)
resultados["modelo"] = resultados["param_clf"].apply(lambda m: type(m).__name__)
resultados["hiperparametros"] = resultados["params"].apply(
    lambda p: {k.replace("clf__", ""): v for k, v in p.items() if k != "clf"})

mejores = (resultados.sort_values("mean_test_AUC", ascending=False).groupby("modelo").head(1)
           .set_index("modelo")[["mean_test_AUC", "std_test_AUC", "mean_test_Accuracy", "std_test_Accuracy",
                                 "mean_train_AUC", "hiperparametros"]])
mejores["brecha_train_CV"] = mejores["mean_train_AUC"] - mejores["mean_test_AUC"]
mejores.insert(0, "puesto", range(1, len(mejores) + 1))
display(mejores.drop(columns="mean_train_AUC").style.format(
    {c: "{:.4f}" for c in ["mean_test_AUC", "std_test_AUC", "mean_test_Accuracy", "std_test_Accuracy", "brecha_train_CV"]}))

# AUC de cada fold para el mejor de cada familia (los folds son los mismos para todos)
folds = pd.DataFrame({m: resultados.loc[resultados.loc[resultados["modelo"] == m, "mean_test_AUC"].idxmax(),
                                        [f"split{k}_test_AUC" for k in range(5)]].values.astype(float)
                      for m in mejores.index}, index=[f"fold {k + 1}" for k in range(5)])
display(folds.style.format("{:.4f}").highlight_max(axis=1, props=f"background-color: {COLOR['secundario']}; color: white"))

# ¿El ganador es significativamente mejor que el segundo? Prueba t pareada sobre los 5 folds
primero, segundo = mejores.index[:2]
t_stat, p_t = stats.ttest_rel(folds[primero], folds[segundo])
print(f"{primero} vs {segundo}: diferencia media = {(folds[primero] - folds[segundo]).mean():+.4f}, "
      f"t pareada p = {p_t:.3f}")

CV_EXTERNO = StratifiedKFold(n_splits=5, shuffle=True, random_state=SEED + 1)
grid_interno = GridSearchCV(pipe, param_grid, cv=CV, scoring="roc_auc", n_jobs=-1)

t0 = time.time()
anidada = cross_validate(grid_interno, X_train, y_train, cv=CV_EXTERNO, scoring=["roc_auc", "accuracy"],
                         return_estimator=True, n_jobs=1)
print(f"Validación anidada en {time.time() - t0:.1f} s")

tabla_anidada = pd.DataFrame({
    "AUC": anidada["test_roc_auc"],
    "Accuracy": anidada["test_accuracy"],
    "modelo_elegido": [type(e.best_params_["clf"]).__name__ for e in anidada["estimator"]],
}, index=[f"fold externo {k + 1}" for k in range(5)])
display(tabla_anidada.style.format({"AUC": "{:.4f}", "Accuracy": "{:.4f}"}))
print(f"AUC anidado     : {anidada['test_roc_auc'].mean():.4f} ± {anidada['test_roc_auc'].std():.4f}")
print(f"Accuracy anidado: {anidada['test_accuracy'].mean():.4f} ± {anidada['test_accuracy'].std():.4f}")
print(f"AUC de la búsqueda simple (Parte 5): {grid.best_score_:.4f}")

modelo_final = grid.best_estimator_
proba_test = modelo_final.predict_proba(X_test)[:, 1]
pred_test = (proba_test >= 0.5).astype(int)

tn, fp, fn, tp = confusion_matrix(y_test, pred_test).ravel()
metricas_test = {
    "AUC": roc_auc_score(y_test, proba_test),
    "Accuracy": accuracy_score(y_test, pred_test),
    "Precision": precision_score(y_test, pred_test),
    "Recall (sensibilidad)": recall_score(y_test, pred_test),
    "Especificidad": tn / (tn + fp),
    "F1": f1_score(y_test, pred_test),
}

# Intervalo de confianza del AUC con bootstrap (2,000 remuestreos del conjunto de prueba)
rng = np.random.default_rng(SEED)
auc_boot = []
for _ in range(2000):
    idx = rng.integers(0, len(y_test), len(y_test))
    if y_test.iloc[idx].nunique() == 2:
        auc_boot.append(roc_auc_score(y_test.iloc[idx], proba_test[idx]))
ic_inf, ic_sup = np.percentile(auc_boot, [2.5, 97.5])

display(pd.Series(metricas_test, name=type(modelo_final.named_steps["clf"]).__name__).to_frame().T
        .style.format("{:.4f}"))
print(f"AUC test = {metricas_test['AUC']:.4f}  (IC 95 % bootstrap: {ic_inf:.4f} – {ic_sup:.4f})")
print(f"Verdaderos negativos = {tn} | Falsos positivos = {fp} | Falsos negativos = {fn} | Verdaderos positivos = {tp}")
print()
print(classification_report(y_test, pred_test, target_names=[CLASE_NOMBRE[0], CLASE_NOMBRE[1]], digits=4))

cm = confusion_matrix(y_test, pred_test)
cm_pct = cm / cm.sum(axis=1, keepdims=True) * 100
etiquetas = [CLASE_NOMBRE[0], CLASE_NOMBRE[1]]
nombres_celda = [["Verdaderos negativos", "Falsos positivos"], ["Falsos negativos", "Verdaderos positivos"]]

fig, ax = plt.subplots(figsize=(6.5, 5.5))
ax.imshow(cm_pct, cmap=CMAP_SEC, vmin=0, vmax=100)
for i in range(2):
    for j in range(2):
        claro = cm_pct[i, j] > 55
        ax.text(j, i - 0.12, f"{cm[i, j]}", ha="center", va="center", fontsize=20, fontweight="semibold",
                color="white" if claro else COLOR["texto"])
        ax.text(j, i + 0.13, f"{cm_pct[i, j]:.1f} % de la fila", ha="center", va="center", fontsize=9,
                color="white" if claro else COLOR["texto_2"])
        ax.text(j, i + 0.3, nombres_celda[i][j], ha="center", va="center", fontsize=8.5,
                color="white" if claro else COLOR["texto_2"])
ax.set_xticks([0, 1], etiquetas)
ax.set_yticks([0, 1], etiquetas)
ax.set_xlabel("Predicción del modelo")
ax.set_ylabel("Diagnóstico real")
ax.grid(False)
for borde in ax.spines.values():
    borde.set_visible(False)
ax.set_title(f"Matriz de confusión en test (umbral 0.5)  ·  Accuracy = {metricas_test['Accuracy']:.3f}")
fig.tight_layout()
fig.savefig(FIG_DIR / "02_matriz_confusion.png")
plt.show()

fpr, tpr, umbrales = roc_curve(y_test, proba_test)
i05 = np.argmin(np.abs(umbrales - 0.5))

fig, ax = plt.subplots(figsize=(6.8, 6.3))
ax.fill_between(fpr, tpr, step="post", color=COLOR["principal"], alpha=0.12)
ax.step(fpr, tpr, where="post", color=COLOR["principal"], linewidth=2.2,
        label=f"{type(modelo_final.named_steps['clf']).__name__} (AUC = {metricas_test['AUC']:.3f}; "
              f"IC 95 % {ic_inf:.3f}–{ic_sup:.3f})")
ax.plot([0, 1], [0, 1], color=COLOR["tenue"], linestyle=":", linewidth=1, label="Azar (AUC = 0.5)")
ax.scatter(fpr[i05], tpr[i05], s=60, color=COLOR["secundario"], zorder=3,
           label=f"Umbral 0.5: sensibilidad {tpr[i05]:.3f}, especificidad {1 - fpr[i05]:.3f}")
ax.set_xlabel("Tasa de falsos positivos (1 − especificidad)")
ax.set_ylabel("Tasa de verdaderos positivos (sensibilidad)")
ax.set_xlim(0, 1)
ax.set_ylim(0, 1.01)
ax.set_aspect("equal")
ax.legend(loc="lower right", fontsize=8.5)
ax.set_title("Curva ROC en el conjunto de prueba")
fig.tight_layout()
fig.savefig(FIG_DIR / "02_curva_roc.png")
plt.show()

proba_oof = cross_val_predict(modelo_final, X_train, y_train, cv=CV, method="predict_proba", n_jobs=-1)[:, 1]
fpr_tr, tpr_tr, umb_tr = roc_curve(y_train, proba_oof)
umbral_youden = float(umb_tr[np.argmax(tpr_tr - fpr_tr)])
print(f"Umbral de Youden (elegido con train, fuera de fold): {umbral_youden:.3f}")

filas = []
for u in sorted({0.3, 0.4, 0.5, 0.6, 0.7, round(umbral_youden, 3)}):
    p = (proba_test >= u).astype(int)
    tn_u, fp_u, fn_u, tp_u = confusion_matrix(y_test, p).ravel()
    filas.append({"umbral": u, "origen": "Youden (train)" if u == round(umbral_youden, 3) else "",
                  "sensibilidad": tp_u / (tp_u + fn_u), "especificidad": tn_u / (tn_u + fp_u),
                  "accuracy": (tp_u + tn_u) / len(p), "falsos_positivos": fp_u, "falsos_negativos": fn_u})
display(pd.DataFrame(filas).set_index("umbral").style.format(
    {"sensibilidad": "{:.3f}", "especificidad": "{:.3f}", "accuracy": "{:.3f}"}))

imp = permutation_importance(modelo_final, X_test, y_test, scoring="roc_auc", n_repeats=30, random_state=SEED, n_jobs=-1)
importancias = (pd.DataFrame({"caida_AUC_media": imp.importances_mean, "caida_AUC_std": imp.importances_std},
                             index=X_test.columns).sort_values("caida_AUC_media"))

fig, ax = plt.subplots(figsize=(8, 5))
colores = np.where(importancias["caida_AUC_media"] > 0.01, COLOR["principal"], COLOR["tenue"])
ax.barh(importancias.index.str.lower(), importancias["caida_AUC_media"], xerr=importancias["caida_AUC_std"],
        color=colores, height=0.65, error_kw=dict(ecolor=COLOR["texto_2"], capsize=2, linewidth=0.8))
ax.axvline(0, color=COLOR["texto_2"], linewidth=0.8)
ax.set_xlabel("Caída del AUC al desordenar la variable")
ax.grid(axis="y", visible=False)
ax.set_title("Importancia por permutación en test  ·  en gris: caída menor que 0.01")
fig.tight_layout()
fig.savefig(FIG_DIR / "02_importancia_permutacion.png")
plt.show()
display(importancias.sort_values("caida_AUC_media", ascending=False).style.format("{:.4f}"))

APP_DIR.mkdir(exist_ok=True)
ruta_modelo = APP_DIR / "model.joblib"
joblib.dump(grid.best_estimator_, ruta_modelo)

metadatos = {
    "modelo": type(modelo_final.named_steps["clf"]).__name__,
    "hiperparametros": {k.replace("clf__", ""): (v if isinstance(v, (int, float, str, type(None))) else str(v))
                        for k, v in grid.best_params_.items() if k != "clf"},
    "columnas_entrada": COLUMNAS_ENTRADA,
    "categorias": {c: sorted(X_train[c].unique().tolist()) for c in COLS_ONEHOT},
    # Rango observado en train (sin los ceros imposibles); sirve para validar entradas en la API
    "rangos_train": {c: [float(X_train.loc[X_train[c] != 0, c].min() if c in COLS_CEROS else X_train[c].min()),
                         float(X_train[c].max())] for c in COLS_ESCALAR + COLS_CEROS},
    "umbral": 0.5,
    "umbral_youden": round(umbral_youden, 4),
    "auc_cv": round(float(grid.best_score_), 4),
    "auc_cv_anidado": round(float(anidada["test_roc_auc"].mean()), 4),
    "metricas_test": {k: round(float(v), 4) for k, v in metricas_test.items()},
    "sklearn_version": sklearn.__version__,
    "fecha_entrenamiento": time.strftime("%Y-%m-%d"),
}
(APP_DIR / "model_metadata.json").write_text(json.dumps(metadatos, indent=2, ensure_ascii=False), encoding="utf-8")

print(f"Modelo guardado   : {ruta_modelo} ({ruta_modelo.stat().st_size / 1e3:.1f} KB)")
print(f"Metadatos guardados: {APP_DIR / 'model_metadata.json'}")
print(json.dumps(metadatos, indent=2, ensure_ascii=False))

modelo_cargado = joblib.load(ruta_modelo)

# Las predicciones del modelo cargado son idénticas a las del modelo en memoria
assert np.allclose(modelo_cargado.predict_proba(X_test)[:, 1], proba_test)
print("El modelo cargado reproduce exactamente las probabilidades de test.")

# Paciente de ejemplo con los datos crudos (incluido un colesterol en 0, que el Pipeline trata como faltante)
paciente = pd.DataFrame([{
    "Age": 58, "Sex": "M", "ChestPainType": "ASY", "RestingBP": 140, "Cholesterol": 0,
    "FastingBS": 1, "RestingECG": "ST", "MaxHR": 110, "ExerciseAngina": "Y", "Oldpeak": 2.0, "ST_Slope": "Flat",
}])[COLUMNAS_ENTRADA]
proba = modelo_cargado.predict_proba(paciente)[0, 1]
print(f"Probabilidad de enfermedad cardíaca: {proba:.3f} -> predicción: {int(proba > 0.5)}")