"""Chapitre 37 — Apprentissage automatique et explicabilité : reproductibilité et bornes de performance des chiffres cités dans le
texte, les tableaux et les solutions. Modèle hédonique de référence : traite.cas1.hedonique (non redéfini). Tolérances déclarées."""
import sys, pathlib, math
import numpy as np, pytest
sys.path.insert(0, str(pathlib.Path(__file__).resolve().parents[1] / "cases"))
import ch37_ml as c37
import ch18_hedonique as c18
from traite import cas1

@pytest.fixture(scope="module")
def res():
    return c37.main(verbose=False)

def test_plis_identiques_au_chapitre_18():
    """Les dix plis aléatoires sont exactement ceux de ch18_hedonique.validation_croisee (même graine, même découpage)."""
    v = cas1.charger()["ventes"]; pa = c37.plis_aleatoires(len(v))
    rng = np.random.default_rng(c18.GRAINE); attendu = np.array_split(rng.permutation(len(v)), 10)
    assert all(np.array_equal(a, b) for a, b in zip(pa, attendu)) and sum(len(f) for f in pa) == 600
    assert c37.GRAINE == c18.GRAINE == 20261031

def test_reproductibilite_et_reference(res):
    """Le modèle hédonique validé ici reproduit le ch. 18 : RMSE 40 888 $, COD 6,83 en validation croisée ; prédiction du sujet 637 953 $."""
    m = res["comparaison"]["metriques"]["Hédonique (référence)"]["aleatoire"]
    assert m["COD"] == pytest.approx(6.83, abs=0.01) and m["rmse_log"] == pytest.approx(0.0839, abs=0.0005) and m["PRD"] == pytest.approx(1.007, abs=0.001)
    assert m["mdape"] == pytest.approx(5.68, abs=0.05) and m["r2_hors"] == pytest.approx(0.846, abs=0.001)
    assert res["modele_ref"]["r2"] == pytest.approx(0.853, abs=0.001) and res["modele_ref"]["sigma"] == pytest.approx(0.083, abs=0.0005)
    s = res["sujet_pred"]
    assert s.loc["Hédonique (référence)", "prediction_cas1"] == pytest.approx(637_953, abs=1)
    assert s.loc["Hédonique (référence)", "ln_pred"] == pytest.approx(13.3627, abs=0.0001)
    # une seconde exécution des modèles stochastiques donne les mêmes prédictions (graine fixée)
    Xm = res["Xm"]; y = res["y"]; f = res["comparaison"]["plis"]["aleatoire"][0]; tr = np.setdiff1d(np.arange(600), f)
    p1 = c37._boosting().fit(Xm.iloc[tr], y[tr]).predict(Xm.iloc[f]); p2 = c37._boosting().fit(Xm.iloc[tr], y[tr]).predict(Xm.iloc[f])
    assert np.allclose(p1, p2) and np.allclose(p1, res["comparaison"]["pred"]["Gradient boosting"]["aleatoire"][f])

