spb/traite-evaluation-qc

1"""Chapitre 19 — tests des fonctions spatiales ajoutées à traite.hedonique et des chiffres cités dans le texte et les solutions.2Tolérances : coefficients ±0,001 (ou relatives 1 %), log-vraisemblances ±0,1, indications en dollars ±1 $ (arrondi au dollar)."""3import sys, pathlib4import numpy as np, pytest5sys.path.insert(0, str(pathlib.Path(__file__).resolve().parents[1] / "cases"))6from traite import hedonique as h7import ch19_spatial as c1989# ----------------------------------------------------------------------------- fonctions génériques (données simulées)1011def _simul(n=400, seed=5):12    rng = np.random.default_rng(seed); coords = rng.uniform(0, 10, (n, 2)); W = h.matrice_poids_knn(coords, 6)13    X = np.column_stack([np.ones(n), rng.normal(size=n), rng.normal(size=n)])14    return rng, coords, W, X1516def test_sem_recupere_lambda_et_beta():17    rng, coords, W, X = _simul(seed=21); n = len(X); beta = np.array([1.0, 0.8, -0.5]); lam = 0.618    u = np.linalg.solve(np.eye(n) - lam * W, rng.normal(0, 0.3, n)); y = X @ beta + u19    r = h.sem_mv(y, X, W)20    assert abs(r["lambda"] - lam) < 0.12 and np.allclose(r["beta"], beta, atol=0.1)21    assert r["loglik"] > h.sar_mv(y, X, W)["loglik"] - 50      # cohérence des échelles de vraisemblance2223def test_sdm_recupere_theta_et_sar_coherent_avec_sar_lag_mv():24    rng, coords, W, X = _simul(seed=8); n = len(X); beta = np.array([1.0, 0.8, -0.5]); theta = np.array([0.4, 0.0]); rho = 0.425    y = np.linalg.solve(np.eye(n) - rho * W, X @ beta + (W @ X[:, 1:]) @ theta + rng.normal(0, 0.3, n))26    r = h.sar_mv(y, X, W, durbin=True)27    assert abs(r["rho"] - rho) < 0.12 and np.allclose(r["beta"][:3], beta, atol=0.12) and abs(r["beta"][3] - 0.4) < 0.1528    # le SAR par raffinement coïncide avec la recherche sur grille de sar_lag_mv (à la résolution de la grille)29    s1 = h.sar_mv(y, X, W); s2 = h.sar_lag_mv(y, X, W)30    assert abs(s1["rho"] - s2["rho"]) < 0.013132def test_tests_lm_detectent_la_bonne_forme():33    rng, coords, W, X = _simul(seed=3); n = len(X); beta = np.array([1.0, 0.8, -0.5])34    y0 = X @ beta + rng.normal(0, 0.3, n)35    t0 = h.tests_lm_spatial(y0, X, W)36    assert t0["LM_erreur"][1] > 0.01 and t0["LM_decalage"][1] > 0.01           # pas de dépendance : tests non significatifs37    u = np.linalg.solve(np.eye(n) - 0.6 * W, rng.normal(0, 0.3, n)); y1 = X @ beta + u38    t1 = h.tests_lm_spatial(y1, X, W)39    assert t1["LM_erreur"][1] < 1e-6 and t1["RLM_erreur"][0] > t1["RLM_decalage"][0]40    assert t1["SARMA"][0] == pytest.approx(t1["LM_erreur"][0] + t1["RLM_decalage"][0])4142def test_effets_sar_identites():43    _, coords, W, X = _simul(); rho = 0.5; beta = np.array([1.0, 0.8, -0.5])44    e = h.effets_sar(rho, beta, W)45    assert np.allclose(e["total"], e["direct"] + e["indirect"])46    # W standardisée en ligne : effet total moyen = β/(1 − ρ) exactement47    assert np.allclose(e["total"], beta / (1 - rho), atol=1e-8)48    assert (np.abs(e["direct"] - beta) < np.abs(e["indirect"]) + 1e-9).all() and e["direct"][1] > beta[1]4950def test_gwr_grande_bande_egale_mco_et_cv_minimal_raisonnable():51    rng, coords, W, X = _simul(seed=2); n = len(X)52    y = X @ np.array([1.0, 0.8, -0.5]) + rng.normal(0, 0.3, n)53    g = h.gwr(y, X, coords, bande=1e4); r = h.mco(X, y)54    assert np.allclose(g["beta"], r["beta"][None, :], atol=1e-4) and abs(g["trace_S"] - 3) < 1e-355    # coefficient variant dans l'espace : la bande choisie par CV est finie et la GWR améliore la CV56    y2 = 1.0 + (0.3 + 0.1 * coords[:, 0]) * X[:, 1] - 0.5 * X[:, 2] + rng.normal(0, 0.3, n)57    sel = h.gwr_bande_cv(y2, X, coords, bornes=(0.5, 50.0), n_grille=12)58    assert sel["bande"] < 20 and sel["cv"] < h.gwr(y2, X, coords, 1e4)["cv"]5960def test_variogramme_bruit_blanc_plat_et_krigeage_exact():61    rng = np.random.default_rng(9); coords = rng.uniform(0, 20, (500, 2)); z = rng.normal(0, 0.1, 500)62    vg = h.variogramme_empirique(coords, z, np.arange(0, 10.1, 1.0))63    assert np.allclose(vg["gamma"], 0.01, atol=0.0025)                              # γ(h) ≈ σ² pour un bruit blanc64    p = {"pepite": 0.0, "palier_partiel": 0.01, "portee": 3.0}65    k = h.krigeage_ordinaire(coords[:5], z[:5], coords[:5], p, voisins=None)66    assert np.allclose(k["prediction"], z[:5], atol=1e-8) and np.allclose(k["variance"], 0.0, atol=1e-8)   # interpolateur exact sans pépite67    assert np.allclose(h.variogramme_exponentiel(np.array([0.0, 1e9]), 0.002, 0.008, 2.0), [0.0, 0.010])6869def test_matrice_spatio_temporelle_unidirectionnelle_et_nilpotente():70    rng = np.random.default_rng(4); coords = rng.uniform(0, 5, (120, 2)); t = np.sort(rng.uniform(0, 36, 120))71    W = h.matrice_poids_spatio_temporelle(coords, t, portee=2.0, fenetre=12)72    assert np.all(np.triu(W) == 0)                                                  # seules les ventes antérieures (indices inférieurs) pèsent73    assert np.abs(np.linalg.eigvals(W)).max() < 1e-10                             # nilpotente : log|I − ρW| = 074    s = W.sum(1); assert np.allclose(s[s > 0], 1.0)75    Wb = h.matrice_poids_spatio_temporelle(coords, t, 2.0, 12, unidirectionnel=False)76    assert (Wb > 0).sum() >= (W > 0).sum()7778def test_coordonnees_km_et_moran_par_distance():79    xy = h.coordonnees_km([45.0, 45.0 + 1 / 110.57], [-75.0, -75.0], 45.0, -75.0)80    assert xy[1, 1] == pytest.approx(1.0, abs=1e-6) and xy[1, 0] == pytest.approx(0.0, abs=1e-9)81    rng = np.random.default_rng(6); coords = rng.uniform(0, 10, (300, 2)); x = coords[:, 0] + rng.normal(0, 0.3, 300)82    cg = h.moran_par_distance(x, coords, [0, 1, 2, 4, 8])83    assert cg[0]["I"] > cg[-1]["I"] and cg[0]["z"] > 38485# ----------------------------------------------------------------------------- chiffres du chapitre (cas 1)8687@pytest.fixture(scope="module")88def res():89    return c19.main(verbose=False)9091def test_moran_residus_avec_et_sans_secteur(res):92    m = res["moran"].set_index("serie")93    assert m.loc["Résidus, modèle de référence", "I_knn"] == pytest.approx(-0.0136, abs=0.001)94    assert abs(m.loc["Résidus, modèle de référence", "z_knn"]) < 1.9695    assert m.loc["Résidus, sans indicatrices de secteur", "I_knn"] == pytest.approx(0.3237, abs=0.001)96    assert m.loc["Résidus, sans indicatrices de secteur", "z_knn"] == pytest.approx(17.11, abs=0.05)97    assert m.loc["Résidus, sans secteur ni distance", "I_knn"] == pytest.approx(0.3595, abs=0.001)98    assert m.loc["ln P (sans modèle)", "I_knn"] == pytest.approx(0.2013, abs=0.001)99    assert m.loc["Résidus, sans indicatrices de secteur", "I_st"] == pytest.approx(0.305, abs=0.001)100    assert res["base"]["modele"]["r2"] == pytest.approx(0.853, abs=0.001) and res["base"]["mco_r"]["r2"] == pytest.approx(0.766, abs=0.001)101102def test_correlogramme(res):103    cg = res["correlogramme"]["sans_secteur"]104    assert cg[0]["I"] == pytest.approx(0.225, abs=0.002) and cg[4]["I"] == pytest.approx(0.261, abs=0.002)105    assert cg[-1]["I"] < -0.25 and all(abs(c["z"]) < 2 for c in res["correlogramme"]["reference"])106107def test_tests_lm(res):108    r = res["lm"]["sans_secteur"]109    assert r["LM_erreur"][0] == pytest.approx(285.18, abs=0.05) and r["LM_decalage"][0] == pytest.approx(116.21, abs=0.05)110    assert r["RLM_erreur"][0] == pytest.approx(170.06, abs=0.05) and r["RLM_decalage"][0] == pytest.approx(1.09, abs=0.01)111    assert r["RLM_decalage"][1] == pytest.approx(0.2975, abs=0.001)112    c = res["lm"]["reference"]113    assert all(p > 0.1 for _, p in c.values())114115def test_sar_sem_sdm(res):116    m = res["modeles"]117    assert m["sar_r"]["rho"] == pytest.approx(0.4385, abs=0.0005) and m["sar_r"]["se_rho"] == pytest.approx(0.0412, abs=0.0005)118    assert m["sem_r"]["lambda"] == pytest.approx(0.6486, abs=0.0005) and m["sem_r"]["se_lambda"] == pytest.approx(0.0447, abs=0.0005)119    assert m["sdm_r"]["rho"] == pytest.approx(0.5089, abs=0.0005)120    assert m["mco_r"]["loglik"] == pytest.approx(510.71, abs=0.05) and m["sar_r"]["loglik"] == pytest.approx(554.23, abs=0.05)121    assert m["sem_r"]["loglik"] == pytest.approx(590.58, abs=0.05) and m["sdm_r"]["loglik"] == pytest.approx(606.72, abs=0.05)122    assert m["mco_c"]["loglik"] == pytest.approx(650.00, abs=0.05)123    assert m["mco_r"]["aic"] == pytest.approx(-995.4, abs=0.1) and m["sar_r"]["aic"] == pytest.approx(-1080.5, abs=0.1)124    assert m["sem_r"]["aic"] == pytest.approx(-1153.2, abs=0.1) and m["sdm_r"]["aic"] == pytest.approx(-1163.4, abs=0.1) and m["mco_c"]["aic"] == pytest.approx(-1266.0, abs=0.1)125    assert m["sar_r"]["moran_residus"] == pytest.approx(0.121, abs=0.002) and m["sem_r"]["moran_residus"] == pytest.approx(-0.041, abs=0.002)126    # modèle de référence : aucun paramètre spatial significatif127    assert abs(m["sar_c"]["rho"] / m["sar_c"]["se_rho"]) < 1.96 and abs(m["sem_c"]["lambda"] / m["sem_c"]["se_lambda"]) < 1.96128    assert m["sar_c"]["rho"] == pytest.approx(-0.080, abs=0.001) and m["sem_c"]["lambda"] == pytest.approx(-0.066, abs=0.001)129    # coefficients cités130    cols = list(res["base"]["Xr"].columns)131    assert m["sar_r"]["beta"][cols.index("ln_S")] == pytest.approx(0.5571, abs=0.0005)132    assert m["sem_r"]["beta"][cols.index("ln_S")] == pytest.approx(0.5416, abs=0.0005)133    assert m["mco_r"]["beta"][cols.index("ln_T")] == pytest.approx(0.0419, abs=0.0005) and m["sem_r"]["beta"][cols.index("ln_T")] == pytest.approx(0.0886, abs=0.0005)134    assert m["mco_r"]["beta"][cols.index("dist_centre")] == pytest.approx(-0.0117, abs=0.0005)135    assert m["sdm_r"]["beta"][len(cols) + cols.index("ln_S") - 1] == pytest.approx(-0.2439, abs=0.0005)   # θ de W·ln_S136137def test_effets_directs_indirects(res):138    e = res["effets"]139    assert e.loc["ln_S", "sar_direct"] == pytest.approx(0.5715, abs=0.0005) and e.loc["ln_S", "sar_indirect"] == pytest.approx(0.4206, abs=0.0005)140    assert e.loc["ln_S", "sar_total"] == pytest.approx(0.9921, abs=0.0005)141    assert e.loc["ln_S", "sar_total"] == pytest.approx(res["modeles"]["sar_r"]["beta"][1] / (1 - res["modeles"]["sar_r"]["rho"]), abs=1e-6)142    assert e.loc["qualite", "sar_total"] == pytest.approx(0.1076, abs=0.0005)143    assert e.loc["ln_S", "sdm_direct"] == pytest.approx(0.5484, abs=0.0005) and e.loc["ln_S", "sdm_indirect"] == pytest.approx(0.0666, abs=0.0005)144    assert e.loc["ln_S", "sdm_total"] == pytest.approx(0.6150, abs=0.0005)145146def test_gwr(res):147    g = res["gwr"]148    assert g["selection"]["bande"] == pytest.approx(4.99, abs=0.03)149    assert g["gwr"]["trace_S"] == pytest.approx(44.8, abs=0.1) and g["gwr"]["aicc"] == pytest.approx(-1233.2, abs=0.2)150    assert g["aicc_mco_sans_secteur"] == pytest.approx(-994.8, abs=0.2) and g["aicc_mco_reference"] == pytest.approx(-1264.9, abs=0.2)151    assert g["gwr"]["r2"] == pytest.approx(0.861, abs=0.001) and g["gwr"]["sigma"] == pytest.approx(0.0829, abs=0.0002)152    q = g["quartiles"]153    assert q.loc["ln_S", "50%"] == pytest.approx(0.5418, abs=0.0005) and q.loc["ln_S", "min"] == pytest.approx(0.5083, abs=0.0005) and q.loc["ln_S", "max"] == pytest.approx(0.5522, abs=0.0005)154    c = g["comparaison_secteurs"]155    assert c.loc["Buckingham", "gwr_localisation_rel_hull"] == pytest.approx(-0.3131, abs=0.001) and c.loc["Buckingham", "reference_localisation"] == pytest.approx(-0.3440, abs=0.001)156    assert c.loc["Buckingham", "generateur_localisation"] == pytest.approx(-0.3418, abs=0.001)157    assert c.loc["Aylmer", "gwr_localisation_rel_hull"] == pytest.approx(0.0034, abs=0.001) and c.loc["Aylmer", "generateur_localisation"] == pytest.approx(0.0622, abs=0.001)158    assert c.loc["Masson-Angers", "gwr_localisation_rel_hull"] == pytest.approx(-0.3190, abs=0.001)159160def test_effets_fixes(res):161    ef = res["effets_fixes"].set_index("specification")162    assert ef.iloc[1]["aic"] == pytest.approx(-1266.0, abs=0.1)                     # référence : meilleur AIC163    assert ef["aic"].idxmin() == ef.index[1]164    g2 = [i for i in ef.index if i.startswith("Grille de 2")][0]; g5 = [i for i in ef.index if i.startswith("Grille de 5")][0]165    assert ef.loc[g2, "r2"] == pytest.approx(0.8721, abs=0.0005) and ef.loc[g2, "aic"] == pytest.approx(-1225.4, abs=0.1) and ef.loc[g2, "parametres"] == 78166    assert ef.loc[g5, "aic"] == pytest.approx(-1252.1, abs=0.1) and ef.loc[g5, "parametres"] == 30167    tri = [i for i in ef.index if "trimestrielles" in i][0]168    assert ef.loc[tri, "aic"] == pytest.approx(-1253.8, abs=0.1)169    idx = res["effets_fixes"].attrs["indice_trimestriel"]170    assert idx["2026Q3"] == pytest.approx(0.199, abs=0.001) and res["effets_fixes"].attrs["tendance_mensuelle"] == pytest.approx(0.0036, abs=0.0001)171172def test_variogramme_krigeage(res):173    vk = res["variogramme"]; p = vk["parametres"]174    assert vk["variogramme"]["gamma"][0] == pytest.approx(0.00699, abs=0.00005) and vk["variogramme"]["gamma"].max() == pytest.approx(0.01702, abs=0.00005)175    assert p["pepite"] == pytest.approx(0.00455, abs=0.00005) and p["palier_partiel"] == pytest.approx(0.00897, abs=0.00005) and p["portee"] == pytest.approx(5.0, abs=1e-6)176    assert vk["variance_residus"] == pytest.approx(0.01067, abs=0.00005) and vk["sigma2_reference"] == pytest.approx(0.00689, abs=0.00005)177    assert vk["sujet"]["prediction"][0] == pytest.approx(-0.0597, abs=0.0005) and vk["sujet"]["variance"][0] == pytest.approx(0.00567, abs=0.00005)178    assert np.all(np.abs(vk["variogramme_reference"]["gamma"] - 0.0068) < 0.0006)   # variogramme plat pour le modèle de référence179180def test_spatio_temporel(res):181    st = res["spatio_temporel"]182    assert st["n_retenues"] == 582 and st["n_sans_voisin"] == 18 and st["max_abs_valeur_propre"] < 1e-10183    assert st["moran_e_st"]["I"] == pytest.approx(0.305, abs=0.001)184    assert st["rho_st"] == pytest.approx(0.3696, abs=0.0005) and st["se_rho_st"] == pytest.approx(0.0532, abs=0.0005)185    assert st["moran_e_apres_lag"] == pytest.approx(0.128, abs=0.002)186    assert st["sem_st"]["lambda"] == pytest.approx(0.7485, abs=0.0005) and st["sem_st"]["aic"] == pytest.approx(-1154.8, abs=0.1)187188def test_validation_croisee(res):189    v = res["validation"].set_index("modele")["rmse_log"]190    attendu = {"MCO, modèle de référence (secteurs)": 0.0841, "MCO, sans indicatrices de secteur": 0.1053, "MCO sans secteur + krigeage des résidus": 0.0878,191               "SEM sans secteur (espérance conditionnelle)": 0.0882, "SAR sans secteur (prédiction conditionnelle)": 0.0950, "GWR sans secteur (validation croisée)": 0.0860}192    for k, x in attendu.items():193        assert v[k] == pytest.approx(x, abs=0.0001), k194    e = res["validation"].set_index("modele")["erreur_abs_moy_pct"]195    assert e["MCO, modèle de référence (secteurs)"] == pytest.approx(0.0681, abs=0.0002) and e["MCO, sans indicatrices de secteur"] == pytest.approx(0.0838, abs=0.0002)196197def test_indications_sujet(res):198    s = res["sujet"]["indication"].round(0).values199    assert list(s) == [637953, 664342, 624799, 633352, 654429, 627519]200    assert res["sujet"]["ln_correction"].iloc[3] == pytest.approx(-0.0520, abs=0.0005) and res["sujet"]["ln_correction"].iloc[5] == pytest.approx(-0.0311, abs=0.0005)201    assert all(sct == "Hull" for sct in res["sujet"].attrs["secteurs_voisins"]) and len(res["sujet"].attrs["voisins_sujet"]) == 8202    assert "G0404" in res["sujet"].attrs["voisins_sujet"] and "G0574" in res["sujet"].attrs["voisins_sujet"]203204def test_sensibilite(res):205    s = res["sensibilite"]; k = s[s["famille"].str.startswith("k")].set_index("parametre")206    assert list(k.index) == [4, 6, 8, 12, 16] and (k["lignes_vides"] == 0).all()207    assert k.loc[4, "rho"] == pytest.approx(0.2359, abs=0.0005) and k.loc[16, "rho"] == pytest.approx(0.6537, abs=0.0005)208    assert k.loc[4, "lambda"] == pytest.approx(0.5066, abs=0.0005) and k.loc[16, "lambda"] == pytest.approx(0.7855, abs=0.0005)209    assert k["rho"].is_monotonic_increasing and k["lambda"].is_monotonic_increasing210    assert (k["I"].between(0.31, 0.34)).all()211    d = s[s["famille"] == "distance inverse"].set_index("parametre")212    assert d.loc[2.0, "rho"] == pytest.approx(0.6360, abs=0.0005) and d.loc[4.0, "lambda"] == pytest.approx(0.8876, abs=0.0005)213214# ----------------------------------------------------------------------------- solutions des exercices215216def test_solutions_exercices():217    from scipy.stats import chi2218    # ex. 19.1 : cinq maisons alignées, contiguïté standardisée219    e = np.array([0.10, 0.06, -0.02, -0.08, -0.06]); W = np.zeros((5, 5))220    for i in range(5):221        for j in (i - 1, i + 1):222            if 0 <= j < 5: W[i, j] = 1223    W = W / W.sum(1, keepdims=True); r = h.moran_i(e, W)224    assert r["I"] == pytest.approx(0.6917, abs=0.0005) and r["E_I"] == pytest.approx(-0.25) and r["z"] == pytest.approx(2.37, abs=0.01)225    # ex. 19.3 : règle de décision226    assert chi2.sf(1.09, 1) == pytest.approx(0.296, abs=0.001) and chi2.ppf(0.95, 1) == pytest.approx(3.841, abs=0.001)227    # ex. 19.4 : effet total = β/(1 − ρ)228    assert 0.5571 / (1 - 0.4385) == pytest.approx(0.992, abs=0.001) and 0.0604 / (1 - 0.4385) == pytest.approx(0.1076, abs=0.0005)229    # ex. 19.5 : poids gaussiens230    d = np.array([0.5, 1, 2, 5])231    assert np.allclose(np.exp(-0.5 * (d / 2) ** 2), [0.969, 0.882, 0.607, 0.044], atol=0.001)232    assert np.allclose(np.exp(-0.5 * (d / 5) ** 2), [0.995, 0.980, 0.923, 0.607], atol=0.001)233    # ex. 19.6 : krigeage ordinaire à trois points234    co = np.array([[1, 0], [0, 2], [-3, 0.0]]); z = np.array([0.08, 0.05, -0.01]); p = {"pepite": 0.0045, "palier_partiel": 0.009, "portee": 5.0}235    k = h.krigeage_ordinaire(co, z, np.array([[0, 0.0]]), p, None)236    assert k["prediction"][0] == pytest.approx(0.0489, abs=0.0001) and k["variance"][0] == pytest.approx(0.00869, abs=0.00001)237    # ex. 19.1 (suite) : signes alternés → I = −1 exactement sur cette configuration238    e2 = np.array([0.10, -0.10, 0.10, -0.10, 0.10]); assert h.moran_i(e2, W)["I"] == pytest.approx(-1.0, abs=1e-9)239    # ex. 19.6 : poids du krigeage vs inverse de la distance240    idw = (1 / np.array([1, 2, 3.0])); idw /= idw.sum(); assert idw[2] == pytest.approx(0.182, abs=0.001) and (idw * z).sum() == pytest.approx(0.0555, abs=0.0001)241242def test_solution_exercice_7(res):243    r = c19.exercice_7(res["base"]); s = r["sem_par_k"].set_index("k")244    assert s.loc[4, "lambda"] == pytest.approx(0.5066, abs=0.0005) and s.loc[16, "lambda"] == pytest.approx(0.7855, abs=0.0005)245    assert s.loc[4, "u_sujet"] == pytest.approx(-0.0189, abs=0.0005) and s.loc[16, "u_sujet"] == pytest.approx(-0.0554, abs=0.0005)246    assert list(s["indication"].round(0)) == [656724, 633352, 628672]247    assert r["st_24"]["sans_voisin"] == 16 and r["st_24"]["rho_st"] == pytest.approx(0.4149, abs=0.0005) and r["st_12"]["rho_st"] == pytest.approx(0.3696, abs=0.0005)248