TeX 52.9%
Python 45.7%
Shell 1.4%
1# -*- coding: utf-8 -*-2"""03 — Modèles d'ajustement et modèles de référence, ré-estimés chaque trimestre avec les seules ventes3antérieures au trimestre (as-of 2024-01, 2024-04, …, 2026-07) — comme le ferait un AVM ré-entraîné régulièrement :4 • grille hédonique globale (coefficients) et régionale rétrécie → ajustements B / C5 • LightGBM (objectif L1 = médiane) sur le prix ajusté dans le temps (avec et sans rôle) → ajustements C-GBM, M116 • hédonique globale directe avec effets de municipalité → M107Sorties : models_<période>.json, gbm_pred_<période>.npy (prédiction f pour TOUTES les ventes)."""8import pickle9from common import *10import lightgbm as lgb1112df = load_parquet("sales.parquet")13meta = load_json("index_meta.json"); regions = meta["regions"]; M_EVAL0 = meta["M_EVAL0"]14rid = {r: i for i, r in enumerate(regions)}15df["rid"] = df["region"].map(rid).astype(int)16idx = np.load(OUT / "index_asof.npy")1718FEATS = ["ln_floor", "ln_land", "age", "age2", "lien1", "lien2", "genre1", "genre2", "genre3", "etages", "neuf"]19def feats(d):20 return pd.DataFrame({21 "ln_floor": d["ln_floor"], "ln_land": d["ln_land"], "age": d["age"], "age2": d["age"] ** 2 / 100,22 "lien1": (d["lien"] == 1).astype(float), "lien2": (d["lien"] >= 2).astype(float),23 "genre1": (d["genre"] == 1).astype(float), "genre2": (d["genre"] == 2).astype(float), "genre3": (d["genre"] == 3).astype(float),24 "etages": d["etages"].astype(float), "neuf": d["neuf"].astype(float)}, index=d.index)2526def time_adjusted(d, m_asof):27 """ln P ramené au niveau du mois m_asof-1 avec l'indice connu au mois m_asof."""28 lvl = idx[m_asof - M_EVAL0]29 return d["ln_price"].to_numpy() - lvl[d["rid"].to_numpy(), d["month"].to_numpy()] + lvl[d["rid"].to_numpy(), m_asof - 1]3031def ols(X, y):32 A = np.column_stack([np.ones(len(X)), X])33 b, *_ = np.linalg.lstsq(A, y, rcond=None)34 return b3536ASOF = list(range(M_EVAL0, int(df["month"].max()) + 1, 3)) # ré-estimation trimestrielle (2024-01, 2024-04, …)37save_json({"asof": ASOF}, "asof.json")38TODO = [int(a) for a in os.environ["CV_ASOF"].split(",")] if os.environ.get("CV_ASOF") else ASOF # sous-ensemble (répartition sur plusieurs nœuds)39for m_asof in TODO:40 period = str(m_asof); start = f"{2021 + m_asof // 12}-{m_asof % 12 + 1:02d}-01"41 tr = df[df["month"] < m_asof].copy()42 y = time_adjusted(tr, m_asof)43 F = feats(tr)44 log(f"[{period}] entraînement sur {len(tr):,} ventes < {start} (indice connu au mois {m_asof})")45 # --- grille hédonique globale : estimateur « within » (effets fixes de municipalité), pour que les β46 # mesurent des différences entre biens voisins et non entre marchés47 def within(frame, cols, by):48 return frame[cols] - frame.groupby(by)[cols].transform("mean")49 # effets mensuels résiduels (base = dernier mois connu m_asof-1) : l'indice régional ne retire pas toute la50 # tendance commune ; sans eux, les modèles globaux sous-estiment les ventes récentes51 months_tr = sorted(tr["month"].unique()); MCOLS = [f"m{m}" for m in months_tr if m != m_asof - 1]52 Fm = pd.concat([F] + [pd.Series((tr["month"] == m).astype(float), index=tr.index, name=f"m{m}") for m in months_tr if m != m_asof - 1], axis=1)53 Fw = within(pd.concat([Fm, tr[["code_mun"]]], axis=1), FEATS + MCOLS, "code_mun")54 yw = y - pd.Series(y, index=tr.index).groupby(tr["code_mun"]).transform("mean").to_numpy()55 bw_all, *_ = np.linalg.lstsq(Fw.to_numpy(), yw, rcond=None)56 bw = bw_all[:len(FEATS)]; bm = bw_all[len(FEATS):]57 beta_global = dict(zip(FEATS, bw.tolist()))58 xb = F.to_numpy() @ bw + Fm[MCOLS].to_numpy() @ bm59 resid = y - xb60 r2w = 1 - np.var(yw - Fw.to_numpy() @ bw_all) / np.var(yw)61 Fw = Fw[FEATS]62 log(f" hédonique within-municipalité : R² within = {r2w:.3f} ; β = { {k: round(v, 3) for k, v in beta_global.items()} }")63 # --- grille régionale (within-municipalité) rétrécie vers la globale (n/(n+300))64 beta_region = {}65 for r in regions:66 sel = (tr["region"] == r).to_numpy()67 if sel.sum() < 100: beta_region[r] = beta_global; continue68 b, *_ = np.linalg.lstsq(Fw.to_numpy()[sel], yw[sel], rcond=None)69 w = sel.sum() / (sel.sum() + 300.0)70 beta_region[r] = {k: float(w * b[j] + (1 - w) * beta_global[k]) for j, k in enumerate(FEATS)}71 # --- M10 : hédonique directe = Xβ + effet de région + MRC + municipalité + unité de voisinage (rétrécis en cascade)72 reg_eff = pd.Series(resid, index=tr.index).groupby(tr["region"]).mean().to_dict()73 tr["res_g"] = resid - tr["region"].map(reg_eff).to_numpy()74 mrc_eff = tr.groupby("mrc")["res_g"].agg(["sum", "count"])75 mrc_eff = (mrc_eff["sum"] / (mrc_eff["count"] + 30)).to_dict()76 tr["res_m"] = tr["res_g"] - tr["mrc"].map(mrc_eff).fillna(0)77 mun_eff = tr.groupby("code_mun")["res_m"].agg(["sum", "count"])78 mun_eff = (mun_eff["sum"] / (mun_eff["count"] + 20)).to_dict()79 tr["res_u"] = tr["res_m"] - tr["code_mun"].map(mun_eff).fillna(0)80 # correction médiane : la MCO vise la moyenne du log ; la distribution est asymétrique (queue basse), on recentre sur la médiane81 uv_tmp = tr.groupby(["code_mun", "unite_voisinage"])["res_u"].transform(lambda v: v.sum() / (len(v) + 10))82 median_shift = float(np.nanmedian(tr["res_u"] - uv_tmp.fillna(0)))83 uv = tr.groupby(["code_mun", "unite_voisinage"])["res_u"].agg(["sum", "count"])84 uv_eff = {f"{a}|{b}": v for (a, b), v in (uv["sum"] / (uv["count"] + 10)).items()}85 # --- LightGBM (prix ajusté dans le temps) avec rôle et sans rôle86 def gbm_frame(d):87 G = feats(d)88 G["lat"] = d["lat"]; G["lng"] = d["lng"]; G["x_km"] = d["x_km"]; G["y_km"] = d["y_km"]89 G["rid"] = d["rid"]; G["btype"] = d["btype"]90 G["ln_role"] = np.log(d["role_vente"]); G["part_terrain"] = d["part_terrain"]; G["role_vintage"] = d["role_vintage"]91 G["month"] = np.minimum(d["month"], m_asof - 1).astype(float) # tendance résiduelle ; hors échantillon = dernier mois connu92 return G93 G = gbm_frame(tr)94 params = dict(objective="regression_l1", learning_rate=0.05, num_leaves=127, min_data_in_leaf=40, feature_fraction=0.8,95 bagging_fraction=0.8, bagging_freq=1, lambda_l2=1.0, verbose=-1, num_threads=int(os.environ.get("CV_THREADS", 20)))96 models = {}97 for name, cols in [("gbm", [c for c in G.columns if c not in ("ln_role", "part_terrain", "role_vintage")]),98 ("gbm_role", list(G.columns))]:99 ds = lgb.Dataset(G[cols], y, categorical_feature=["rid", "btype"] if "btype" in cols else ["rid"])100 # validation interne : derniers 10 % chronologiques pour fixer le nombre d'arbres101 n = len(G); cut = int(n * 0.9)102 dtr = lgb.Dataset(G[cols].iloc[:cut], y[:cut], categorical_feature=["rid", "btype"])103 dva = lgb.Dataset(G[cols].iloc[cut:], y[cut:], reference=dtr)104 bst = lgb.train(params, dtr, num_boost_round=3000, valid_sets=[dva], callbacks=[lgb.early_stopping(100, verbose=False)])105 best = max(bst.best_iteration, 300) # plancher : l'arrêt précoce sur la tranche finale (mois partiels) peut être trop court106 bst = lgb.train(params, ds, num_boost_round=best)107 models[name] = bst108 Gall = gbm_frame(df); Gall["month"] = float(m_asof - 1) # f évalué au niveau du dernier mois connu pour toutes les ventes109 pred_all = bst.predict(Gall[cols])110 np.save(OUT / f"{name}_pred_{period}.npy", pred_all.astype(np.float32))111 rr = y - bst.predict(G[cols]); r2 = 1 - np.var(rr) / np.var(y)112 log(f" {name}: {best} arbres, R² (in-sample) = {r2:.3f}, résidu médian = {np.median(rr) * 100:.2f} %")113 bst.save_model(str(OUT / f"{name}_{period}.txt"))114 save_json({"period": period, "m_asof": m_asof, "start": start, "n_train": len(tr), "feats": FEATS,115 "beta_global": beta_global, "intercept_global": 0.0, "region_fe": reg_eff,116 "beta_region": beta_region, "mrc_eff": mrc_eff, "mun_eff": mun_eff, "uv_eff": uv_eff,117 "r2_within": float(r2w), "sd_resid_u": float(tr["res_u"].std()), "median_shift": median_shift},118 f"models_{period}.json")119log("modèles terminés")120