def test_comparaison_aleatoire(res):
    """Sur un jeu engendré par une loi semi-log additive, l'hédonique bien spécifié domine ; le boosting est le meilleur des modèles d'arbres."""
    M = {k: v["aleatoire"] for k, v in res["comparaison"]["metriques"].items()}
    attendu = {"Hédonique (référence)": (0.0839, 6.83, 5.68), "Arbre (profondeur 3)": (0.1468, 11.63, 9.85), "Forêt aléatoire": (0.1067, 8.64, 7.52),
               "Gradient boosting": (0.0931, 7.42, 6.22), "XGBoost": (0.0921, 7.33, 6.09)}
    for nom, (rmse, cod, mdape) in attendu.items():
        assert M[nom]["rmse_log"] == pytest.approx(rmse, abs=0.0015), nom
        assert M[nom]["COD"] == pytest.approx(cod, abs=0.12), nom
        assert M[nom]["mdape"] == pytest.approx(mdape, abs=0.15), nom
    assert M["Hédonique (référence)"]["rmse_log"] < M["XGBoost"]["rmse_log"] < M["Gradient boosting"]["rmse_log"] < M["Forêt aléatoire"]["rmse_log"] < M["Arbre (profondeur 3)"]["rmse_log"]
    # bornes : tous les modèles d'ensemble restent dans la plage IAAO 5-10 ; l'arbre seul en sort
    assert 5.0 < M["Gradient boosting"]["COD"] < 10.0 and 5.0 < M["Forêt aléatoire"]["COD"] < 10.0 and M["Arbre (profondeur 3)"]["COD"] > 10.0
    # régressivité des arbres : PRB plus négatif que l'hédonique, PRD plus grand
    assert M["Forêt aléatoire"]["PRB"] < -0.15 and M["Forêt aléatoire"]["PRD"] > M["Hédonique (référence)"]["PRD"]
    assert M["Forêt aléatoire"]["PRB"] == pytest.approx(-0.196, abs=0.01) and M["Gradient boosting"]["PRB"] == pytest.approx(-0.088, abs=0.01)
    # niveau : médianes des ratios ≈ 1 pour tous
    assert all(0.99 < m["ratio_median"] < 1.015 for m in M.values())
    # RMSE en dollars
    assert M["Hédonique (référence)"]["rmse_dollars"] == pytest.approx(40_889, abs=30) and M["Gradient boosting"]["rmse_dollars"] == pytest.approx(45_644, abs=300)

def test_validation_spatiale_et_temporelle(res):
    M = res["comparaison"]["metriques"]
    assert M["Hédonique (référence)"]["spatial"]["COD"] == pytest.approx(7.00, abs=0.1) and M["Gradient boosting"]["spatial"]["COD"] == pytest.approx(7.62, abs=0.12)
    assert M["Forêt aléatoire"]["spatial"]["COD"] == pytest.approx(8.96, abs=0.12) and M["XGBoost"]["spatial"]["COD"] == pytest.approx(7.53, abs=0.12)
    # la validation spatiale est toujours plus sévère que la validation aléatoire
    for nom in M: assert M[nom]["spatial"]["rmse_log"] >= M[nom]["aleatoire"]["rmse_log"] - 1e-9, nom
    bloc = res["comparaison"]["plis"]["bloc_spatial"]; assert len(np.unique(bloc)) == 10 and np.bincount(bloc).min() >= 9
    # temporel : 101 ventes de 2026 ; les arbres sous-prédisent (pas d'extrapolation de la tendance), l'hédonique presque pas
    T = {k: v["temporel"] for k, v in M.items()}
    assert T["Hédonique (référence)"]["n"] == 101
    assert T["Hédonique (référence)"]["biais_log"] == pytest.approx(-0.0102, abs=0.001) and T["Gradient boosting"]["biais_log"] == pytest.approx(-0.034, abs=0.004)
    assert T["Forêt aléatoire"]["biais_log"] == pytest.approx(-0.060, abs=0.005) and T["Arbre (profondeur 3)"]["biais_log"] == pytest.approx(-0.085, abs=0.006)
    assert T["Gradient boosting"]["ratio_median"] == pytest.approx(0.966, abs=0.005) and T["Hédonique (référence)"]["ratio_median"] == pytest.approx(0.988, abs=0.003)
    d = res["detendance"]
    assert d["metriques_detendance"]["biais_log"] == pytest.approx(-0.018, abs=0.004) and d["metriques_detendance"]["COD"] == pytest.approx(7.44, abs=0.15)
    assert d["metriques_feuille1"]["biais_log"] == pytest.approx(-0.134, abs=0.01) and d["n_entrainement_apres_mois_17"] == 14
    assert d["beta_mois_entrainement"] == pytest.approx(0.00342, abs=0.00005)
    # le biais du boosting brut s'aggrave au fil des mois de 2026
    b = res["comparaison"]["biais_mensuel_2026"]["Gradient boosting"]; assert b[25.0] < b[19.0] and b[25.0] == pytest.approx(-0.076, abs=0.01)

