spb/traite-evaluation-qc

1#!/usr/bin/env python32"""Chapitre 20 — figures : densités bootstrap et énumération exhaustive ; résidus des candidats et élargissement ;3fourchettes par méthode ; sensibilité de la conclusion aux poids et variance de la combinaison selon la corrélation.4Usage (racine du dépôt) : PYTHONPATH=analysis/python python3 analysis/python/figures/ch20_figures.py"""5from __future__ import annotations6import sys, pathlib7import numpy as np8import matplotlib.pyplot as plt9from matplotlib.ticker import FuncFormatter10from scipy.stats import gaussian_kde11sys.path.insert(0, str(pathlib.Path(__file__).resolve().parents[1]))12sys.path.insert(0, str(pathlib.Path(__file__).resolve().parents[1] / "cases"))13from traite.style import setup, palette, savefig, LARGEUR, fmt_dollars14import ch20_reconciliation as c201516def fmt_k(x, pos=None):17    return f"{x/1000:,.0f}".replace(",", " ") + " k$"181920def fig_bootstrap(r):21    b, bs, ex = r["base"], r["bootstrap"], r["exhaustif"]22    fig, axes = plt.subplots(1, 2, figsize=(LARGEUR, 3.3))23    ax = axes[0]24    grille = np.linspace(520000, 700000, 400)25    for i, (k, lab) in enumerate([("VA", "candidats"), ("VB", "ajustements"), ("VC", "les deux")]):26        d = gaussian_kde(bs[k])(grille)27        ax.plot(grille, d / d.max(), color=palette.cat[i], lw=1.6 if k != "VB" else 1.2, ls="-" if k != "VB" else "--", label=lab)28    ax.axvline(b["V_comp"], color=palette.ink2, lw=0.8); ax.text(b["V_comp"] - 2000, 1.04, "comparaison\n(5 retenus)", ha="right", va="bottom", fontsize=7, color=palette.ink2)29    ax.axvline(b["V_hed"], color=palette.ink2, lw=0.8, ls=":"); ax.text(b["V_hed"] + 2000, 0.04, "hédonique", ha="left", va="bottom", fontsize=7, color=palette.ink2)30    ax.set_ylim(0, 1.22); ax.set_yticks([]); ax.set_xlabel("Valeur pondérée par comparaison"); ax.set_ylabel("Densité (échelle relative)")31    ax.xaxis.set_major_formatter(FuncFormatter(fmt_k)); ax.legend(loc="upper right", fontsize=7, frameon=False, title="rééchantillonnage", title_fontsize=7)32    ax.set_title("(a) Bootstrap, B = 2 000", loc="left", fontsize=8.5)33    ax = axes[1]34    ax.hist(ex["valeurs"], bins=36, color=palette.cat[0], alpha=0.85, edgecolor="white", lw=0.3)35    ax.axvline(b["V_comp"], color=palette.cat[1], lw=1.4); ax.text(b["V_comp"] - 1500, ax.get_ylim()[1] * 0.95, "sélection\ndu ch. 16", ha="right", va="top", fontsize=7, color=palette.cat[1])36    ax.axvline(b["V_hed"], color=palette.ink2, lw=0.8, ls=":"); ax.text(b["V_hed"] + 1500, ax.get_ylim()[1] * 0.95, "hédonique", ha="left", va="top", fontsize=7, color=palette.ink2)37    ax.axvline(ex["moyenne"], color=palette.ink, lw=0.8, ls="--"); ax.text(ex["moyenne"] + 1500, ax.get_ylim()[1] * 0.55, "moyenne des\n4 368 sous-ensembles", ha="left", va="top", fontsize=7)38    ax.set_xlabel("Valeur pondérée du sous-ensemble de 5 candidats"); ax.set_ylabel("Nombre de sous-ensembles")39    ax.xaxis.set_major_formatter(FuncFormatter(fmt_k)); ax.set_title("(b) Tous les sous-ensembles de 5 parmi 16", loc="left", fontsize=8.5)40    fig.tight_layout(w_pad=2.0)41    savefig(fig, "ch20/ch20-fig-bootstrap-densite.pdf")424344def fig_residus_elargissement(r):45    b, el = r["base"], r["elargissement"]46    a = b["adj16"]47    fig, axes = plt.subplots(1, 2, figsize=(LARGEUR, 3.2))48    ax = axes[0]49    rang = np.arange(1, len(a) + 1); res = (np.exp(a["residu_log"].values) - 1) * 10050    ret = a["retenu"].values51    ax.axhline(0, color=palette.ink2, lw=0.8)52    ax.bar(rang[ret], res[ret], color=palette.cat[1], width=0.7, label="cinq comparables retenus")53    ax.bar(rang[~ret], res[~ret], color=palette.cat[0], width=0.7, alpha=0.75, label="onze autres candidats")54    ax.axhline(res[ret].mean(), color=palette.cat[1], lw=1.0, ls="--"); ax.text(16.4, res[ret].mean(), f"moy. {res[ret].mean():.1f} %".replace(".", ","), fontsize=7, va="center", ha="left", color=palette.cat[1])55    ax.axhline(res[~ret].mean(), color=palette.cat[0], lw=1.0, ls="--"); ax.text(16.4, res[~ret].mean() + 1.0, f"moy. {res[~ret].mean():.1f} %".replace(".", ","), fontsize=7, va="bottom", ha="left", color=palette.cat[0])56    ax.set_xticks(rang); ax.set_xticklabels(a["id_vente"].values, rotation=90, fontsize=6.5)57    ax.set_xlim(0.3, 18.5); ax.set_ylabel("Résidu du modèle à la vente (%)"); ax.set_xlabel("Candidats, par distance d'attributs croissante")58    ax.set_ylim(-22, 14); ax.legend(loc="upper left", fontsize=7, frameon=False); ax.set_title("(a) Résidus des 16 candidats", loc="left", fontsize=8.5)59    ax = axes[1]60    ax.plot(el["k"], el["valeur_ponderee"], marker="o", color=palette.cat[0], lw=1.6, label="valeur pondérée (inverse du brut)")61    ax.plot(el["k"], el["moyenne"], marker="s", ms=3.5, color=palette.cat[2], lw=1.0, ls="--", label="moyenne simple")62    ax.axhline(b["V_hed"], color=palette.ink2, lw=0.8, ls=":"); ax.text(16.2, b["V_hed"] - 1500, "hédonique", fontsize=7, ha="right", va="top", color=palette.ink2)63    ax.set_xlabel("Nombre k de comparables retenus (par distance croissante)"); ax.set_ylabel("Indication par comparaison")64    ax.yaxis.set_major_formatter(FuncFormatter(fmt_k)); ax.set_xticks(el["k"]); ax.set_ylim(560000, 650000); ax.legend(loc="lower right", fontsize=7, frameon=False)65    ax.set_title("(b) Élargissement de l'échantillon", loc="left", fontsize=8.5)66    fig.tight_layout(w_pad=2.0)67    savefig(fig, "ch20/ch20-fig-residus-elargissement.pdf")686970def fig_fourchettes(r):71    b, bs, co, cb, vt = r["base"], r["bootstrap"], r["cout"], r["combinaison"], r["verite"]72    ip = b["ip_hed"]73    z90 = 1.644974    lignes = [75        ("Prix de vente probable\n(intervalle de prédiction à 95 %, ch. 18)", b["V_hed"], ip["bas"], ip["haut"], palette.cat[3], 0.5),76        ("Coût (hypothèses α, vie économique :\nétendue de la grille)", co["valeur"], co["min_grille"], co["max_grille"], palette.cat[2], 1.0),77        ("Comparaison, 5 comparables\n(bootstrap conjoint, 90 %)", b["V_comp"], bs["resume"]["conjoint"]["q05"], bs["resume"]["conjoint"]["q95"], palette.cat[1], 1.0),78        ("Comparaison élargie, 16 candidats\n(point)", r["elargissement"]["valeur_ponderee"].iloc[-1], np.nan, np.nan, palette.cat[1], 1.0),79        ("Hédonique (estimation + spécification,\n± 1,645 écart-type)", b["V_hed"], b["V_hed"] - z90 * cb["s"][1], b["V_hed"] + z90 * cb["s"][1], palette.cat[0], 1.0),80        ("Conclusion : valeur la plus probable,\nfourchette conclue et bootstrap 90 %", cb["conclusion"], cb["fourchette_conclue"][0], cb["fourchette_conclue"][1], palette.brand, 1.0),81    ]82    fig, ax = plt.subplots(figsize=(LARGEUR, 3.6))83    y = np.arange(len(lignes))[::-1]84    for yi, (lab, pt, lo, hi, col, alpha) in zip(y, lignes):85        if not np.isnan(lo):86            ax.plot([lo, hi], [yi, yi], color=col, lw=2.2 if alpha == 1.0 else 3.5, alpha=alpha, solid_capstyle="butt")87            ax.plot([lo, lo], [yi - 0.12, yi + 0.12], color=col, lw=1.0, alpha=alpha); ax.plot([hi, hi], [yi - 0.12, yi + 0.12], color=col, lw=1.0, alpha=alpha)88        ax.plot(pt, yi, "o", color=col, ms=6, zorder=5, mec="white", mew=0.8)89        ax.text(pt, yi + 0.2, fmt_dollars(pt), ha="center", va="bottom", fontsize=7)90    # bande bootstrap 90 % de la conclusion (plus épaisse, sous la fourchette conclue)91    yc = y[-1]92    ax.plot([cb["fourchette90"][0], cb["fourchette90"][1]], [yc, yc], color=palette.brand, lw=6, alpha=0.35, solid_capstyle="butt", zorder=1)93    ax.axvline(cb["conclusion"], color=palette.brand, lw=0.8, ls="--", zorder=0)94    ax.axvline(vt["V_true"], color=palette.ink2, lw=0.8, ls=":", zorder=0)95    ax.text(vt["V_true"] + 3000, -0.55, "espérance de la loi\nde génération", fontsize=7, color=palette.ink2, ha="left", va="center")96    ax.set_ylim(-0.9, len(lignes) - 0.3); ax.set_yticks(y); ax.set_yticklabels([l[0] for l in lignes], fontsize=7.2)97    ax.set_xlim(270000, 770000); ax.xaxis.set_major_formatter(FuncFormatter(fmt_k)); ax.set_xlabel("Dollars (synthétiques), 1er septembre 2026")98    ax.grid(axis="y", visible=False)99    fig.tight_layout()100    savefig(fig, "ch20/ch20-fig-fourchettes-methodes.pdf")101102103def fig_sensibilite_poids(r):104    cb = r["combinaison"]105    fig, axes = plt.subplots(1, 2, figsize=(LARGEUR, 3.2))106    ax = axes[0]107    for i, (wc, courbe) in enumerate(cb["courbes"].items()):108        lab = "poids du coût = " + (f"{wc*100:.0f} %")109        ax.plot(cb["wh"], courbe, color=palette.cat[i], lw=1.6 if wc == 0.05 else 1.1, ls="-" if wc == 0.05 else "--", label=lab)110    tab = cb["tableau"]111    marqueurs = {"raisonnée (retenue)": ("o", palette.brand), "inverse de la variance (indépendance)": ("s", palette.cat[0]), "variance minimale avec corrélation (MCG)": ("^", palette.cat[0]), "égale (un tiers chacune)": ("D", palette.cat[2])}112    for _, row in tab.iterrows():113        if row["schema"] in marqueurs:114            m, c = marqueurs[row["schema"]]115            ax.plot(row["w_hed"], row["valeur"], m, color=c, ms=6, mec="white", mew=0.6, zorder=5)116    ax.annotate("raisonnée\n(0,40 / 0,55 / 0,05)", (0.55, tab.loc[tab["schema"] == "raisonnée (retenue)", "valeur"].iloc[0]), xytext=(0.62, 572000), fontsize=7, arrowprops=dict(arrowstyle="-", color=palette.ink2, lw=0.6))117    ax.annotate("inverse de la variance\net MCG", (0.70, tab.loc[tab["schema"] == "inverse de la variance (indépendance)", "valeur"].iloc[0]), xytext=(0.42, 628000), fontsize=7, arrowprops=dict(arrowstyle="-", color=palette.ink2, lw=0.6))118    ax.annotate("égale", (1 / 3, tab.loc[tab["schema"] == "égale (un tiers chacune)", "valeur"].iloc[0]), xytext=(0.12, 540000), fontsize=7, arrowprops=dict(arrowstyle="-", color=palette.ink2, lw=0.6))119    ax.set_xlabel("Poids de l'indication hédonique $w_h$ (le reste à la comparaison)"); ax.set_ylabel("Valeur combinée")120    ax.yaxis.set_major_formatter(FuncFormatter(fmt_k)); ax.set_xlim(0, 1); ax.legend(loc="lower right", fontsize=7, frameon=False)121    ax.set_title("(a) Conclusion selon les poids", loc="left", fontsize=8.5)122    ax = axes[1]123    s = cb["s"]124    for i, (rho, var) in enumerate(cb["var_deux"].items()):125        sd = np.sqrt(np.maximum(var, 0)) / 1000126        ax.plot(cb["wh"], sd, color=palette.seq[min(2 + i, len(palette.seq) - 1)], lw=1.4, label=f"ρ = {rho:.1f}".replace(".", ","))127        w1 = c20.poids_optimaux_deux(s[0], s[1], rho)128        if 0 <= 1 - w1 <= 1:129            ax.plot(1 - w1, np.sqrt(w1 ** 2 * s[0] ** 2 + (1 - w1) ** 2 * s[1] ** 2 + 2 * w1 * (1 - w1) * rho * s[0] * s[1]) / 1000, "o", color=palette.seq[min(2 + i, len(palette.seq) - 1)], ms=5, mec="white", mew=0.6)130    ax.set_xlabel("Poids de l'indication hédonique $w_h$ ($w_c = 1 - w_h$)"); ax.set_ylabel("Écart-type de la combinaison (k$)")131    ax.set_xlim(0, 1); ax.legend(loc="upper center", fontsize=7, frameon=False, ncol=2)132    ax.set_title("(b) Variance de deux indications corrélées", loc="left", fontsize=8.5)133    fig.tight_layout(w_pad=2.0)134    savefig(fig, "ch20/ch20-fig-sensibilite-poids.pdf")135136137def main():138    setup()139    r = c20.tout()140    fig_bootstrap(r); fig_residus_elargissement(r); fig_fourchettes(r); fig_sensibilite_poids(r)141    print("figures ch20 écrites dans assets/figures/ch20/")142143144if __name__ == "__main__":145    main()146