def test_un_secteur_exclu(res):
    S = res["secteurs"]
    assert S["metriques"]["Hédonique (référence)"]["COD"] == pytest.approx(10.44, abs=0.1) and S["metriques"]["Gradient boosting"]["COD"] == pytest.approx(8.59, abs=0.15)
    assert S["metriques"]["Forêt aléatoire"]["COD"] == pytest.approx(10.21, abs=0.15)
    bh = S["biais_secteur"]["Hédonique (référence)"]
    assert bh["Buckingham"] == pytest.approx(0.234, abs=0.005) and bh["Masson-Angers"] == pytest.approx(0.261, abs=0.005) and abs(bh["Gatineau"]) < 0.02
    # le biais de l'hédonique sur un secteur exclu est (au bruit près) la prime de secteur omise
    assert bh["Buckingham"] == pytest.approx(-c37.VRAIS["sect"]["Buckingham"], abs=0.06) and bh["Masson-Angers"] == pytest.approx(-c37.VRAIS["sect"]["Masson-Angers"], abs=0.06)
    bg = S["biais_secteur"]["Gradient boosting"]
    assert abs(bg["Buckingham"]) < 0.02 and bg["Masson-Angers"] == pytest.approx(0.072, abs=0.01) and bg["Gatineau"] == pytest.approx(-0.068, abs=0.01)

def test_scenario_non_lineaire(res):
    """Perturbation déclarée : le boosting bat l'hédonique mal spécifié, et l'hédonique respécifié reprend l'avantage."""
    sc = res["scenario"]; M = sc["metriques"]
    assert sc["n_prime"] == 76 and sc["n_mur"] == 235
    assert M["Hédonique (référence)"]["COD"] == pytest.approx(8.15, abs=0.1) and M["Gradient boosting"]["COD"] == pytest.approx(7.46, abs=0.12)
    assert M["XGBoost"]["COD"] == pytest.approx(7.40, abs=0.12) and M["Hédonique respécifié"]["COD"] == pytest.approx(6.85, abs=0.1) and M["Forêt aléatoire"]["COD"] == pytest.approx(9.17, abs=0.15)
    assert M["Gradient boosting"]["rmse_log"] < M["Hédonique (référence)"]["rmse_log"] and M["Hédonique respécifié"]["rmse_log"] < M["Gradient boosting"]["rmse_log"]
    assert M["Hédonique (référence)"]["rmse_log"] == pytest.approx(0.1003, abs=0.001) and M["Gradient boosting"]["rmse_log"] == pytest.approx(0.0944, abs=0.0015)
    c = sc["coef_respecifie"]
    assert c["prime"] == pytest.approx(0.129, abs=0.001) and c["mur"] == pytest.approx(-0.092, abs=0.001) and c["sigma"] == pytest.approx(0.0831, abs=0.0005)
    assert abs(c["prime"] - c37.SCENARIO["prime_qualite_age"]) < 1.96 * c["se_prime"] and abs(c["mur"] - c37.SCENARIO["mur"]) < 1.96 * c["se_mur"]
    # l'interaction âge × qualité devient la première interaction SHAP dans le scénario, et n'apparaît pas en tête dans le cas de base
    assert res["shap"]["scenario_interactions"].index[0] == ("age", "qualite") and res["shap"]["scenario_interactions"].iloc[0] == pytest.approx(0.0068, abs=0.001)
    assert res["shap"]["interactions"].index[0] != ("age", "qualite") and res["shap"]["interactions"].iloc[0] < 0.004

def test_shap(res):
    sh = res["shap"]; Xm = res["Xm"]
    assert sh["efficacite_max_err"] < 1e-9 and sh["efficacite_sujet_err"] < 1e-9           # axiome d'efficacité (TreeSHAP exact)
    assert sh["base"] == pytest.approx(13.0732, abs=0.0005) and sh["base"] == pytest.approx(res["y"].mean(), abs=0.001)
    assert math.exp(sh["base"]) == pytest.approx(476_010, rel=0.001)
    imp = sh["importance"]
    assert list(imp.index[:3]) == ["superficie", "age", "mois"] and imp["superficie"] == pytest.approx(0.111, abs=0.002) and imp["age"] == pytest.approx(0.0432, abs=0.002)
    assert imp[[f"sect_{s}" for s in cas1.SECTEURS]].max() < 0.002                             # secteurs absorbés par lat/lon/distance
    s = sh["sujet"]
    assert s["mois"] == pytest.approx(0.0986, abs=0.002) and s["superficie"] == pytest.approx(0.0871, abs=0.002) and s["age"] == pytest.approx(0.0322, abs=0.002)
    assert s["terrain"] == pytest.approx(-0.0127, abs=0.002) and s.sum() + sh["base"] == pytest.approx(res["sujet_pred"].loc["Gradient boosting", "ln_pred"], abs=1e-9)
    # symétrie/efficacité : la somme des |φ| moyens n'excède pas la dispersion de y
    assert imp.sum() < np.std(res["y"]) * 2
    # importance par permutation : même tête de classement
    p = res["permutation"]["Gradient boosting"]
    assert p.idxmax() == "superficie" and p["superficie"] == pytest.approx(0.0348, abs=0.003) and p["age"] == pytest.approx(0.0061, abs=0.0015)
    assert res["permutation"]["Gradient boosting (EQM de base)"] == pytest.approx(0.00867, abs=0.0002)

def test_sujet_et_pdp(res):
    s = res["sujet_pred"]
    assert s.loc["Gradient boosting", "prediction"] == pytest.approx(605_825, rel=0.003) and s.loc["XGBoost", "prediction"] == pytest.approx(613_976, rel=0.003)
    assert s.loc["Forêt aléatoire", "prediction"] == pytest.approx(586_818, rel=0.003) and s.loc["Arbre (profondeur 3)", "prediction"] == pytest.approx(549_573, rel=0.003)
    assert s.loc["Hédonique (référence)", "prediction"] == pytest.approx(638_011, abs=50)
    # tous les modèles d'arbres prédisent le sujet sous l'hédonique (compression vers la moyenne)
    assert (s["prediction"].drop("Hédonique (référence)") < s.loc["Hédonique (référence)", "prediction"]).all()
    assert s.loc["Gradient boosting", "ln_pred"] - s.loc["Gradient boosting", "ln_pred"] == 0 and math.exp(s.loc["Gradient boosting", "ln_pred"]) == pytest.approx(604_595, rel=0.003)
    pdp = res["pdp"]
    # la PDP de la superficie du boosting suit la loi de génération : corrélation > 0,99 et écart maximal < 6 points (aux extrémités du support)
    for nom in ("Gradient boosting", "Hédonique (référence)"):
        a = pdp["superficie"][nom] - pdp["superficie"][nom].mean(); b = pdp["superficie"]["vraie"] - pdp["superficie"]["vraie"].mean()
        assert np.corrcoef(a, b)[0, 1] > 0.99 and np.abs(a - b).max() < 0.06, nom
        assert np.abs((a - b)[(pdp["superficie"]["grille"] >= 100) & (pdp["superficie"]["grille"] <= 220)]).max() < 0.04, nom
    # scénario : la PDP de l'âge du boosting reproduit le mur (creux entre 25 et 45 ans) que le polynôme hédonique ne peut pas suivre
    g = pdp["age_scenario"]["grille"]; gb = pdp["age_scenario"]["Gradient boosting"]; hd = pdp["age_scenario"]["Hédonique (référence)"]; vr = pdp["age_scenario"]["vraie"]
    m = (g >= 25) & (g <= 45)
    assert np.abs((gb - gb.mean()) - (vr - vr.mean())).max() < np.abs((hd - hd.mean()) - (vr - vr.mean())).max()
    assert np.mean(gb[m]) < np.mean(gb[(g > 46) & (g < 60)])

def test_courbes_et_sensibilite(res):
    c = res["courbes"]
    assert c[4]["optimum"] < c[2]["optimum"] <= c[1]["optimum"]                   # plus l'arbre est profond, plus tôt il surapprend
    assert c[4]["entrainement"][-1] < 0.02 and c[4]["validation"][-1] > c[4]["rmse_min"]  # erreur d'entraînement quasi nulle, validation remonte
    assert c[1]["rmse_min"] == pytest.approx(0.0887, abs=0.001) and c[2]["rmse_min"] == pytest.approx(0.0923, abs=0.001) and c[4]["rmse_min"] == pytest.approx(0.0966, abs=0.0015)
    assert c[2]["optimum"] == pytest.approx(730, abs=60) and c[4]["optimum"] == pytest.approx(367, abs=60)
    s = res["sensibilite"].set_index(["modele", "profondeur", "arbres", "nu", "feuille_min", "sous_ech"])
    g = s.loc["Gradient boosting"]
    assert g.loc[(1, 1200, 0.03, 5, 0.8), "rmse_log"] == pytest.approx(0.0887, abs=0.001)      # souches additives : le meilleur boosting
    assert g.loc[(2, 1200, 0.03, 5, 0.8), "rmse_log"] == pytest.approx(0.0931, abs=0.001)
    assert g.loc[(4, 1200, 0.03, 5, 0.8), "rmse_log"] == pytest.approx(0.0984, abs=0.0015) and g.loc[(2, 1200, 0.10, 5, 0.8), "rmse_log"] == pytest.approx(0.0988, abs=0.0015)
    assert g.loc[(2, 300, 0.03, 5, 0.8), "rmse_log"] == pytest.approx(0.0977, abs=0.0015) and g.loc[(2, 3000, 0.03, 5, 0.8), "rmse_log"] == pytest.approx(0.0962, abs=0.0015)
    f = s.loc["Forêt aléatoire"]
    assert f["rmse_log"].between(0.105, 0.117).all() and f.loc[("max", 100, slice(None), 2, "bootstrap"), "rmse_log"].iloc[0] == pytest.approx(0.1065, abs=0.001)

def test_solutions_exercices(res):
    e = res["exercices"]
    # 37.1 : meilleure division des huit ventes
    b = e["1"]["best"]; assert b["seuil"] == 113.5 and b["k"] == 2 and b["gain"] == pytest.approx(0.338, abs=0.001) and e["1"]["sse0"] == pytest.approx(0.462, abs=0.001)
    assert b["c_gauche"] == pytest.approx(12.611, abs=0.001) and b["c_droite"] == pytest.approx(13.0857, abs=0.001)
    assert e["1"]["S"] == [70.0, 99.0, 128.0, 156.0, 160.0, 170.0, 221.0, 241.0]
    # 37.2 : variance d'une moyenne d'arbres corrélés
    assert e["2"]["B1"] == pytest.approx(0.0144) and e["2"]["B100"] == pytest.approx(0.004421, abs=1e-6) and e["2"]["limite"] == pytest.approx(0.00432) and e["2"]["ecart_type_limite"] == pytest.approx(0.0657, abs=0.0001)
    # 37.3 : boosting à la main
    assert e["3"]["F"][-1] == pytest.approx(12.8542, abs=1e-4) and e["3"]["residus"][-1] == pytest.approx(0.1458, abs=1e-4) and e["3"]["m_1pct"] == pytest.approx(43.7, abs=0.05)
    # 37.4 : Shapley à trois joueurs
    phi = e["4"]["phi"]; assert phi["S"] == pytest.approx(70.833, abs=0.001) and phi["Q"] == pytest.approx(35.833, abs=0.001) and phi["L"] == pytest.approx(23.333, abs=0.001) and e["4"]["somme"] == pytest.approx(130.0)
    # 37.5 : dix prédictions d'un AVM
    r = e["5"]; assert r["rmse_log"] == pytest.approx(0.0503, abs=0.0001) and r["mdape"] == pytest.approx(4.51, abs=0.01) and r["COD"] == pytest.approx(4.42, abs=0.01)
    assert r["PRD"] == pytest.approx(1.006, abs=0.001) and r["PRB"] == pytest.approx(-0.034, abs=0.001) and r["ratio_median"] == pytest.approx(0.971, abs=0.001) and r["mape"] == pytest.approx(4.82, abs=0.01)
    # 37.6 et 37.7 : chiffres cohérents avec les tableaux
    assert e["6"]["Hédonique (référence)"]["secteur"] == pytest.approx(10.44, abs=0.1) and e["7"]["Hédonique respécifié"] == pytest.approx(6.85, abs=0.1)