SPB Git

Ch. 13 livré et validé (compilation isolée, tests, fusion des registres)

Co-Authored-By: Claude Fable 5.1 <noreply@anthropic.com>
Simon-Pierre Boucher committed (Oct 31, 2026) parent 5e249f7 Browse files

65 changed files +6,394 −93

65 changed files
modified analysis/python/cases/ch14_nettoyage.py +18 −1
@@ -476,6 +476,22 @@ def comparer_imputations(P: dict) -> pd.DataFrame:
476 476 "age": c.loc["age", "coef"], "prediction_sujet": pred, **{f"coef_{k}": v for k, v in ind.items()}})
477 477 return pd.DataFrame(rows).set_index("methode")
478 478
479 +def variante_mnar_terrain(seed: int = GRAINE) -> pd.DataFrame:
480 + """Exercice 14.7 : les 20 terrains manquants ne sont plus tirés au hasard (MCAR) mais sont les 20 plus petits terrains du jeu
481 + parmi les lignes non altérées autrement (MNAR). Retourne le tableau de `comparer_imputations` pour cette variante (+ e.-t. des indicatrices)."""
482 + g = generer_brut(seed); b = g["brut"].copy(); ver = g["verite"]
483 + v = cas1.charger()["ventes"].set_index("id_vente")
484 + ids = ver.loc[ver["anomalie"] == "terrain_manquant", "id_vente_origine"]
485 + m = b["id_vente_origine"].isin(ids); b.loc[m, "superficie_terrain_m2"] = v.loc[b.loc[m, "id_vente_origine"], "superficie_terrain_m2"].values
486 + cand = b[~b["id_inscription"].isin(set(ver["id_inscription"]))].nsmallest(20, "superficie_terrain_m2")
487 + b.loc[cand.index, "superficie_terrain_m2"] = np.nan
488 + P = pipeline(b); tab = comparer_imputations(P)
489 + # erreur-type de l'indicatrice de terrain manquant
490 + base = geocodage_repli(P["regles"]["df"], []); d = imputer(base, None, "indicateur")
491 + S = scores_robustes(d, None); E = exclusions_regression(S["garde"], None); mi = _modele_indicateurs(E["df"])
492 + tab.loc["indicateur", "se_manq_terrain"] = mi["coef"].loc["manq_terrain", "se_hc1"]
493 + return tab
494 +
479 495 def sensibilite_seuils(b: pd.DataFrame, verite: pd.DataFrame) -> pd.DataFrame:
480 496 """Grille seuil MAD × seuil t : exclusions, σ, β_ln S, et détection des anomalies « à exclure » (ventes liées non déclarées)."""
481 497 a_detecter = set(verite.loc[verite["anomalie"] == "vente_liee_non_declaree", "id_inscription"])
@@ -519,7 +535,8 @@ def main(ecrire: bool = True, verbose: bool = True) -> dict:
519 535 A = appariement_role(b, role, verite); P0("== Appariement :\n", A["taux"].to_string(), "\n", A["pct"].to_string(), "\nambigus :", A["ambigus"], "; faux appariements :", A.get("faux_appariements"), A.get("faux_par_palier"))
520 536 I = comparer_imputations(P); P0("== Imputations :\n", I.round(4).to_string())
521 537 T = sensibilite_seuils(b, verite); P0("== Sensibilité aux seuils :\n", T.round(4).to_string(index=False))
522 − return {"brut": b, "verite": verite, "role": role, "pipeline": P, "modeles": C, "geocodage": G, "appariement": A, "imputations": I, "sensibilite": T, "masquage": M}
538 + V = variante_mnar_terrain(); P0("== Variante MNAR (ex. 14.7) :\n", V.round(4).to_string())
539 + return {"brut": b, "verite": verite, "role": role, "pipeline": P, "modeles": C, "geocodage": G, "appariement": A, "imputations": I, "sensibilite": T, "masquage": M, "mnar": V}
523 540
524 541 if __name__ == "__main__":
525 542 main(ecrire="--sans-ecriture" not in sys.argv)
modified analysis/python/cases/ch20_reconciliation.py +22 −3
@@ -145,6 +145,19 @@ def exhaustif(b: dict, P_aj16: np.ndarray, brut16: np.ndarray, k: int = 5) -> di
145 145 "rang_selection": int((vals < b["V_comp"]).sum()), "part_sous_hed": float((vals < b["V_hed"]).mean())}
146 146
147 147
148 +def elargissement(b: dict, P_aj16: np.ndarray, brut16: np.ndarray) -> pd.DataFrame:
149 + """Indication par comparaison lorsque l'on retient les k premiers candidats (k = 5 … 16) : valeur pondérée (inverse du brut),
150 + moyenne simple, résidu moyen (exp(e) − 1) des comparables retenus. Le modèle est fixe (arbitrage biais–variance du ch. 16)."""
151 + res = b["adj16"]["residu_log"].values
152 + rows = []
153 + for k in range(5, len(P_aj16) + 1):
154 + i = np.arange(k)
155 + rows.append({"k": k, "valeur_ponderee": valeur_ponderee(P_aj16[i], brut16[i]), "moyenne": float(P_aj16[i].mean()),
156 + "residu_moyen_pct": float(np.mean(np.exp(res[i]) - 1)), "brut_moyen": float(brut16[i].mean()),
157 + "distance_max": float(b["adj16"]["distance"].values[k - 1])})
158 + return pd.DataFrame(rows)
159 +
160 +
148 161 # ---------------------------------------------------------------- 3. indication par le coût (hypothèses déclarées)
149 162 def cout(b: dict, alpha: float = ALPHA_TERRAIN, vie: float = VIE_ECONOMIQUE, age_eff: float = AGE_EFFECTIF) -> dict:
150 163 devis, sujet = c21.charger()
@@ -208,6 +221,8 @@ def combinaison(b: dict, bs: dict, co: dict, graine: int = GRAINE) -> dict:
208 221 arrondi = ARRONDI_CONCLUSION
209 222 conclusion = round(V_raisonnee / arrondi) * arrondi
210 223 fourchette90 = (math.floor(q[1] / arrondi) * arrondi, math.ceil(q[3] / arrondi) * arrondi)
224 + # fourchette de valeur conclue : encadrée par les deux indications crédibles (comparaison, hédonique), arrondies vers l'extérieur
225 + fourchette_conclue = (math.floor(min(V[0], V[1]) / arrondi) * arrondi, math.ceil(max(V[0], V[1]) / arrondi) * arrondi)
211 226 # courbes de sensibilité aux poids : w_hed de 0 à 1 − w_cout, pour trois poids du coût
212 227 wh = np.linspace(0, 1, 101)
213 228 courbes = {wc: np.array([ (1 - h - wc) * V[0] + h * V[1] + wc * V[2] if h <= 1 - wc + 1e-12 else np.nan for h in wh]) for wc in (0.0, 0.05, 0.10)}
@@ -216,7 +231,8 @@ def combinaison(b: dict, bs: dict, co: dict, graine: int = GRAINE) -> dict:
216 231 w_opt_deux = {r: poids_optimaux_deux(s[0], s[1], r) for r in (0.0, rho, 0.3, 0.6, 0.9)}
217 232 return {"V": V, "s": s, "sd_spec_hed": sd_spec, "sd_pred_boot": sd_pred_boot, "rho": rho, "Sigma": Sigma, "schemas": schemas, "tableau": tab,
218 233 "Vb": Vb, "cout_b": cout_b, "q": {"q025": q[0], "q05": q[1], "q50": q[2], "q95": q[3], "q975": q[4]}, "V_raisonnee": V_raisonnee,
219 − "conclusion": conclusion, "fourchette90": fourchette90, "wh": wh, "courbes": courbes, "var_deux": var_deux, "w_opt_deux": w_opt_deux,
234 + "conclusion": conclusion, "fourchette90": fourchette90, "fourchette_conclue": fourchette_conclue,
235 + "largeur_relative_conclue": (fourchette_conclue[1] - fourchette_conclue[0]) / conclusion, "wh": wh, "courbes": courbes, "var_deux": var_deux, "w_opt_deux": w_opt_deux,
220 236 "ecart_conclusion_comp": conclusion / b["V_comp"] - 1, "ecart_conclusion_hed": conclusion / b["V_hed"] - 1}
221 237
222 238
@@ -283,11 +299,12 @@ def tout(verbose: bool = False) -> dict:
283 299 b = base()
284 300 bs = bootstrap(b)
285 301 ex = exhaustif(b, bs["P_aj16"], bs["brut16"])
302 + el = elargissement(b, bs["P_aj16"], bs["brut16"])
286 303 co = cout(b)
287 304 cb = combinaison(b, bs, co)
288 305 vt = verite(b)
289 306 exo = exercices(b, bs, co, cb)
290 − r = {"base": b, "bootstrap": bs, "exhaustif": ex, "cout": co, "combinaison": cb, "verite": vt, "exercices": exo}
307 + r = {"base": b, "bootstrap": bs, "exhaustif": ex, "elargissement": el, "cout": co, "combinaison": cb, "verite": vt, "exercices": exo}
291 308 _CACHE["r"] = r
292 309 if verbose:
293 310 imprimer(r)
@@ -314,6 +331,8 @@ def imprimer(r: dict) -> None:
314 331 print(f"n = {ex['n']} ; moyenne {f(ex['moyenne'])} ; é.-t. {f(ex['ecart_type'])} ; min {f(ex['min'])} ; max {f(ex['max'])} ; q5 {f(ex['q05'])} ; q50 {f(ex['q50'])} ; q95 {f(ex['q95'])}")
315 332 print(f"rang de la sélection du ch. 16 : {ex['rang_selection']} sous-ensembles sur {ex['n']} donnent moins ; part sous la prédiction hédonique {ex['part_sous_hed']:.3f}")
316 333 print(ex["par_n_retenus"].round(0).to_string())
334 + print("\n=== Élargissement k = 5 … 16 (modèle fixe) ===")
335 + print(r["elargissement"].round(4).to_string(index=False))
317 336 print("\n=== Coût (hypothèses : α = %.2f, vie = %d ans, âge effectif = %d) ===" % (co["alpha"], co["vie"], co["age_eff"]))
318 337 print(f"C_N {f(co['CN'])} ; V_L {f(co['V_L'])} ; D âge-vie {f(co['D'])} ({co['taux_D']:.3f}) ; contribution bâtiment {f(co['contribution_batiment'])} ; valeur {f(co['valeur'])}")
319 338 print(f"C_N + V_L = {f(co['cout_neuf_plus_terrain'])} ; extraction du marché à la prédiction hédonique : D = {f(co['extraction_hed']['depreciation'])} (taux {co['extraction_hed']['taux']:.3f}) ; à la comparaison : D = {f(co['extraction_comp']['depreciation'])} (taux {co['extraction_comp']['taux']:.3f})")
@@ -325,7 +344,7 @@ def imprimer(r: dict) -> None:
325 344 print(cb["tableau"].round(3).to_string(index=False))
326 345 print(f"poids optimal de la comparaison à deux indications selon ρ : { {round(k,3): round(v,3) for k, v in cb['w_opt_deux'].items()} }")
327 346 q = cb["q"]
328 − print(f"conclusion raisonnée {f(cb['V_raisonnee'])} → arrondie {f(cb['conclusion'])} ; bootstrap de la combinaison : q2.5 {f(q['q025'])} q5 {f(q['q05'])} q50 {f(q['q50'])} q95 {f(q['q95'])} q97.5 {f(q['q975'])} ; fourchette 90 % arrondie {f(cb['fourchette90'][0])}–{f(cb['fourchette90'][1])}")
347 + print(f"conclusion raisonnée {f(cb['V_raisonnee'])} → arrondie {f(cb['conclusion'])} ; bootstrap de la combinaison : q2.5 {f(q['q025'])} q5 {f(q['q05'])} q50 {f(q['q50'])} q95 {f(q['q95'])} q97.5 {f(q['q975'])} ; fourchette 90 % arrondie {f(cb['fourchette90'][0])}–{f(cb['fourchette90'][1])} ; fourchette conclue {f(cb['fourchette_conclue'][0])}–{f(cb['fourchette_conclue'][1])} ({cb['largeur_relative_conclue']:.3f})")
329 348 print(f"écart conclusion / comparaison {cb['ecart_conclusion_comp']:.4f} ; / hédonique {cb['ecart_conclusion_hed']:.4f}")
330 349 print("\n=== Loi de génération ===")
331 350 print(f"c_hat {vt['c_hat']:.4f} (e.-t. {vt['se_c']:.4f}, σ_e {vt['sigma_e']:.4f}) ; V_vraie {f(vt['V_true'])} [{f(vt['V_true_bas'])} ; {f(vt['V_true_haut'])}] ; écart comparaison {vt['ecart_comp']:.4f} ; hédonique {vt['ecart_hed']:.4f} ; conclusion {cb['conclusion']/vt['V_true']-1:.4f}")
added analysis/python/cases/ch34_tal.py +440 −0
@@ -0,0 +1,440 @@
1 +#!/usr/bin/env python3
2 +"""Chapitre 34 — TAL, encadrement locatif et scénarios : calculs du chapitre (cas 3, 12 logements, Montréal,
3 +DONNÉES SYNTHÉTIQUES, sauf l'IPC du Québec qui est réel : Statistique Canada, tableau 18-10-0004-01).
4 +
5 +Contenu :
6 + 1. Pourcentage de base de la fixation (Règlement sur les critères de fixation de loyer, art. 3.1, tel que remplacé par
7 + D. 1455-2025) recalculé à partir des moyennes annuelles de l'IPC du Québec (vérification des 4,5 % de 2025 et des
8 + 3,1 % de 2026 publiés par le TAL) et historique « contrefactuel » de la formule depuis 2006.
9 + 2. Calcul pas à pas d'une augmentation de loyer selon la méthode en vigueur depuis le 1er janvier 2026 (art. 3, 3.1, 4,
10 + 5 et 15.1) pour un logement du cas 3, et reconstitution pédagogique de la structure de l'ancienne méthode.
11 + 3. Trajectoire du loyer contractuel sous contrainte : forme fermée de l'espérance du rapport loyer/marché sous rotation
12 + aléatoire, vérifiée par simulation, et comparaison avec la projection déterministe de traite.cas3.projection_dcf.
13 + 4. Scénarios : fonction locale `projection_scenario` qui reprend la logique de traite.cas3.projection_dcf avec des
14 + paramètres variés (rotation, pourcentage de base, croissance du marché, rattrapage partiel à la relocation,
15 + récupération de 5 % des dépenses d'immobilisation, mode déterministe ou espérance) SANS modifier traite/ ; valeur
16 + DCF de chaque scénario par traite.cas3.valeur_dcf.
17 +Usage (racine du dépôt) : PYTHONPATH=analysis/python python3 analysis/python/cases/ch34_tal.py
18 +Les fonctions sont importées par analysis/python/figures/ch34_*.py et analysis/python/tests/test_ch34.py."""
19 +from __future__ import annotations
20 +import sys, csv, copy, math, pathlib
21 +import numpy as np, pandas as pd
22 +sys.path.insert(0, str(pathlib.Path(__file__).resolve().parents[1]))
23 +from traite import cas3, finance as fin
24 +
25 +ROOT = pathlib.Path(__file__).resolve().parents[3]
26 +IPC_RAW = ROOT / "data/raw/statcan_18100004/18100004.csv"
27 +IPC_PROC = ROOT / "data/processed/ipc_quebec_mensuel.csv"
28 +VECTEUR_IPC_QC = "v41691783" # IPC Québec, Ensemble, 2002 = 100, non désaisonnalisé (tableau 18-10-0004-01)
29 +PCT_TAL_PUBLIES = {2025: 0.045, 2026: 0.031} # TAL, « Pourcentages applicables », page lue le 2026-10-31
30 +GRAINE = 20261031
31 +
32 +# ============================================================================ 1. IPC du Québec et pourcentage de base (art. 3.1)
33 +
34 +def ipc_quebec_mensuel() -> pd.DataFrame:
35 + """Série mensuelle de l'IPC du Québec (Ensemble, 2002 = 100). Extraite une fois du fichier brut de Statistique Canada
36 + (170 Mo) vers data/processed/ipc_quebec_mensuel.csv ; les appels suivants lisent le fichier traité."""
37 + if IPC_PROC.exists():
38 + return pd.read_csv(IPC_PROC)
39 + if not IPC_RAW.exists():
40 + raise FileNotFoundError(f"ni {IPC_PROC} ni {IPC_RAW} : télécharger le tableau 18-10-0004-01 (voir data/raw/SOURCES.md)")
41 + rows = []
42 + with IPC_RAW.open(encoding="utf-8-sig", newline="") as f:
43 + for r in csv.reader(f, delimiter=";"):
44 + if len(r) > 10 and r[8] == VECTEUR_IPC_QC and r[10] not in ("", ".."):
45 + rows.append({"periode": r[0], "ipc": float(r[10])})
46 + df = pd.DataFrame(rows).sort_values("periode").reset_index(drop=True)
47 + IPC_PROC.parent.mkdir(parents=True, exist_ok=True)
48 + df.to_csv(IPC_PROC, index=False)
49 + return df
50 +
51 +def ipc_moyennes_annuelles(df: pd.DataFrame | None = None) -> pd.Series:
52 + """Moyenne annuelle de l'IPC (années civiles complètes seulement, 12 observations)."""
53 + df = ipc_quebec_mensuel() if df is None else df
54 + d = df.assign(an=df["periode"].str[:4].astype(int))
55 + g = d.groupby("an")["ipc"].agg(["mean", "count"])
56 + return g.loc[g["count"] == 12, "mean"]
57 +
58 +def pourcentage_base(annee_fixation: int, moy: pd.Series | None = None) -> dict:
59 + """Art. 3.1 : { [(A−B)/B] + [(B−C)/C] + [(C−D)/D] } / 3, A = IPC moyen de la période de référence (année civile précédant le
60 + terme du bail pour un bail se terminant du 1er avril au 31 décembre, art. 1), B, C, D = les trois périodes de 12 mois
61 + antérieures ; résultat négatif réputé nul (3e al.)."""
62 + moy = ipc_moyennes_annuelles() if moy is None else moy
63 + A, B, C, D = (float(moy.loc[annee_fixation - k]) for k in (1, 2, 3, 4))
64 + v1, v2, v3 = (A - B) / B, (B - C) / C, (C - D) / D
65 + brut = (v1 + v2 + v3) / 3
66 + return {"annee": annee_fixation, "A": A, "B": B, "C": C, "D": D, "var_A_B": v1, "var_B_C": v2, "var_C_D": v3,
67 + "pct_brut": brut, "pct_base": max(brut, 0.0), "pct_publie": PCT_TAL_PUBLIES.get(annee_fixation)}
68 +
69 +def historique_pourcentage_base(debut: int = 2006, fin: int = 2026) -> pd.DataFrame:
70 + """Pour chaque année de fixation : inflation de la période de référence (A−B)/B et pourcentage de base qu'aurait donné la
71 + formule de l'art. 3.1 (contrefactuel avant 2026, car la formule n'existait pas ; valeurs réelles 2025-2026 à comparer aux
72 + pourcentages publiés par le TAL)."""
73 + moy = ipc_moyennes_annuelles()
74 + rows = []
75 + for an in range(debut, fin + 1):
76 + if an - 4 in moy.index and an - 1 in moy.index:
77 + p = pourcentage_base(an, moy)
78 + rows.append({"annee": an, "inflation_ref": p["var_A_B"], "pct_formule": p["pct_base"], "pct_publie": p["pct_publie"]})
79 + return pd.DataFrame(rows)
80 +
81 +# ============================================================================ 2. Fixation du loyer : nouvelle méthode (2026) et ancienne
82 +
83 +def arrondi_art_15_1(x: float) -> float:
84 + """Art. 15.1 : fraction < 0,50 $ → dollar inférieur ; ≥ 0,50 $ → dollar supérieur (arrondi « demi vers le haut », non bancaire)."""
85 + return math.floor(x + 0.5)
86 +
87 +def part_attribuable(loyer_mensuel: float, revenus_annuels: float) -> float:
88 + """Art. 5 : part attribuable au logement = loyer au terme du bail / revenus de l'immeuble (ici loyer annuel / revenus annuels,
89 + loyers résidentiels et autres revenus compris ; aucune part pour un local non résidentiel)."""
90 + return loyer_mensuel * 12 / revenus_annuels
91 +
92 +def taux_net(variation: float, pct_base: float) -> float:
93 + """Art. 4, 2e al. : si la variation est positive, on en soustrait la valeur absolue du pourcentage de base (brut, avant
94 + plancher à zéro) ; si le résultat est négatif, il est réputé nul ; une variation négative ne donne rien."""
95 + if variation <= 0:
96 + return 0.0
97 + return max(variation - abs(pct_base), 0.0)
98 +
99 +def fixation_nouvelle(loyer: float, pct_base: float, part: float, taxes_mun: tuple[float, float], taxes_scol: tuple[float, float],
100 + assurances: tuple[float, float], immobilisations: float = 0.0, nouveaux_services: float = 0.0,
101 + pct_immobilisations: float = 0.05) -> dict:
102 + """Augmentation mensuelle selon les art. 3, 3.1, 4, 5 et 15.1 du Règlement (version en vigueur le 1er janvier 2026).
103 + taxes_mun = (exigibles pendant la période de référence, exigibles l'année suivante) ; taxes_scol et assurances =
104 + (année précédant la période de référence, période de référence). Montants annuels pour l'immeuble entier.
105 + Reconstitution pédagogique : le formulaire officiel du TAL (Outil de calcul) fait foi."""
106 + base = loyer * pct_base
107 + v_mun = (taxes_mun[1] - taxes_mun[0]) / taxes_mun[0]
108 + v_scol = (taxes_scol[1] - taxes_scol[0]) / taxes_scol[0]
109 + v_ass = (assurances[1] - assurances[0]) / assurances[0]
110 + n_mun, n_scol, n_ass = taux_net(v_mun, pct_base), taux_net(v_scol, pct_base), taux_net(v_ass, pct_base)
111 + m_mun = n_mun * taxes_mun[0] / 12 * part # art. 4, 3e al., 1° : montant mensuel de la période de référence
112 + m_scol = n_scol * taxes_scol[0] / 12 * part # art. 4, 3e al., 2° : montant mensuel de la période précédente
113 + m_ass = n_ass * assurances[0] / 12 * part
114 + m_immo = pct_immobilisations * immobilisations * part / 12
115 + m_serv = nouveaux_services * part / 12
116 + total = base + m_mun + m_scol + m_ass + m_immo + m_serv
117 + nouveau = arrondi_art_15_1(loyer + total)
118 + return {"loyer": loyer, "pct_base": pct_base, "part": part, "base": base,
119 + "var_taxes_mun": v_mun, "net_taxes_mun": n_mun, "m_taxes_mun": m_mun,
120 + "var_taxes_scol": v_scol, "net_taxes_scol": n_scol, "m_taxes_scol": m_scol,
121 + "var_assurances": v_ass, "net_assurances": n_ass, "m_assurances": m_ass,
122 + "immobilisations": immobilisations, "m_immobilisations": m_immo, "m_services": m_serv,
123 + "total_brut": total, "augmentation": nouveau - loyer, "nouveau_loyer": nouveau, "taux_effectif": (nouveau - loyer) / loyer}
124 +
125 +PCT_ANCIENNE_2026 = {"electricite": 0.029, "gaz": 0.137, "mazout": 0.049, "entretien": 0.057, "services_immeuble": 0.034,
126 + "gestion": 0.072, "revenu_net": 0.057, "immobilisations": 0.04}
127 +PCT_ANCIENNE_2025 = {"electricite": 0.029, "gaz": -0.058, "mazout": -0.029, "entretien": 0.069, "services_immeuble": 0.045,
128 + "gestion": 0.082, "revenu_net": 0.069, "immobilisations": 0.047}
129 +
130 +def fixation_ancienne(loyer: float, part: float, revenus: float, depenses: dict, pct: dict, taxes_mun: tuple[float, float],
131 + taxes_scol: tuple[float, float], assurances: tuple[float, float], immobilisations: float = 0.0,
132 + gestion_forfait: float = 0.05) -> dict:
133 + """Reconstitution PÉDAGOGIQUE de la structure de l'ancienne méthode (avis donnés avant le 1er janvier 2026 ; formulaire de
134 + calcul du TAL) : chaque poste de dépenses d'exploitation est ajusté par le pourcentage de sa composante (énergie,
135 + entretien, services, gestion forfaitaire de 5 % des revenus), les taxes et assurances passent par leur variation réelle,
136 + le revenu net (revenus − dépenses d'exploitation) est majoré du pourcentage « revenu net » et les dépenses d'immobilisation
137 + du pourcentage « dépenses d'immobilisation » ; la somme est répartie selon la part attribuable. Le détail exact des lignes
138 + du formulaire officiel (et ses exclusions) fait foi ; `depenses` : électricité, gaz, mazout, entretien, services_immeuble."""
139 + gestion = gestion_forfait * revenus
140 + aj = {k: depenses.get(k, 0.0) * pct.get(k, 0.0) for k in ("electricite", "gaz", "mazout", "entretien", "services_immeuble")}
141 + aj["gestion"] = gestion * pct["gestion"]
142 + aj["taxes_mun"] = taxes_mun[1] - taxes_mun[0]
143 + aj["taxes_scol"] = taxes_scol[1] - taxes_scol[0]
144 + aj["assurances"] = assurances[1] - assurances[0]
145 + dep_expl = sum(depenses.get(k, 0.0) for k in ("electricite", "gaz", "mazout", "entretien", "services_immeuble")) + gestion + taxes_mun[0] + taxes_scol[0] + assurances[0]
146 + revenu_net = revenus - dep_expl
147 + aj["revenu_net"] = revenu_net * pct["revenu_net"]
148 + aj["immobilisations"] = immobilisations * pct["immobilisations"]
149 + total_immeuble = sum(aj.values())
150 + mensuel = total_immeuble * part / 12
151 + nouveau = arrondi_art_15_1(loyer + mensuel)
152 + return {"ajustements_immeuble": aj, "depenses_exploitation": dep_expl, "revenu_net": revenu_net, "total_immeuble": total_immeuble,
153 + "mensuel_brut": mensuel, "augmentation": nouveau - loyer, "nouveau_loyer": nouveau, "taux_effectif": (nouveau - loyer) / loyer}
154 +
155 +def donnees_immeuble_cas3(d: dict | None = None) -> dict:
156 + """Montants de l'immeuble du cas 3 utiles à la fixation : revenus (loyers contractuels + autres revenus), taxes 2024-2026,
157 + assurances 2024-2025, dépenses 2025 par composante (données synthétiques)."""
158 + d = cas3.charger() if d is None else d
159 + etat, dep, hyp = d["etat"], d["dep"], d["hyp"]
160 + p = dep.pivot(index="poste", columns="annee", values="montant").astype(float)
161 + rev = cas3.revenus(etat, hyp)
162 + taxes_mun_2026 = hyp["valeur_role_2026"] * hyp["taux_taxe_fonciere"]
163 + return {"revenus": rev["rbp_contrat"] + rev["autres_revenus"], "rbp_contrat": rev["rbp_contrat"], "autres_revenus": rev["autres_revenus"],
164 + "taxes_mun": {2024: p.loc["Taxes municipales", 2024], 2025: p.loc["Taxes municipales", 2025], 2026: taxes_mun_2026},
165 + "taxes_scol": {2024: p.loc["Taxe scolaire", 2024], 2025: p.loc["Taxe scolaire", 2025]},
166 + "assurances": {2024: p.loc["Assurances", 2024], 2025: p.loc["Assurances", 2025]},
167 + "depenses_2025": {"electricite": p.loc["Électricité (aires communes)", 2025], "gaz": p.loc["Chauffage (aires communes, gaz)", 2025], "mazout": 0.0,
168 + "entretien": p.loc["Entretien et réparations", 2025],
169 + "services_immeuble": p.loc["Conciergerie", 2025] + p.loc["Déneigement et paysagement", 2025]},
170 + "etat": etat}
171 +
172 +def exemple_fixation(numero: int = 103, d: dict | None = None) -> dict:
173 + """Exemple intégral : logement `numero` du cas 3, pourcentages 2026 (3,1 % ; 5 %).
174 + Variante A : données du dossier (taxes municipales 2025 → 2026 normalisées, scolaire et assurances 2024 → 2025, aucune
175 + dépense d'immobilisation en 2025). Variante B (illustrative) : mêmes loyers, mais hausse de taxes municipales de 8 %,
176 + hausse d'assurances de 10 % et dépenses d'immobilisation de 85 000 $ (les travaux de maçonnerie et de balcons programmés
177 + en 2027 dans le dossier, supposés réalisés pendant la période de référence)."""
178 + im = donnees_immeuble_cas3(d)
179 + etat = im["etat"]; row = etat[etat["logement"] == numero].iloc[0]
180 + loyer = float(row["loyer_mensuel_contrat"]); marche = float(row["loyer_mensuel_marche"])
181 + part = part_attribuable(loyer, im["revenus"])
182 + pb = PCT_TAL_PUBLIES[2026]
183 + A = fixation_nouvelle(loyer, pb, part, (im["taxes_mun"][2025], im["taxes_mun"][2026]), (im["taxes_scol"][2024], im["taxes_scol"][2025]),
184 + (im["assurances"][2024], im["assurances"][2025]), immobilisations=0.0)
185 + B = fixation_nouvelle(loyer, pb, part, (im["taxes_mun"][2025], im["taxes_mun"][2025] * 1.08), (im["taxes_scol"][2024], im["taxes_scol"][2025]),
186 + (im["assurances"][2024], im["assurances"][2024] * 1.10), immobilisations=85_000.0)
187 + anc = fixation_ancienne(loyer, part, im["revenus"], im["depenses_2025"], PCT_ANCIENNE_2026, (im["taxes_mun"][2025], im["taxes_mun"][2026]),
188 + (im["taxes_scol"][2024], im["taxes_scol"][2025]), (im["assurances"][2024], im["assurances"][2025]), immobilisations=0.0)
189 + anc_B = fixation_ancienne(loyer, part, im["revenus"], im["depenses_2025"], PCT_ANCIENNE_2026, (im["taxes_mun"][2025], im["taxes_mun"][2025] * 1.08),
190 + (im["taxes_scol"][2024], im["taxes_scol"][2025]), (im["assurances"][2024], im["assurances"][2024] * 1.10), immobilisations=85_000.0)
191 + return {"numero": numero, "loyer": loyer, "marche": marche, "ecart_marche": marche - loyer, "ecart_marche_pct": marche / loyer - 1,
192 + "part": part, "revenus": im["revenus"], "immeuble": im, "A": A, "B": B, "ancienne_A": anc, "ancienne_B": anc_B,
193 + "annees_rattrapage_base": math.log(marche / loyer) / math.log(1 + pb)}
194 +
195 +# ============================================================================ 3. Trajectoire espérée sous rotation aléatoire
196 +
197 +def ratio_espere(q0: float, rho: float, g_c: float, g_m: float, t: int | np.ndarray) -> np.ndarray:
198 + """E[q_t] = q0·a^t + ρ·(1 − a^t)/(1 − a), a = (1 − ρ)·δ, δ = (1 + g_c)/(1 + g_m) : rapport espéré loyer/marché d'un logement dont
199 + le locataire part chaque année avec probabilité ρ (i.i.d.), le logement étant alors reloué au marché, puis indexé à g_c."""
200 + delta = (1 + g_c) / (1 + g_m); a = (1 - rho) * delta
201 + t = np.asarray(t, dtype=float)
202 + return q0 * a ** t + rho * (1 - a ** t) / (1 - a)
203 +
204 +def ratio_limite(rho: float, g_c: float, g_m: float) -> float:
205 + """q_∞ = ρ / (1 − (1 − ρ)·δ)."""
206 + delta = (1 + g_c) / (1 + g_m)
207 + return rho / (1 - (1 - rho) * delta)
208 +
209 +def ratio_sans_relocation(q0: float, g_c: float, g_m: float, t) -> np.ndarray:
210 + delta = (1 + g_c) / (1 + g_m)
211 + return q0 * delta ** np.asarray(t, dtype=float)
212 +
213 +def simulation_ratio(q0: float, rho: float, g_c: float, g_m: float, T: int, n: int = 20_000, graine: int = GRAINE) -> np.ndarray:
214 + """Vérification Monte-Carlo de la forme fermée : n trajectoires, départ i.i.d. de probabilité ρ par année."""
215 + rng = np.random.default_rng(graine); delta = (1 + g_c) / (1 + g_m)
216 + q = np.full(n, q0, dtype=float); out = [q.mean()]
217 + for _ in range(T):
218 + q = q * delta
219 + depart = rng.random(n) < rho
220 + q[depart] = 1.0
221 + out.append(q.mean())
222 + return np.array(out)
223 +
224 +def trajectoire_cas3_esperee(d: dict | None = None, rho: float | None = None, g_c: float | None = None, g_m: float | None = None, T: int = 10) -> pd.DataFrame:
225 + """Espérance du RBP contractuel du cas 3 sous rotation aléatoire (forme fermée appliquée logement par logement) et
226 + comparaison avec la projection déterministe « les plus sous le marché d'abord » de traite.cas3.projection_dcf."""
227 + d = cas3.charger() if d is None else d
228 + hyp = d["hyp"]; etat = d["etat"]
229 + rho = hyp["rotation_annuelle"] if rho is None else rho; g_c = hyp["croissance_loyers_contrat"] if g_c is None else g_c
230 + g_m = hyp["croissance_loyers_marche"] if g_m is None else g_m
231 + L = etat["loyer_mensuel_contrat"].astype(float).values; M = etat["loyer_mensuel_marche"].astype(float).values
232 + t = np.arange(0, T + 1)
233 + marche_t = np.outer((1 + g_m) ** t, M) # (T+1) × n
234 + q_t = np.column_stack([ratio_espere(L[i] / M[i], rho, g_c, g_m, t) for i in range(len(L))])
235 + rbp_esp = (q_t * marche_t).sum(axis=1) * 12
236 + rbp_marche = marche_t.sum(axis=1) * 12
237 + rbp_sans = (np.column_stack([ratio_sans_relocation(L[i] / M[i], g_c, g_m, t) for i in range(len(L))]) * marche_t).sum(axis=1) * 12
238 + proj = projection_scenario(etat, d["dep"], hyp, rotation=rho, g_c=g_c, g_m=g_m, horizon=T)
239 + det = proj.loc[proj["annee"] <= T + 1, "rbp"].values # années 1..T+1 ↔ t = 0..T
240 + return pd.DataFrame({"t": t, "exercice": int(hyp["date_evaluation"][:4]) + t, "rbp_marche": rbp_marche, "rbp_esperance": rbp_esp,
241 + "rbp_deterministe": det, "rbp_sans_relocation": rbp_sans,
242 + "ratio_esperance": rbp_esp / rbp_marche, "ratio_deterministe": det / rbp_marche, "ratio_sans": rbp_sans / rbp_marche})
243 +
244 +# ============================================================================ 4. Scénarios : projection locale paramétrée et valeur DCF
245 +
246 +def projection_scenario(etat: pd.DataFrame, dep: pd.DataFrame, hyp: dict, rotation: float | None = None, g_c: float | None = None,
247 + g_m: float | None = None, fraction_rattrapage: float = 1.0, mode: str = "deterministe",
248 + recuperation_capex: bool = False, pct_recuperation: float = 0.05, horizon: int | None = None) -> pd.DataFrame:
249 + """Reprend la logique de traite.cas3.projection_dcf (indexation des baux à g_c, relocation d'une fraction `rotation` des
250 + logements au loyer de marché croissant à g_m, dépenses +g_d, réserve constante, CAPEX aux années prévues) avec :
251 + - mode « deterministe » : k = round(rotation × n) logements, les plus sous le marché d'abord (identique au module commun
252 + lorsque les paramètres sont ceux du dossier) ;
253 + - mode « esperance » : chaque loyer devient (1 − ρ)·L(1+g_c) + ρ·[L(1+g_c) + φ·(M − L(1+g_c))], soit l'espérance sous départ
254 + i.i.d. (rotation fractionnaire, sans arrondi) ;
255 + - fraction_rattrapage φ ∈ [0, 1] : à la relocation, le nouveau loyer n'absorbe qu'une fraction φ de l'écart au marché
256 + (φ < 1 représente un recours en fixation du nouveau locataire, art. 1896 et 1950 C.c.Q., ou un rattrapage partiel) ;
257 + - recuperation_capex : l'année qui suit une dépense d'immobilisation, les loyers sont majorés de pct_recuperation × CAPEX
258 + (art. 3, par. 5° du Règlement), réparti selon la part attribuable (proportionnellement aux loyers).
259 + Colonnes identiques à projection_dcf, plus `ecart_marche` et `logements_reloues`."""
260 + T = horizon or hyp["horizon_ans"]
261 + base = cas3.rne(etat, dep, hyp)
262 + loyers = etat["loyer_mensuel_contrat"].astype(float).values.copy(); marche = etat["loyer_mensuel_marche"].astype(float).values.copy()
263 + rot = hyp["rotation_annuelle"] if rotation is None else rotation
264 + gc = hyp["croissance_loyers_contrat"] if g_c is None else g_c
265 + gm = hyp["croissance_loyers_marche"] if g_m is None else g_m
266 + g_d = hyp["croissance_depenses"]; v = hyp["taux_inoccupation_pertes"]; autres = hyp["autres_revenus_annuels"]
267 + de0 = base["depenses_hors_reserve"]; res0 = base["reserve"]
268 + annee0 = int(hyp["date_evaluation"][:4]); rows = []; capex_prec = 0.0
269 + for t in range(1, T + 2):
270 + k = 0
271 + if t > 1:
272 + loyers *= (1 + gc); marche *= (1 + gm)
273 + if mode == "deterministe":
274 + k = int(round(rot * len(loyers)))
275 + idx = np.argsort(loyers / marche)[:k]
276 + loyers[idx] = loyers[idx] + fraction_rattrapage * (marche[idx] - loyers[idx])
277 + elif mode == "esperance":
278 + loyers = (1 - rot) * loyers + rot * (loyers + fraction_rattrapage * (marche - loyers))
279 + k = rot * len(loyers)
280 + else:
281 + raise ValueError(mode)
282 + if recuperation_capex and capex_prec > 0:
283 + loyers *= 1 + pct_recuperation * capex_prec / (loyers.sum() * 12)
284 + rbp = loyers.sum() * 12; rbe = fin.revenu_brut_effectif(rbp, v, autres * (1 + g_d) ** (t - 1))
285 + de = de0 * (1 + g_d) ** (t - 1); res = res0; noi = rbe - de - res
286 + capex = float(hyp["capex"].get(str(annee0 + t - 1), 0.0)) if t <= T else 0.0
287 + rows.append({"annee": t, "exercice": annee0 + t - 1, "rbp": rbp, "rbp_marche": marche.sum() * 12, "rbe": rbe, "depenses": de, "reserve": res,
288 + "rne": noi, "capex": capex, "flux_avant_dette": noi - capex, "ecart_marche": marche.sum() * 12 - rbp, "logements_reloues": k})
289 + capex_prec = capex
290 + return pd.DataFrame(rows)
291 +
292 +def valeur_scenario(d: dict, **kw) -> dict:
293 + """Projection locale + valeur DCF (traite.cas3.valeur_dcf, hypothèses r, R_T, c_v, T du dossier) + mesures de l'écart au marché."""
294 + hyp = d["hyp"]; T = hyp["horizon_ans"]; r = hyp["taux_actualisation"]
295 + proj = projection_scenario(d["etat"], d["dep"], hyp, **kw)
296 + val = cas3.valeur_dcf(proj, hyp)
297 + p = proj[proj["annee"] <= T]
298 + va_ecart = float(sum(e * fin.facteur_va(r, int(t)) for e, t in zip(p["ecart_marche"], p["annee"])))
299 + return {**val, "projection": proj, "rne_1": float(proj.loc[proj["annee"] == 1, "rne"].iloc[0]), "rbp_1": float(proj.loc[proj["annee"] == 1, "rbp"].iloc[0]),
300 + "rbp_T1": float(proj.loc[proj["annee"] == T + 1, "rbp"].iloc[0]), "ecart_T1_pct": float(proj.loc[proj["annee"] == T + 1, "ecart_marche"].iloc[0] / proj.loc[proj["annee"] == T + 1, "rbp_marche"].iloc[0]),
301 + "va_ecart_marche": va_ecart, "croissance_rbp": (float(proj.loc[proj["annee"] == T + 1, "rbp"].iloc[0]) / float(proj.loc[proj["annee"] == 1, "rbp"].iloc[0])) ** (1 / T) - 1}
302 +
303 +def grille_scenarios(d: dict, rotations=(0.05, 0.15, 0.25), bases=(0.02, 0.031, 0.04), **kw) -> pd.DataFrame:
304 + """Valeur DCF pour chaque combinaison rotation × pourcentage de base (autres hypothèses du dossier inchangées)."""
305 + V0 = d["dcf"]["valeur"] if "dcf" in d else None
306 + rows = []
307 + for rho in rotations:
308 + for b in bases:
309 + v = valeur_scenario(d, rotation=rho, g_c=b, **kw)
310 + rows.append({"rotation": rho, "base": b, "logements_par_an": int(round(rho * len(d["etat"]))) if kw.get("mode", "deterministe") == "deterministe" else rho * len(d["etat"]),
311 + "valeur": v["valeur"], "rne_T1": v["rne_T1"], "rbp_T1": v["rbp_T1"], "ecart_T1_pct": v["ecart_T1_pct"],
312 + "va_ecart_marche": v["va_ecart_marche"], "part_terminale": v["part_terminale"], "croissance_rbp": v["croissance_rbp"]})
313 + g = pd.DataFrame(rows)
314 + if V0 is not None:
315 + g["ecart_pct_base"] = g["valeur"] / V0 - 1
316 + return g
317 +
318 +def scenarios_complementaires(d: dict) -> pd.DataFrame:
319 + """Scénarios narratifs autour du cas de base (rotation 15 %, base 3,1 %, marché 3,5 %, relocation intégrale au marché)."""
320 + V0 = d["dcf"]["valeur"]
321 + defs = [
322 + ("Base (dossier)", {}),
323 + ("Base, mode espérance (rotation fractionnaire)", {"mode": "esperance"}),
324 + ("Indexation seule (aucune relocation)", {"rotation": 0.0}),
325 + ("Encadrement serré (ch. 30 : base 2 %, rotation 8 %, marché 3 %)", {"g_c": 0.02, "rotation": 0.08, "g_m": 0.03}),
326 + ("Rattrapage partiel à la relocation (φ = 0,5)", {"fraction_rattrapage": 0.5}),
327 + ("Encadrement renforcé (base 2 %, rotation 8 %, φ = 0,5, marché 3 %)", {"g_c": 0.02, "rotation": 0.08, "fraction_rattrapage": 0.5, "g_m": 0.03}),
328 + ("Récupération de 5 % des dépenses d'immobilisation", {"recuperation_capex": True}),
329 + ("Rotation 25 % et base 4 %", {"rotation": 0.25, "g_c": 0.04}),
330 + ("Tout au marché dès l'année 1 (borne supérieure théorique)", {"rotation": 1.0}),
331 + ]
332 + rows = []
333 + for nom, kw in defs:
334 + v = valeur_scenario(d, **kw)
335 + rows.append({"scenario": nom, **{k: kw.get(k) for k in ("rotation", "g_c", "g_m", "fraction_rattrapage", "mode", "recuperation_capex")},
336 + "rne_1": v["rne_1"], "rne_T1": v["rne_T1"], "rbp_T1": v["rbp_T1"], "ecart_T1_pct": v["ecart_T1_pct"], "va_ecart_marche": v["va_ecart_marche"],
337 + "valeur": v["valeur"], "ecart_pct_base": v["valeur"] / V0 - 1, "part_terminale": v["part_terminale"]})
338 + return pd.DataFrame(rows)
339 +
340 +def elasticites(d: dict, h_rot: float = 0.05, h_base: float = 0.005) -> dict:
341 + """Variation de valeur par point de rotation et par point de pourcentage de base autour de la base (différences centrées,
342 + mode espérance pour éviter les sauts d'arrondi du nombre de logements)."""
343 + hyp = d["hyp"]; rho = hyp["rotation_annuelle"]; b = hyp["croissance_loyers_contrat"]
344 + vb = valeur_scenario(d, mode="esperance")["valeur"]
345 + v_rp = valeur_scenario(d, mode="esperance", rotation=rho + h_rot)["valeur"]; v_rm = valeur_scenario(d, mode="esperance", rotation=rho - h_rot)["valeur"]
346 + v_bp = valeur_scenario(d, mode="esperance", g_c=b + h_base)["valeur"]; v_bm = valeur_scenario(d, mode="esperance", g_c=b - h_base)["valeur"]
347 + return {"valeur_esperance": vb, "dV_par_point_rotation": (v_rp - v_rm) / (2 * h_rot * 100), "dV_par_point_base": (v_bp - v_bm) / (2 * h_base * 100)}
348 +
349 +# ============================================================================ 5. Exercices
350 +
351 +def exercices() -> dict:
352 + out = {}
353 + # Ex. 1 : pourcentage de base avec des IPC fictifs (A=161,0 ; B=157,2 ; C=153,6 ; D=147,0)
354 + A, B, C, D = 161.0, 157.2, 153.6, 147.0
355 + out["ex1"] = {"v1": (A - B) / B, "v2": (B - C) / C, "v3": (C - D) / D, "pct": ((A - B) / B + (B - C) / C + (C - D) / D) / 3}
356 + # Ex. 2 : taxes municipales 24 000 → 26 400 (+10 %), base 3,1 %, loyer 1 100 sur revenus 180 000
357 + part = part_attribuable(1100, 180_000)
358 + out["ex2"] = {"part": part, "variation": 0.10, "net": taux_net(0.10, 0.031), "montant": taux_net(0.10, 0.031) * 24_000 / 12 * part,
359 + "montant_si_6pct": taux_net(0.06, 0.031) * 24_000 / 12 * part, "montant_si_2pct": taux_net(0.02, 0.031) * 24_000 / 12 * part}
360 + # Ex. 3 : fixation complète logement 201 du cas 3 avec données du dossier (variante A) et travaux 40 000 $ (fenêtres)
361 + ex = exemple_fixation(201)
362 + im = ex["immeuble"]
363 + f3 = fixation_nouvelle(ex["loyer"], 0.031, ex["part"], (im["taxes_mun"][2025], im["taxes_mun"][2025] * 1.05), (im["taxes_scol"][2024], im["taxes_scol"][2025]),
364 + (im["assurances"][2024], im["assurances"][2025]), immobilisations=40_000.0)
365 + out["ex3"] = {"loyer": ex["loyer"], "marche": ex["marche"], "part": ex["part"], "fixation": f3}
366 + # Ex. 4 : forme fermée, q0 = 0,70, ρ = 0,10, g_c = 0,031, g_m = 0,035
367 + out["ex4"] = {"q": ratio_espere(0.70, 0.10, 0.031, 0.035, np.array([0, 1, 5, 10])), "q_inf": ratio_limite(0.10, 0.031, 0.035),
368 + "q_inf_20": ratio_limite(0.20, 0.031, 0.035), "survie_10": 0.9 ** 10, "temps_moyen": 1 / 0.10}
369 + # Ex. 5 : scénarios par script
370 + d = cas3.tout()
371 + out["ex5"] = {"base": d["dcf"]["valeur"], "rot10_base25": valeur_scenario(d, rotation=0.10, g_c=0.025)["valeur"],
372 + "rot10_base25_phi05": valeur_scenario(d, rotation=0.10, g_c=0.025, fraction_rattrapage=0.5)["valeur"],
373 + "rot10_base25_recup": valeur_scenario(d, rotation=0.10, g_c=0.025, recuperation_capex=True)["valeur"]}
374 + # Ex. 6 : éviction — indemnité (art. 1965) pour le logement 103 (bail depuis 2008, 18 ans → plafond 24 mois) vs gain annuel au marché
375 + loyer103 = 728.0; marche103 = 1047.0
376 + indemnite = min(max(18, 3), 24) * loyer103
377 + out["ex6"] = {"mois": min(max(18, 3), 24), "indemnite": indemnite, "gain_annuel": (marche103 - loyer103) * 12, "demenagement_hyp": 1_500.0,
378 + "recuperation_ans": (indemnite + 1_500.0) / ((marche103 - loyer103) * 12)}
379 + # Ex. 7 : dérivée de q_inf par rapport à ρ
380 + delta = 1.031 / 1.035
381 + out["ex7"] = {"delta": delta, "dq_drho_15": (1 - delta) / (1 - (1 - 0.15) * delta) ** 2, "dq_drho_5": (1 - delta) / (1 - (1 - 0.05) * delta) ** 2,
382 + "q_inf_5": ratio_limite(0.05, 0.031, 0.035), "q_inf_15": ratio_limite(0.15, 0.031, 0.035), "q_inf_25": ratio_limite(0.25, 0.031, 0.035)}
383 + return out
384 +
385 +# ============================================================================ impression
386 +
387 +def _d(x): return f"{x:,.0f}".replace(",", " ") + " $"
388 +def _d2(x): return f"{x:,.2f}".replace(",", " ").replace(".", ",") + " $"
389 +def _p(x, d=2): return f"{x*100:.{d}f}".replace(".", ",") + " %"
390 +
391 +def main():
392 + print("=== 1. Pourcentage de base (art. 3.1) à partir de l'IPC du Québec (StatCan 18-10-0004-01, v41691783) ===")
393 + for an in (2025, 2026):
394 + p = pourcentage_base(an)
395 + print(f"{an} : A={p['A']:.3f} B={p['B']:.3f} C={p['C']:.3f} D={p['D']:.3f} ; variations {_p(p['var_A_B'])}, {_p(p['var_B_C'])}, {_p(p['var_C_D'])} ; "
396 + f"moyenne {_p(p['pct_base'], 3)} ; publié TAL {_p(p['pct_publie'], 1)}")
397 + h = historique_pourcentage_base()
398 + print(h.to_string(index=False, float_format=lambda x: f"{x:.4f}"))
399 +
400 + print("\n=== 2. Exemple de fixation 2026, logement 103 du cas 3 ===")
401 + ex = exemple_fixation(103)
402 + im = ex["immeuble"]
403 + print(f"loyer {_d(ex['loyer'])}/mois ; marché {_d(ex['marche'])} ; écart {_d(ex['ecart_marche'])} ({_p(ex['ecart_marche_pct'],1)}) ; revenus immeuble {_d(ex['revenus'])} ; part {_p(ex['part'],3)}")
404 + print(f"taxes mun. 2025 {_d(im['taxes_mun'][2025])} → 2026 {_d(im['taxes_mun'][2026])} ; scolaire 2024 {_d(im['taxes_scol'][2024])} → 2025 {_d(im['taxes_scol'][2025])} ; "
405 + f"assurances 2024 {_d(im['assurances'][2024])} → 2025 {_d(im['assurances'][2025])}")
406 + for nom in ("A", "B"):
407 + f = ex[nom]
408 + print(f"-- variante {nom} : base {_d2(f['base'])} ; taxes mun. var {_p(f['var_taxes_mun'],1)} net {_p(f['net_taxes_mun'],1)} → {_d2(f['m_taxes_mun'])} ; "
409 + f"scol var {_p(f['var_taxes_scol'],1)} net {_p(f['net_taxes_scol'],1)} → {_d2(f['m_taxes_scol'])} ; assur var {_p(f['var_assurances'],1)} net {_p(f['net_assurances'],1)} → {_d2(f['m_assurances'])} ; "
410 + f"immob {_d(f['immobilisations'])} → {_d2(f['m_immobilisations'])} ; total {_d2(f['total_brut'])} ; nouveau loyer {_d(f['nouveau_loyer'])} (+{_d(f['augmentation'])}, {_p(f['taux_effectif'])})")
411 + for nom in ("ancienne_A", "ancienne_B"):
412 + a = ex[nom]
413 + print(f"-- {nom} (reconstitution) : ajustements immeuble " + ", ".join(f"{k} {_d(v)}" for k, v in a["ajustements_immeuble"].items()) +
414 + f" ; dépenses expl. {_d(a['depenses_exploitation'])} ; revenu net {_d(a['revenu_net'])} ; total {_d(a['total_immeuble'])} ; mensuel {_d2(a['mensuel_brut'])} ; nouveau loyer {_d(a['nouveau_loyer'])} ({_p(a['taux_effectif'])})")
415 + print(f"années de rattrapage du marché à 3,1 %/an si le marché était figé : {ex['annees_rattrapage_base']:.1f}")
416 +
417 + print("\n=== 3. Trajectoire espérée (forme fermée) et vérification Monte-Carlo ===")
418 + for rho in (0.05, 0.15, 0.25):
419 + q = ratio_espere(0.695, rho, 0.031, 0.035, np.array([0, 1, 2, 5, 10]))
420 + print(f"ρ={rho:.2f} : q_t(t=0,1,2,5,10) = " + ", ".join(f"{x:.4f}" for x in q) + f" ; q_inf = {ratio_limite(rho, 0.031, 0.035):.4f} ; survie 10 ans {(1-rho)**10:.3f}")
421 + sim = simulation_ratio(0.695, 0.15, 0.031, 0.035, 10)
422 + print("Monte-Carlo ρ=0,15 (t=0..10) : " + ", ".join(f"{x:.4f}" for x in sim))
423 + d = cas3.tout()
424 + tr = trajectoire_cas3_esperee(d)
425 + print(tr.to_string(index=False, float_format=lambda x: f"{x:,.4f}"))
426 +
427 + print("\n=== 4. Scénarios : grille rotation × base ===")
428 + g = grille_scenarios(d)
429 + print(g.to_string(index=False, float_format=lambda x: f"{x:,.4f}"))
430 + print("\n=== 4b. Scénarios complémentaires ===")
431 + s = scenarios_complementaires(d)
432 + print(s[["scenario", "rne_1", "rne_T1", "rbp_T1", "ecart_T1_pct", "va_ecart_marche", "valeur", "ecart_pct_base", "part_terminale"]].to_string(index=False, float_format=lambda x: f"{x:,.4f}"))
433 + print("\n=== 4c. Élasticités (mode espérance) ===")
434 + print(elasticites(d))
435 + print("\n=== 5. Exercices ===")
436 + for k, v in exercices().items():
437 + print(k, v)
438 +
439 +if __name__ == "__main__":
440 + main()
added analysis/python/cases/ch35_indices.py +369 −0
@@ -0,0 +1,369 @@
1 +#!/usr/bin/env python3
2 +"""Chapitre 35 — Indices immobiliers et ventes répétées : exemple intégral sur le cas 1 (Gatineau, DONNÉES SYNTHÉTIQUES, graine 20261031).
3 +
4 +Construit et compare, sur les 600 ventes du cas 1 (data/synthetic/cas1_gatineau_ventes.csv, 15 trimestres 2023T1-2026T3) :
5 + * l'indice de la médiane et de la moyenne géométrique des prix (biais de composition, décomposition exacte) ;
6 + * l'indice hédonique à indicatrices temporelles (traite.hedonique.indice_dummies_temps, modèle de référence traite.cas1.matrice
7 + avec les indicatrices trimestrielles du ch. 18), l'indice chaîné par périodes adjacentes et l'indice d'imputation géométrique
8 + (égal à l'indice à indicatrices quand les coefficients sont communs) ;
9 + * un jeu de REVENTES SYNTHÉTIQUES dérivé des 600 ventes (chaque vente observée est la première vente ; la revente est simulée
10 + avec la loi de génération du jeu : tendance 4,43 %/an, bruit de transaction σ = 0,085 à chaque vente, marche aléatoire
11 + idiosyncrasique δ² = 0,0004 par mois, 15 % de reventes rapides « rénovées » avec prime) ;
12 + * les estimateurs de ventes répétées de Bailey-Muth-Nourse (traite.hedonique.ventes_repetees_bmn) et de Case-Shiller pondéré en
13 + trois étapes (traite.hedonique.ventes_repetees_case_shiller), avec erreurs-types ;
14 + * l'indice SPAR (prix / valeur au rôle, Bourassa, Hoesli et Sun 2006) construit avec valeur_role_2025 ;
15 + * la comparaison à la vérité de génération, l'ajustement de conditions du marché du comparable G0404 selon l'indice,
16 + la sensibilité (millésimes et révision, filtre des reventes rapides, taille de l'échantillon, fenêtre d'estimation) et les
17 + chiffres des solutions d'exercices.
18 +Usage (racine du dépôt) : PYTHONPATH=analysis/python python3 analysis/python/cases/ch35_indices.py"""
19 +from __future__ import annotations
20 +import sys, math, pathlib
21 +import numpy as np, pandas as pd
22 +sys.path.insert(0, str(pathlib.Path(__file__).resolve().parents[1]))
23 +sys.path.insert(0, str(pathlib.Path(__file__).resolve().parent))
24 +from traite import cas1
25 +from traite import hedonique as h
26 +import ch18_hedonique as c18
27 +
28 +ROOT = pathlib.Path(__file__).resolve().parents[3]
29 +GRAINE = 20261031
30 +T = 15 # trimestres 2023T1 … 2026T3
31 +LABELS = [f"{2023 + k // 4}T{k % 4 + 1}" for k in range(T)]
32 +BETA_VRAI = c18.VRAIS["mois"] # ln(1,0443)/12 : tendance mensuelle de génération
33 +SIGMA = c18.VRAIS["sigma"] # bruit de transaction de génération
34 +DELTA2 = 0.0004 # variance mensuelle de la marche aléatoire idiosyncrasique des reventes (hypothèse)
35 +P_REVENTE = 0.60 # probabilité qu'une vente observée soit suivie d'une revente dans la fenêtre
36 +PART_FLIP = 0.15 # part des reventes rapides « rénovées »
37 +PRIME_FLIP = (0.12, 0.04) # prime moyenne et écart-type (log) des reventes rapides
38 +FIN_FENETRE = pd.Timestamp("2026-08-31")
39 +
40 +# ----------------------------------------------------------------------------- données et vérité
41 +def charger() -> pd.DataFrame:
42 + v = cas1.charger()["ventes"].copy()
43 + v["trim"] = c18.trimestre(v["date_vente"])
44 + v["mois"] = ((v["date_vente"] - cas1.REF_DATE).dt.days / 30.4375).astype(float)
45 + return v
46 +
47 +def mois_mi_trimestre() -> np.ndarray:
48 + mid = pd.to_datetime([f"{2023 + k // 4}-{3 * (k % 4) + 2:02d}-15" for k in range(T)])
49 + return ((mid - cas1.REF_DATE).days / 30.4375).to_numpy()
50 +
51 +def verite(beta: float = BETA_VRAI) -> np.ndarray:
52 + """Indice vrai (loi de génération) au milieu de chaque trimestre, base 2023T1 = 100."""
53 + m = mois_mi_trimestre()
54 + return 100 * np.exp(beta * (m - m[0]))
55 +
56 +# ----------------------------------------------------------------------------- médiane, moyenne, composition
57 +def indices_naifs(v: pd.DataFrame, td: dict) -> pd.DataFrame:
58 + """Indice de la médiane, de la moyenne arithmétique, de la moyenne géométrique et du prix médian au m² par trimestre ;
59 + décomposition exacte de la moyenne géométrique : mean(ln P)_t = α + x̄_t'β + τ_t (les résidus MCO somment à zéro par trimestre)."""
60 + g = v.groupby("trim")
61 + out = pd.DataFrame({"trimestre": LABELS, "n": g.size().reindex(range(T)).to_numpy(),
62 + "mediane": g["prix_vente"].median().to_numpy(), "moyenne": g["prix_vente"].mean().to_numpy(),
63 + "moy_geo": np.exp(g["prix_vente"].apply(lambda s: np.log(s).mean())).to_numpy(),
64 + "mediane_m2": g.apply(lambda d: (d["prix_vente"] / d["superficie_habitable_m2"]).median()).to_numpy(),
65 + "superficie_mediane": g["superficie_habitable_m2"].median().to_numpy(),
66 + "part_aylmer": g["secteur"].apply(lambda s: (s == "Aylmer").mean()).to_numpy(),
67 + "part_hull": g["secteur"].apply(lambda s: (s == "Hull").mean()).to_numpy(),
68 + "age_median": g["age"].median().to_numpy()})
69 + for c in ("mediane", "moyenne", "moy_geo", "mediane_m2"):
70 + out[f"indice_{c}"] = 100 * out[c] / out[c].iloc[0]
71 + # composition : x̄_t'β̂ (hors indicatrices temporelles) relativement à 2023T1
72 + X = td["X"]; beta = pd.Series(td["mco"]["beta"], index=X.columns)
73 + cols = [c for c in X.columns if not c.startswith("trim_")]
74 + xb = (X[cols].values @ beta[cols].values)
75 + comp = pd.Series(xb, index=v.index).groupby(v["trim"]).mean().to_numpy()
76 + out["composition_log"] = comp - comp[0]
77 + out["indice_composition"] = 100 * np.exp(out["composition_log"])
78 + out["indice_moy_geo_corrigee"] = 100 * np.exp(np.log(out["moy_geo"] / out["moy_geo"].iloc[0]) - out["composition_log"])
79 + return out
80 +
81 +# ----------------------------------------------------------------------------- hédonique : indicatrices, chaînage, imputation
82 +def indice_td(v: pd.DataFrame) -> dict:
83 + """Indice hédonique à indicatrices trimestrielles (2023T1 = référence) sur le modèle de référence du cas 1."""
84 + X = c18.matrice_trimestres(v); y = np.log(v["prix_vente"].astype(float)).values
85 + r = h.mco(X.values, y)
86 + cols = [c for c in X.columns if c.startswith("trim_")]
87 + ix = [list(X.columns).index(c) for c in cols]
88 + tau = r["beta"][ix]; se = r["se_hc1"][ix]
89 + idx = h.indice_dummies_temps(tau)
90 + return {"X": X, "y": y, "mco": r, "tau": np.concatenate([[0.0], tau]), "se": np.concatenate([[0.0], se]), "indice": idx,
91 + "bas": np.concatenate([[100.0], 100 * np.exp(tau - 1.96 * se)]), "haut": np.concatenate([[100.0], 100 * np.exp(tau + 1.96 * se)]),
92 + "sigma": float(np.sqrt(r["sigma2"])), "r2_ajuste": r["r2_ajuste"], "k": r["k"]}
93 +
94 +def indice_chaine(v: pd.DataFrame) -> dict:
95 + """Indice hédonique chaîné : pour chaque paire de trimestres adjacents (t−1, t), MCO sur les ventes des deux trimestres avec une
96 + indicatrice de t (coefficients des caractéristiques propres à la paire) ; I_t = I_{t−1} exp(δ_t)."""
97 + delta = np.zeros(T); se = np.zeros(T)
98 + for t in range(1, T):
99 + sub = v[v["trim"].isin([t - 1, t])]
100 + X = cas1.matrice(sub).drop(columns=["mois"])
101 + X = X[[c for c in X.columns if X[c].std() > 0 or c == "const"]]
102 + X["d_t"] = (sub["trim"] == t).astype(float)
103 + r = h.mco(X.values, np.log(sub["prix_vente"].astype(float)).values)
104 + j = list(X.columns).index("d_t"); delta[t] = r["beta"][j]; se[t] = r["se_hc1"][j]
105 + ln_idx = np.cumsum(delta)
106 + return {"delta": delta, "se_delta": se, "indice": 100 * np.exp(ln_idx), "se_cumul": np.sqrt(np.cumsum(se ** 2))}
107 +
108 +def indice_imputation(v: pd.DataFrame, td: dict) -> dict:
109 + """Imputation géométrique (Laspeyres : panier de la période de base ; Paasche : panier courant ; Fisher = moyenne géométrique)
110 + avec le modèle à indicatrices (coefficients communs) : chaque vente de la période 0 (resp. t) est « revendue » à la période t (resp. 0)
111 + au prix prédit. Avec des coefficients communs, les trois indices coïncident avec exp(τ_t) ; la fonction le vérifie numériquement."""
112 + X = td["X"]; beta = td["mco"]["beta"]; cols = list(X.columns)
113 + tcols = [c for c in cols if c.startswith("trim_")]
114 + def predire_a(Xsub: pd.DataFrame, t: int) -> np.ndarray:
115 + Z = Xsub.copy()
116 + for c in tcols: Z[c] = 0.0
117 + if t > 0: Z[f"trim_{t}"] = 1.0
118 + return Z.values @ beta
119 + lasp = np.zeros(T); paas = np.zeros(T)
120 + X0 = X[v["trim"] == 0]
121 + for t in range(T):
122 + Xt = X[v["trim"] == t]
123 + lasp[t] = np.mean(predire_a(X0, t) - predire_a(X0, 0))
124 + paas[t] = np.mean(predire_a(Xt, t) - predire_a(Xt, 0))
125 + fisher = 0.5 * (lasp + paas)
126 + return {"laspeyres": 100 * np.exp(lasp), "paasche": 100 * np.exp(paas), "fisher": 100 * np.exp(fisher),
127 + "ecart_max_vs_td": float(np.max(np.abs(100 * np.exp(fisher) - td["indice"])))}
128 +
129 +# ----------------------------------------------------------------------------- reventes synthétiques
130 +def generer_reventes(v: pd.DataFrame, graine: int = GRAINE, p_revente: float = P_REVENTE, part_flip: float = PART_FLIP,
131 + delta2: float = DELTA2, beta: float = BETA_VRAI, sigma: float = SIGMA) -> pd.DataFrame:
132 + """Paires de reventes synthétiques. Pour chaque vente observée (première vente, prix P1 à la date t1), avec probabilité p_revente :
133 + intervalle Δ (mois) tiré d'une loi déclarée (reventes « rapides » : 3 à 9 mois, prime log N(0,08 ; 0,03) ; autres : 6 + Exp(16), ≥ 6) ;
134 + ln P2 = ln P1 + β Δ + prime·flip + w + ε2 − ε1', avec ε1', ε2 ~ N(0, σ²) (bruit de transaction à chaque vente) et w ~ N(0, δ² Δ).
135 + Seules les reventes datées au plus tard le 31 août 2026 sont conservées ; prix arrondis à 500 $."""
136 + rng = np.random.default_rng(graine)
137 + rows = []
138 + for _, s in v.iterrows():
139 + if rng.random() > p_revente:
140 + continue
141 + flip = rng.random() < part_flip
142 + if flip:
143 + gap = int(rng.integers(3, 10)); prime = rng.normal(*PRIME_FLIP)
144 + else:
145 + gap = int(max(6, round(6 + rng.exponential(16.0)))); prime = 0.0
146 + d2 = s["date_vente"] + pd.DateOffset(months=gap)
147 + if d2 > FIN_FENETRE:
148 + continue
149 + w = rng.normal(0.0, math.sqrt(delta2 * gap)); e1 = rng.normal(0.0, sigma); e2 = rng.normal(0.0, sigma)
150 + lnp2 = math.log(s["prix_vente"]) + beta * gap + prime + w + e2 - e1
151 + p2 = 500 * round(math.exp(lnp2) / 500)
152 + rows.append({"id_vente": s["id_vente"], "secteur": s["secteur"], "date_1": s["date_vente"], "date_2": d2, "prix_1": float(s["prix_vente"]),
153 + "prix_2": float(p2), "t_1": int(s["trim"]), "t_2": int(c18.trimestre(pd.Series([d2]))[0]), "intervalle_mois": gap,
154 + "flip": bool(flip), "prime": prime})
155 + df = pd.DataFrame(rows)
156 + df["ln_ratio"] = np.log(df["prix_2"] / df["prix_1"])
157 + df["id_paire"] = [f"R{i + 1:04d}" for i in range(len(df))]
158 + return df
159 +
160 +def resume_reventes(p: pd.DataFrame) -> dict:
161 + return {"n_paires": len(p), "n_flips": int(p["flip"].sum()), "intervalle_median": float(p["intervalle_mois"].median()),
162 + "intervalle_moyen": float(p["intervalle_mois"].mean()), "intervalle_min": int(p["intervalle_mois"].min()), "intervalle_max": int(p["intervalle_mois"].max()),
163 + "part_moins_12": float((p["intervalle_mois"] < 12).mean()), "variation_mediane": float(np.exp(p["ln_ratio"].median()) - 1),
164 + "variation_moyenne_log": float(p["ln_ratio"].mean()), "sd_ln_ratio": float(p["ln_ratio"].std(ddof=1)),
165 + "paires_par_annee_1": p["date_1"].dt.year.value_counts().sort_index().to_dict(),
166 + "paires_par_annee_2": p["date_2"].dt.year.value_counts().sort_index().to_dict(),
167 + "meme_trimestre": int((p["t_1"] == p["t_2"]).sum())}
168 +
169 +# ----------------------------------------------------------------------------- ventes répétées : BMN et Case–Shiller avec erreurs-types
170 +def _se_mco_sans_constante(D: np.ndarray, y: np.ndarray, beta: np.ndarray, w: np.ndarray | None = None) -> np.ndarray:
171 + if w is None: w = np.ones(len(y))
172 + Dw = D * w[:, None]; yw = y * w
173 + res = yw - Dw @ beta; n, k = D.shape
174 + s2 = res @ res / (n - k)
175 + return np.sqrt(np.diag(s2 * np.linalg.pinv(Dw.T @ Dw)))
176 +
177 +def bmn(p: pd.DataFrame) -> dict:
178 + r = h.ventes_repetees_bmn(p["prix_1"], p["prix_2"], p["t_1"], p["t_2"], T)
179 + se = _se_mco_sans_constante(r["D"], r["y"], r["beta"])
180 + r["se"] = np.concatenate([[0.0], se]); r["bas"] = 100 * np.exp(np.concatenate([[0.0], r["beta"]]) - 1.96 * r["se"])
181 + r["haut"] = 100 * np.exp(np.concatenate([[0.0], r["beta"]]) + 1.96 * r["se"]); r["sigma"] = float(np.sqrt(r["residus"] @ r["residus"] / (len(r["y"]) - T + 1)))
182 + # vérification de la forme fermée : β = (D'D)^{-1} D'y
183 + r["beta_ferme"] = np.linalg.solve(r["D"].T @ r["D"], r["D"].T @ r["y"])
184 + return r
185 +
186 +def case_shiller(p: pd.DataFrame) -> dict:
187 + r = h.ventes_repetees_case_shiller(p["prix_1"], p["prix_2"], p["t_1"], p["t_2"], T)
188 + e1 = h.ventes_repetees_bmn(p["prix_1"], p["prix_2"], p["t_1"], p["t_2"], T)
189 + gap_t = (p["t_2"] - p["t_1"]).to_numpy(float)
190 + se = _se_mco_sans_constante(e1["D"], e1["y"], r["beta"], r["poids"])
191 + r["se"] = np.concatenate([[0.0], se]); r["bas"] = 100 * np.exp(np.concatenate([[0.0], r["beta"]]) - 1.96 * r["se"])
192 + r["haut"] = 100 * np.exp(np.concatenate([[0.0], r["beta"]]) + 1.96 * r["se"])
193 + r["gap_trimestres"] = gap_t; r["residus_bmn2"] = e1["residus"] ** 2
194 + a = r["coef_variance"] # var(ε) = a0 + a1·(intervalle en trimestres)
195 + r["variance_a_gap"] = {g: float(a[0] + a[1] * g) for g in (1, 2, 4, 8, 12)}
196 + r["poids_relatif_a_gap"] = {g: float(math.sqrt((a[0] + a[1] * 1) / (a[0] + a[1] * g))) for g in (1, 2, 4, 8, 12)}
197 + # pente théorique : var = 2σ² + δ²·(3 mois par trimestre)·gap, en unités de trimestre
198 + r["variance_theorique"] = (2 * SIGMA ** 2, DELTA2 * 3)
199 + return r
200 +
201 +# ----------------------------------------------------------------------------- SPAR
202 +def spar(v: pd.DataFrame) -> dict:
203 + """Indice SPAR (sale price appraisal ratio) : arithmétique [(ΣP_t/ΣA_t)/(ΣP_0/ΣA_0)] et géométrique [exp(mean ln(P/A)_t − mean ln(P/A)_0)],
204 + A = valeur au rôle 2025 (date de référence 2024-07-01, fictive)."""
205 + g = v.groupby("trim")
206 + ratio_arith = (g["prix_vente"].sum() / g["valeur_role_2025"].sum()).to_numpy()
207 + lr = np.log(v["prix_vente"] / v["valeur_role_2025"])
208 + ratio_geo = lr.groupby(v["trim"]).mean().to_numpy(); sd_geo = lr.groupby(v["trim"]).std(ddof=1).to_numpy(); n = g.size().to_numpy()
209 + se_geo = np.sqrt(sd_geo ** 2 / n + sd_geo[0] ** 2 / n[0]); se_geo[0] = 0.0
210 + return {"ratio_arith": ratio_arith, "ratio_geo": np.exp(ratio_geo), "indice_arith": 100 * ratio_arith / ratio_arith[0],
211 + "indice_geo": 100 * np.exp(ratio_geo - ratio_geo[0]), "se_log": se_geo,
212 + "bas": 100 * np.exp(ratio_geo - ratio_geo[0] - 1.96 * se_geo), "haut": 100 * np.exp(ratio_geo - ratio_geo[0] + 1.96 * se_geo),
213 + "sd_ln_ratio": float(lr.std(ddof=1)), "ratio_median_global": float(np.median(v["prix_vente"] / v["valeur_role_2025"]))}
214 +
215 +# ----------------------------------------------------------------------------- comparaison
216 +def precision(idx: np.ndarray, vrai: np.ndarray) -> dict:
217 + e = 100 * (idx / vrai - 1)
218 + return {"EAM_pct": float(np.mean(np.abs(e[1:]))), "RMSE_pct": float(np.sqrt(np.mean(e[1:] ** 2))), "max_abs_pct": float(np.max(np.abs(e))),
219 + "croissance_totale_pct": float(idx[-1] - 100), "ecart_final_pct": float(e[-1]),
220 + "variations_trim_sd_pct": float(np.std(np.diff(np.log(idx))) * 100)}
221 +
222 +def tableau_indices(vrai, naifs, td, ch, b, cs, sp) -> pd.DataFrame:
223 + return pd.DataFrame({"trimestre": LABELS, "n_ventes": naifs["n"], "vrai": vrai, "mediane": naifs["indice_mediane"].to_numpy(),
224 + "moy_geo": naifs["indice_moy_geo"].to_numpy(), "hedonique": td["indice"], "hed_bas": td["bas"], "hed_haut": td["haut"],
225 + "chaine": ch["indice"], "bmn": b["indice"], "bmn_bas": b["bas"], "bmn_haut": b["haut"], "cs": cs["indice"],
226 + "cs_bas": cs["bas"], "cs_haut": cs["haut"], "spar": sp["indice_geo"], "spar_arith": sp["indice_arith"]})
227 +
228 +def ajustement_comparable(v: pd.DataFrame, tab: pd.DataFrame, ref: dict, id_vente: str = "G0404", date_eval: str = "2026-09-01") -> dict:
229 + """Ajustement de conditions du marché d'un comparable du ch. 16 selon l'indice retenu : A = P × (I_eval / I_vente − 1)."""
230 + c = v[v["id_vente"] == id_vente].iloc[0]; P = float(c["prix_vente"]); t1 = int(c["trim"]); t2 = T - 1
231 + dm = (pd.Timestamp(date_eval) - c["date_vente"]).days / 30.4375
232 + bm = ref["beta"][list(ref["X"].columns).index("mois")]
233 + out = {"prix": P, "trimestre_vente": LABELS[t1], "mois_ecoules": float(dm), "tendance": P * (math.exp(bm * dm) - 1)}
234 + for k in ("mediane", "hedonique", "chaine", "bmn", "cs", "spar"):
235 + out[k] = P * (tab[k].iloc[t2] / tab[k].iloc[t1] - 1)
236 + out["vrai"] = P * (math.exp(BETA_VRAI * dm) - 1)
237 + return out
238 +
239 +# ----------------------------------------------------------------------------- sensibilité
240 +def sensibilite(v: pd.DataFrame, p: pd.DataFrame, vrai: np.ndarray, td: dict, graine: int = GRAINE) -> dict:
241 + out = {}
242 + # (a) millésimes : révision des indices de ventes répétées quand de nouvelles paires s'ajoutent
243 + mill = {}
244 + for fin, lab in ((9, "2025T2"), (11, "2025T4"), (T - 1, "2026T3")):
245 + sub = p[p["t_2"] <= fin]
246 + rb = h.ventes_repetees_bmn(sub["prix_1"], sub["prix_2"], sub["t_1"], sub["t_2"], fin + 1)
247 + rc = h.ventes_repetees_case_shiller(sub["prix_1"], sub["prix_2"], sub["t_1"], sub["t_2"], fin + 1)
248 + mill[lab] = {"n": len(sub), "bmn": rb["indice"], "cs": rc["indice"]}
249 + out["millesimes"] = mill
250 + rev = {}
251 + for t, lab in ((9, "2025T2"), (7, "2024T4"), (5, "2024T2")):
252 + rev[lab] = {"bmn_v2025T2": float(mill["2025T2"]["bmn"][t]), "bmn_v2025T4": float(mill["2025T4"]["bmn"][t]), "bmn_v2026T3": float(mill["2026T3"]["bmn"][t]),
253 + "cs_v2025T2": float(mill["2025T2"]["cs"][t]), "cs_v2026T3": float(mill["2026T3"]["cs"][t])}
254 + out["revisions"] = rev
255 + # (b) filtre des reventes rapides
256 + filt = {}
257 + for seuil, lab in ((0, "toutes"), (12, "intervalle ≥ 12 mois"), (18, "intervalle ≥ 18 mois")):
258 + sub = p[p["intervalle_mois"] >= seuil] if seuil else p
259 + rb = h.ventes_repetees_bmn(sub["prix_1"], sub["prix_2"], sub["t_1"], sub["t_2"], T)
260 + rc = h.ventes_repetees_case_shiller(sub["prix_1"], sub["prix_2"], sub["t_1"], sub["t_2"], T)
261 + filt[lab] = {"n": len(sub), "n_flips": int(sub["flip"].sum()), "bmn": precision(rb["indice"], vrai), "cs": precision(rc["indice"], vrai),
262 + "bmn_indice": rb["indice"], "cs_indice": rc["indice"]}
263 + sans_flip = p[~p["flip"]]
264 + rb = h.ventes_repetees_bmn(sans_flip["prix_1"], sans_flip["prix_2"], sans_flip["t_1"], sans_flip["t_2"], T)
265 + filt["sans les reventes rénovées (oracle)"] = {"n": len(sans_flip), "n_flips": 0, "bmn": precision(rb["indice"], vrai), "bmn_indice": rb["indice"]}
266 + out["filtre"] = filt
267 + # (c) taille de l'échantillon : sous-échantillons aléatoires de paires, erreur-type du dernier trimestre et EAM (moyenne de 200 tirages)
268 + rng = np.random.default_rng(graine); taille = {}
269 + for frac in (1.0, 0.5, 0.25):
270 + eam, fin = [], []
271 + for _ in range(200 if frac < 1 else 1):
272 + sub = p.sample(frac=frac, random_state=int(rng.integers(1e9))) if frac < 1 else p
273 + try:
274 + rb = h.ventes_repetees_bmn(sub["prix_1"], sub["prix_2"], sub["t_1"], sub["t_2"], T)
275 + except Exception:
276 + continue
277 + eam.append(precision(rb["indice"], vrai)["EAM_pct"]); fin.append(rb["indice"][-1])
278 + taille[frac] = {"n": int(round(frac * len(p))), "EAM_moyen_pct": float(np.mean(eam)), "sd_indice_final": float(np.std(fin)) if len(fin) > 1 else float("nan")}
279 + out["taille"] = taille
280 + # (d) fenêtre d'estimation de l'indice hédonique : huit derniers trimestres seulement, rechaîné sur 2024T4
281 + sub = v[v["trim"] >= 7]
282 + X = cas1.matrice(sub).drop(columns=["mois"]); X = X[[c for c in X.columns if X[c].std() > 0 or c == "const"]]
283 + for k in range(8, T): X[f"trim_{k}"] = (sub["trim"] == k).astype(float)
284 + r = h.mco(X.values, np.log(sub["prix_vente"].astype(float)).values)
285 + tau = np.array([r["beta"][list(X.columns).index(f"trim_{k}")] for k in range(8, T)])
286 + idx_court = td["indice"][7] * np.exp(np.concatenate([[0.0], tau]))
287 + out["fenetre"] = {"n": len(sub), "indice_court": idx_court, "indice_long": td["indice"][7:], "ecart_max_pct": float(np.max(np.abs(idx_court / td["indice"][7:] - 1)) * 100),
288 + "sigma_court": float(np.sqrt(r["sigma2"]))}
289 + return out
290 +
291 +# ----------------------------------------------------------------------------- exercices
292 +def exercices(tab: pd.DataFrame, b: dict, cs: dict, sp: dict, ref: dict) -> dict:
293 + e = {}
294 + # 35.1 : médiane et composition (deux trimestres fictifs)
295 + e["1"] = {"mediane_t0": 450_000, "mediane_t1": 480_000, "var_mediane_pct": (480_000 / 450_000 - 1) * 100,
296 + "composition_pct": (math.exp(0.539 * math.log(150 / 140)) - 1) * 100,
297 + "variation_pure_pct": (480_000 / 450_000 / math.exp(0.539 * math.log(150 / 140)) - 1) * 100}
298 + # 35.2 : indice à indicatrices
299 + tau = np.array([0.0, 0.031, 0.058, 0.102]); e["2"] = {"indice": 100 * np.exp(tau), "variation_t3_t2_pct": (math.exp(0.102 - 0.058) - 1) * 100,
300 + "ic_t3": (100 * math.exp(0.102 - 1.96 * 0.021), 100 * math.exp(0.102 + 1.96 * 0.021))}
301 + # 35.3 : BMN à la main sur quatre paires et trois périodes
302 + D = np.array([[1.0, 0.0], [-1.0, 1.0], [0.0, 1.0], [1.0, 0.0]]); y = np.log(np.array([462_000 / 430_000, 512_000 / 470_000, 540_000 / 495_000, 455_000 / 436_000]))
303 + beta = np.linalg.solve(D.T @ D, D.T @ y)
304 + e["3"] = {"y": y, "DtD": D.T @ D, "Dty": D.T @ y, "beta": beta, "indice": 100 * np.exp(np.concatenate([[0.0], beta])), "residus": y - D @ beta}
305 + # 35.4 : Case–Shiller : variance et poids
306 + a0, a1 = 0.0145, 0.0012; gaps = np.array([1, 4, 8, 12])
307 + var = a0 + a1 * gaps; e["4"] = {"var": var, "poids": 1 / np.sqrt(var), "poids_rel": np.sqrt(var[0] / var), "sd_pct": 100 * np.sqrt(var)}
308 + # 35.5 : SPAR sur dix ventes fictives
309 + P = np.array([512_000, 468_500, 601_000, 445_000, 530_000, 487_500, 575_000, 455_000, 498_000, 522_000.0]); A = np.array([470_000, 452_000, 545_000, 425_000, 480_000, 468_000, 510_000, 440_000, 455_000, 471_000.0])
310 + base_ratio = 1.045; r_ar = P.sum() / A.sum(); r_geo = math.exp(np.mean(np.log(P / A)))
311 + e["5"] = {"ratio_arith": r_ar, "ratio_geo": r_geo, "spar_arith": 100 * r_ar / base_ratio, "spar_geo": 100 * r_geo / base_ratio, "ratios": P / A}
312 + # 35.6 : ajustement de G0404 selon l'indice
313 + e["6"] = {"P": 506_500.0, "I_2025T2": {k: float(tab[k].iloc[9]) for k in ("hedonique", "bmn", "cs", "mediane")},
314 + "I_2026T3": {k: float(tab[k].iloc[14]) for k in ("hedonique", "bmn", "cs", "mediane")}}
315 + e["6"]["ajustements"] = {k: 506_500.0 * (e["6"]["I_2026T3"][k] / e["6"]["I_2025T2"][k] - 1) for k in e["6"]["I_2025T2"]}
316 + bm = ref["beta"][list(ref["X"].columns).index("mois")]; e["6"]["tendance"] = 506_500.0 * (math.exp(bm * 15.2) - 1)
317 + # 35.7 : révision et biais de sélection (chiffres du script : voir sensibilité) — valeurs conceptuelles
318 + # 35.8 : script (pas de chiffre fixe)
319 + return e
320 +
321 +# ----------------------------------------------------------------------------- principal
322 +def main(verbose: bool = True) -> dict:
323 + v = charger(); ref = cas1.hedonique(v.drop(columns=["trim", "mois"]))
324 + vrai = verite()
325 + td = indice_td(v); ch = indice_chaine(v); imp = indice_imputation(v, td)
326 + naifs = indices_naifs(v, td)
327 + p = generer_reventes(v); rp = resume_reventes(p)
328 + b = bmn(p); cs = case_shiller(p); sp = spar(v)
329 + tab = tableau_indices(vrai, naifs, td, ch, b, cs, sp)
330 + prec = {k: precision(tab[k].to_numpy(), vrai) for k in ("mediane", "moy_geo", "hedonique", "chaine", "bmn", "cs", "spar", "spar_arith")}
331 + prec["hedonique"]["largeur_ic_finale_pct"] = float((td["haut"][-1] / td["bas"][-1] - 1) * 100)
332 + prec["bmn"]["largeur_ic_finale_pct"] = float((b["haut"][-1] / b["bas"][-1] - 1) * 100)
333 + prec["cs"]["largeur_ic_finale_pct"] = float((cs["haut"][-1] / cs["bas"][-1] - 1) * 100)
334 + prec["spar"]["largeur_ic_finale_pct"] = float((sp["haut"][-1] / sp["bas"][-1] - 1) * 100)
335 + # corrélation entre l'écart médiane − hédonique et la composition
336 + ecart_med = np.log(naifs["indice_mediane"].to_numpy() / td["indice"]); corr_comp = float(np.corrcoef(ecart_med, naifs["composition_log"])[0, 1])
337 + corr_surf = float(np.corrcoef(ecart_med, np.log(naifs["superficie_mediane"]))[0, 1])
338 + adj = ajustement_comparable(v, tab, ref)
339 + sens = sensibilite(v, p, vrai, td)
340 + ex = exercices(tab, b, cs, sp, ref)
341 + out = {"ventes": v, "ref": ref, "vrai": vrai, "td": td, "chaine": ch, "imputation": imp, "naifs": naifs, "paires": p, "resume_paires": rp,
342 + "bmn": b, "cs": cs, "spar": sp, "tableau": tab, "precision": prec, "corr_ecart_composition": corr_comp, "corr_ecart_superficie": corr_surf,
343 + "ajustement": adj, "sensibilite": sens, "exercices": ex, "beta_vrai": BETA_VRAI}
344 + if verbose:
345 + pd.set_option("display.width", 240); pd.set_option("display.max_columns", 40)
346 + print(f"== Vérité de génération : β = {BETA_VRAI:.6f}/mois ({(math.exp(12 * BETA_VRAI) - 1) * 100:.2f} %/an) ; indice final {vrai[-1]:.2f}")
347 + print("== Médiane, moyennes et composition :\n", naifs.round(3).to_string(index=False))
348 + print(f" corrélation(ln médiane/hédonique, composition) = {corr_comp:.3f} ; avec ln superficie médiane = {corr_surf:.3f}")
349 + print(f"== Hédonique à indicatrices : k = {td['k']}, σ = {td['sigma']:.4f}, R² aj. = {td['r2_ajuste']:.4f}")
350 + print(f"== Imputation géométrique : écart max à l'indice à indicatrices = {imp['ecart_max_vs_td']:.2e}")
351 + print("== Reventes synthétiques :", {k: (round(x, 4) if isinstance(x, float) else x) for k, x in rp.items()})
352 + print(f"== BMN : σ = {b['sigma']:.4f} ; écart max forme fermée vs lstsq = {np.max(np.abs(b['beta'] - b['beta_ferme'])):.1e}")
353 + print(f"== Case–Shiller : var(ε) = {cs['coef_variance'][0]:.5f} + {cs['coef_variance'][1]:.5f} × intervalle (trimestres) ; théorie {cs['variance_theorique']}")
354 + print(" poids relatifs :", {k: round(x, 3) for k, x in cs["poids_relatif_a_gap"].items()})
355 + print(f"== SPAR : sd ln(P/A) = {sp['sd_ln_ratio']:.4f} ; ratio médian global {sp['ratio_median_global']:.4f}")
356 + print("== Indices par trimestre :\n", tab.round(2).to_string(index=False))
357 + print("== Précision vs vérité :")
358 + for k, d in prec.items(): print(f" {k:11s}", {kk: round(x, 2) for kk, x in d.items()})
359 + print("== Ajustement de G0404 :", {k: (round(x, 1) if isinstance(x, float) else x) for k, x in adj.items()})
360 + print("== Sensibilité — révisions :", sens["revisions"])
361 + for lab, d in sens["filtre"].items(): print(f" filtre {lab:36s} n={d['n']:3d} flips={d['n_flips']:2d} BMN EAM={d['bmn']['EAM_pct']:.2f} final={d['bmn']['croissance_totale_pct']:.2f}", f"CS EAM={d['cs']['EAM_pct']:.2f}" if "cs" in d else "")
362 + print(" taille :", sens["taille"])
363 + print(" fenêtre hédonique :", {k: (np.round(x, 2) if isinstance(x, np.ndarray) else round(x, 3)) for k, x in sens["fenetre"].items()})
364 + print("== Exercices :")
365 + for k, e in ex.items(): print(" ", k, {kk: (np.round(vv, 4) if isinstance(vv, (np.ndarray, float)) else vv) for kk, vv in e.items()})
366 + return out
367 +
368 +if __name__ == "__main__":
369 + main()
added analysis/python/figures/ch20_figures.py +145 −0
@@ -0,0 +1,145 @@
1 +#!/usr/bin/env python3
2 +"""Chapitre 20 — figures : densités bootstrap et énumération exhaustive ; résidus des candidats et élargissement ;
3 +fourchettes par méthode ; sensibilité de la conclusion aux poids et variance de la combinaison selon la corrélation.
4 +Usage (racine du dépôt) : PYTHONPATH=analysis/python python3 analysis/python/figures/ch20_figures.py"""
5 +from __future__ import annotations
6 +import sys, pathlib
7 +import numpy as np
8 +import matplotlib.pyplot as plt
9 +from matplotlib.ticker import FuncFormatter
10 +from scipy.stats import gaussian_kde
11 +sys.path.insert(0, str(pathlib.Path(__file__).resolve().parents[1]))
12 +sys.path.insert(0, str(pathlib.Path(__file__).resolve().parents[1] / "cases"))
13 +from traite.style import setup, palette, savefig, LARGEUR, fmt_dollars
14 +import ch20_reconciliation as c20
15 +
16 +def fmt_k(x, pos=None):
17 + return f"{x/1000:,.0f}".replace(",", " ") + " k$"
18 +
19 +
20 +def 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")
42 +
43 +
44 +def 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) * 100
50 + ret = a["retenu"].values
51 + 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")
68 +
69 +
70 +def 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.6449
74 + 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")
101 +
102 +
103 +def 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)) / 1000
126 + 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")
135 +
136 +
137 +def 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/")
142 +
143 +
144 +if __name__ == "__main__":
145 + main()
added analysis/python/figures/ch34_calcul_tal.py +36 −0
@@ -0,0 +1,36 @@
1 +#!/usr/bin/env python3
2 +"""Ch. 34, figure « exemple de calcul TAL » : composantes de l'augmentation mensuelle du logement 103 du cas 3 (DONNÉES
3 +SYNTHÉTIQUES) selon la méthode en vigueur depuis le 1er janvier 2026 (pourcentages 2026), variante A (données du dossier) et
4 +variante B (illustrative : taxes +8 %, assurances +10 %, immobilisations 85 000 $), avec l'écart au loyer de marché en repère.
5 +Sortie : ch34-fig-calcul-tal.pdf."""
6 +from __future__ import annotations
7 +import sys, pathlib
8 +sys.path.insert(0, str(pathlib.Path(__file__).resolve().parents[1]))
9 +import numpy as np
10 +import matplotlib.pyplot as plt
11 +from traite.style import setup, palette, savefig, LARGEUR
12 +from cases import ch34_tal as c
13 +
14 +def main():
15 + setup()
16 + ex = c.exemple_fixation(103)
17 + comp = [("Pourcentage de base\n(3,1 % × 728 $)", "base"), ("Taxes municipales\n(au-delà de 3,1 %)", "m_taxes_mun"), ("Taxes scolaires\n(au-delà de 3,1 %)", "m_taxes_scol"),
18 + ("Assurances\n(au-delà de 3,1 %)", "m_assurances"), ("Dépenses d'immobilisation\n(5 % × part)", "m_immobilisations"), ("Total avant arrondi", "total_brut")]
19 + fig, ax = plt.subplots(figsize=(LARGEUR, 3.2))
20 + ys = np.arange(len(comp)); h = 0.36
21 + for k, (nom, col, lab) in enumerate((("A", palette.cat[0], "variante A — données du dossier (taxes et assurances en baisse, aucune immobilisation)"),
22 + ("B", palette.cat[1], "variante B — illustrative (taxes +8 %, assurances +10 %, travaux de 85 000 $)"))):
23 + vals = [ex[nom][key] for _, key in comp]
24 + bars = ax.barh(ys + (0.5 - k) * h, vals, h, color=col, label=lab)
25 + for bar, v in zip(bars, vals):
26 + ax.text(v + 0.5, bar.get_y() + bar.get_height() / 2, f"{v:.2f} $".replace(".", ","), va="center", fontsize=7, color=palette.ink)
27 + ax.set_yticks(ys); ax.set_yticklabels([n for n, _ in comp], fontsize=7.5); ax.invert_yaxis()
28 + ax.set_xlabel("Augmentation mensuelle (dollars par mois)")
29 + ax.set_xlim(0, 62)
30 + ax.text(61, 0.05, f"Repère : écart au loyer de marché du logement 103 = {ex['ecart_marche']:.0f} $/mois ({ex['ecart_marche_pct']*100:.1f} % du loyer), hors échelle.".replace(".", ","),
31 + ha="right", va="top", fontsize=7, color=palette.ink2)
32 + ax.legend(loc="lower right", fontsize=7)
33 + savefig(fig, "ch34/ch34-fig-calcul-tal.pdf")
34 +
35 +if __name__ == "__main__":
36 + main()
added analysis/python/figures/ch34_ipc_base.py +31 −0
@@ -0,0 +1,31 @@
1 +#!/usr/bin/env python3
2 +"""Ch. 34, figure « IPC et pourcentage de base » (DONNÉES RÉELLES : Statistique Canada, tableau 18-10-0004-01, IPC du Québec,
3 +Ensemble, 2002 = 100, extrait dans data/processed/ipc_quebec_mensuel.csv) : inflation de la période de référence et pourcentage
4 +de base que donne la formule de l'art. 3.1 du Règlement (moyenne des trois dernières variations annuelles), 2006-2026 ; les
5 +pourcentages publiés par le TAL pour 2025 (4,5 %) et 2026 (3,1 %) sont marqués. Avant 2026, la formule est contrefactuelle.
6 +Sortie : ch34-fig-ipc-base.pdf."""
7 +from __future__ import annotations
8 +import sys, pathlib
9 +sys.path.insert(0, str(pathlib.Path(__file__).resolve().parents[1]))
10 +import matplotlib.pyplot as plt
11 +from matplotlib.ticker import FuncFormatter
12 +from traite.style import setup, palette, savefig, LARGEUR
13 +from cases import ch34_tal as c
14 +
15 +def main():
16 + setup()
17 + h = c.historique_pourcentage_base(2006, 2026)
18 + fig, ax = plt.subplots(figsize=(LARGEUR, 3.0))
19 + ax.bar(h["annee"], h["inflation_ref"] * 100, color=palette.cat[2], alpha=0.75, width=0.7, label="inflation de la période de référence (IPC Québec, moyenne annuelle)")
20 + ax.plot(h["annee"], h["pct_formule"] * 100, color=palette.cat[0], marker="o", ms=3.5, label="pourcentage de base selon la formule de l'art. 3.1 (contrefactuel avant 2026)")
21 + pub = h.dropna(subset=["pct_publie"])
22 + ax.plot(pub["annee"], pub["pct_publie"] * 100, ls="none", marker="D", ms=6, mfc="white", mec=palette.cat[1], mew=1.4, label="pourcentage de base publié par le TAL (2025 : 4,5 % ; 2026 : 3,1 %)")
23 + ax.axvspan(2025.5, 2026.5, color=palette.grid, alpha=0.6, lw=0)
24 + ax.text(2026, 6.6, "nouvelle\nméthode", ha="center", fontsize=7, color=palette.ink2)
25 + ax.set_xlabel("Année de fixation (période de référence = année civile précédente)"); ax.set_ylabel("Pourcentage (%)")
26 + ax.yaxis.set_major_formatter(FuncFormatter(lambda v, q: f"{v:.0f}")); ax.set_xticks(range(2006, 2027, 2)); ax.set_ylim(0, 7.2)
27 + ax.legend(fontsize=7, loc="upper left")
28 + savefig(fig, "ch34/ch34-fig-ipc-base.pdf")
29 +
30 +if __name__ == "__main__":
31 + main()
added analysis/python/figures/ch34_rotation.py +38 −0
@@ -0,0 +1,38 @@
1 +#!/usr/bin/env python3
2 +"""Ch. 34, figure « rotation espérée » : (gauche) rapport espéré loyer contractuel / loyer de marché d'un logement partant de
3 +q0 = 0,695 (logement 103 du cas 3, DONNÉES SYNTHÉTIQUES) pour une rotation de 5, 15 et 25 %, forme fermée (trait) et simulation
4 +Monte-Carlo (points), avec la limite q_∞ ; (droite) probabilité que le locataire initial soit encore en place après t années.
5 +Sortie : ch34-fig-rotation-esperee.pdf."""
6 +from __future__ import annotations
7 +import sys, pathlib
8 +sys.path.insert(0, str(pathlib.Path(__file__).resolve().parents[1]))
9 +import numpy as np
10 +import matplotlib.pyplot as plt
11 +from matplotlib.ticker import FuncFormatter
12 +from traite.style import setup, palette, savefig, LARGEUR
13 +from cases import ch34_tal as c
14 +
15 +def main():
16 + setup()
17 + q0, gc, gm, T = 0.695, 0.031, 0.035, 15
18 + t = np.arange(0, T + 1)
19 + fig, (a1, a2) = plt.subplots(1, 2, figsize=(LARGEUR, 3.0), gridspec_kw={"width_ratios": [1.6, 1]})
20 + for rho, col in ((0.05, palette.cat[2]), (0.15, palette.cat[0]), (0.25, palette.cat[1])):
21 + q = c.ratio_espere(q0, rho, gc, gm, t); sim = c.simulation_ratio(q0, rho, gc, gm, T, n=20_000)
22 + a1.plot(t, q, color=col, label=f"rotation {rho*100:.0f} % : forme fermée")
23 + a1.plot(t, sim, color=col, ls="none", marker="o", ms=3, mfc="white", mew=1.0)
24 + qi = c.ratio_limite(rho, gc, gm)
25 + a1.axhline(qi, color=col, lw=0.7, ls=":")
26 + a1.text(T + 0.2, qi, f"q∞ = {qi:.3f}".replace(".", ","), fontsize=7, color=palette.ink, va="center")
27 + a1.plot(t, c.ratio_sans_relocation(q0, gc, gm, t), color=palette.ink2, ls="--", lw=1.0, label="sans relocation (indexation seule)")
28 + a1.plot([], [], color=palette.ink2, ls="none", marker="o", ms=3, mfc="white", label="simulation (20 000 trajectoires)")
29 + a1.set_xlabel("Années depuis l'évaluation"); a1.set_ylabel("Loyer contractuel / loyer de marché (espérance)")
30 + a1.set_ylim(0.6, 1.02); a1.set_xlim(0, T + 3.5); a1.legend(fontsize=7, loc="lower right")
31 + for rho, col in ((0.05, palette.cat[2]), (0.15, palette.cat[0]), (0.25, palette.cat[1])):
32 + a2.plot(t, (1 - rho) ** t, color=col, label=f"rotation {rho*100:.0f} %")
33 + a2.set_xlabel("Années depuis l'évaluation"); a2.set_ylabel("Probabilité que le locataire initial soit en place")
34 + a2.yaxis.set_major_formatter(FuncFormatter(lambda v, q: f"{v*100:.0f} %")); a2.set_ylim(0, 1.0); a2.legend(fontsize=7)
35 + savefig(fig, "ch34/ch34-fig-rotation-esperee.pdf")
36 +
37 +if __name__ == "__main__":
38 + main()
added analysis/python/figures/ch34_trajectoires.py +33 −0
@@ -0,0 +1,33 @@
1 +#!/usr/bin/env python3
2 +"""Ch. 34, figure « trajectoires RBP » : revenu brut potentiel contractuel du cas 3 (DONNÉES SYNTHÉTIQUES) selon la rotation
3 +(5, 15, 25 % ; pourcentage de base 3,1 %), avec le RBP au marché (+3,5 %/an), l'indexation seule et le scénario d'encadrement
4 +renforcé (base 2 %, rotation 8 %, rattrapage partiel φ = 0,5, marché 3 %). Sortie : ch34-fig-trajectoires-rbp.pdf."""
5 +from __future__ import annotations
6 +import sys, pathlib
7 +sys.path.insert(0, str(pathlib.Path(__file__).resolve().parents[1]))
8 +import matplotlib.pyplot as plt
9 +from matplotlib.ticker import FuncFormatter
10 +from traite import cas3
11 +from traite.style import setup, palette, savefig, LARGEUR, fmt_milliers
12 +from cases import ch34_tal as c
13 +
14 +def main():
15 + setup()
16 + d = cas3.tout(); hyp = d["hyp"]
17 + fig, ax = plt.subplots(figsize=(LARGEUR, 3.3))
18 + base = c.valeur_scenario(d)["projection"]; x = base["exercice"].values
19 + ax.plot(x, base["rbp_marche"] / 1000, color=palette.ink2, lw=1.2, ls="--", label="RBP au loyer de marché (+3,5 %/an)")
20 + for rho, col, lab in ((0.05, palette.cat[2], "rotation 5 % (1 logement/an)"), (0.15, palette.cat[0], "rotation 15 % (2 logements/an) — base du dossier"), (0.25, palette.cat[1], "rotation 25 % (3 logements/an)")):
21 + p = c.valeur_scenario(d, rotation=rho)["projection"]
22 + ax.plot(x, p["rbp"] / 1000, color=col, marker="o", ms=3, label=f"RBP contractuel, base 3,1 %, {lab}")
23 + p0 = c.valeur_scenario(d, rotation=0.0)["projection"]
24 + ax.plot(x, p0["rbp"] / 1000, color=palette.cat[0], ls=":", lw=1.2, label="Indexation seule (aucune relocation)")
25 + pr = c.valeur_scenario(d, g_c=0.02, rotation=0.08, fraction_rattrapage=0.5, g_m=0.03)["projection"]
26 + ax.plot(x, pr["rbp"] / 1000, color=palette.cat[3], marker="s", ms=3, label="Encadrement renforcé (base 2 %, rotation 8 %, φ = 0,5)")
27 + ax.set_xticks(x); ax.set_xlabel("Exercice"); ax.set_ylabel("Revenu brut potentiel (milliers de $)")
28 + ax.yaxis.set_major_formatter(FuncFormatter(lambda v, q: fmt_milliers(v)))
29 + ax.legend(loc="upper left", fontsize=7.2)
30 + savefig(fig, "ch34/ch34-fig-trajectoires-rbp.pdf")
31 +
32 +if __name__ == "__main__":
33 + main()
added analysis/python/figures/ch34_valeur_scenarios.py +41 −0
@@ -0,0 +1,41 @@
1 +#!/usr/bin/env python3
2 +"""Ch. 34, figure « valeur par scénario » : valeur DCF du cas 3 (DONNÉES SYNTHÉTIQUES) pour la grille rotation (5, 15, 25 %) ×
3 +pourcentage de base (2 ; 3,1 ; 4 %), avec repères : valeur de base (2 208 566 $), capitalisation directe (2 037 503 $), indexation
4 +seule et encadrement renforcé. Sortie : ch34-fig-valeur-scenarios.pdf."""
5 +from __future__ import annotations
6 +import sys, pathlib
7 +sys.path.insert(0, str(pathlib.Path(__file__).resolve().parents[1]))
8 +import numpy as np
9 +import matplotlib.pyplot as plt
10 +from matplotlib.ticker import FuncFormatter
11 +from traite import cas3
12 +from traite.style import setup, palette, savefig, LARGEUR, fmt_milliers
13 +from cases import ch34_tal as c
14 +
15 +def main():
16 + setup()
17 + d = cas3.tout()
18 + g = c.grille_scenarios(d); s = c.scenarios_complementaires(d).set_index("scenario")
19 + rotations = sorted(g["rotation"].unique()); bases = sorted(g["base"].unique())
20 + fig, ax = plt.subplots(figsize=(LARGEUR, 3.3))
21 + w = 0.26; xs = np.arange(len(rotations))
22 + for j, b in enumerate(bases):
23 + vals = [g[(g["rotation"] == r) & (g["base"] == b)]["valeur"].iloc[0] / 1e6 for r in rotations]
24 + bars = ax.bar(xs + (j - 1) * w, vals, w, color=palette.cat[j], label=f"pourcentage de base {b*100:.1f} %".replace(".", ","))
25 + for bar, v in zip(bars, vals):
26 + ax.text(bar.get_x() + bar.get_width() / 2, v + 0.012, f"{v:.3f}".replace(".", ","), ha="center", va="bottom", fontsize=7, color=palette.ink)
27 + V0 = d["dcf"]["valeur"] / 1e6; Vcap = d["valeur_capitalisation"] / 1e6
28 + ax.axhline(V0, color=palette.ink2, lw=0.9, ls="--"); ax.text(2.42, V0 + 0.01, "base du dossier", fontsize=7, color=palette.ink2, ha="right")
29 + ax.axhline(Vcap, color=palette.ink2, lw=0.9, ls=":"); ax.text(2.42, Vcap - 0.045, "capitalisation directe (4,415 %)", fontsize=7, color=palette.ink2, ha="right")
30 + v_idx = s.loc["Indexation seule (aucune relocation)", "valeur"] / 1e6
31 + v_ren = s.loc["Encadrement renforcé (base 2 %, rotation 8 %, φ = 0,5, marché 3 %)", "valeur"] / 1e6
32 + ax.axhline(v_idx, color=palette.cat[3], lw=0.9, ls="-."); ax.text(-0.42, v_idx + 0.01, "indexation seule, sans relocation", fontsize=7, color=palette.cat[3])
33 + ax.axhline(v_ren, color=palette.cat[3], lw=0.9, ls="-"); ax.text(-0.42, v_ren + 0.01, "encadrement renforcé", fontsize=7, color=palette.cat[3])
34 + ax.set_xticks(xs); ax.set_xticklabels([f"rotation {r*100:.0f} %\n({int(round(r*12))} logement{'s' if round(r*12) > 1 else ''}/an)" for r in rotations])
35 + ax.set_ylabel("Valeur DCF (millions de $)"); ax.set_ylim(1.6, 2.45)
36 + ax.yaxis.set_major_formatter(FuncFormatter(lambda v, q: f"{v:.1f}".replace(".", ",")))
37 + ax.legend(loc="upper left", fontsize=7.5, ncol=3)
38 + savefig(fig, "ch34/ch34-fig-valeur-scenarios.pdf")
39 +
40 +if __name__ == "__main__":
41 + main()
modified analysis/python/tests/test_ch14.py +11 −0
@@ -183,3 +183,14 @@ def test_solutions_exercices():
183 183 w_acc = [math.log2(a / b) for a, b in zip(m, u)]; w_dis = [math.log2((1 - a) / (1 - b)) for a, b in zip(m, u)]
184 184 assert sum(w_acc) == pytest.approx(16.35, abs=0.05) and (w_acc[0] + w_acc[1] + w_dis[2]) == pytest.approx(11.2, abs=0.05)
185 185 assert (w_dis[0] + w_acc[1] + w_acc[2]) == pytest.approx(0.78, abs=0.05)
186 +
187 +def test_variante_mnar(R):
188 + """Exercice 14.7 : terrains manquants = les 20 plus petits (MNAR)."""
189 + V = R["mnar"]
190 + assert V.loc["suppression", "n"] == 536 and V.loc["suppression", "ln_T"] == pytest.approx(0.1007, abs=5e-4)
191 + assert V.loc["mediane", "ln_T"] == pytest.approx(0.0977, abs=5e-4) and V.loc["regression", "ln_T"] == pytest.approx(0.0972, abs=5e-4)
192 + assert V.loc["indicateur", "ln_T"] == pytest.approx(0.0978, abs=5e-4)
193 + assert V.loc["suppression", "prediction_sujet"] == pytest.approx(638976, abs=5) and V.loc["mediane", "prediction_sujet"] == pytest.approx(636195, abs=5)
194 + assert V.loc["regression", "prediction_sujet"] == pytest.approx(636224, abs=5) and V.loc["indicateur", "prediction_sujet"] == pytest.approx(638201, abs=5)
195 + assert V.loc["indicateur", "coef_manq_terrain"] == pytest.approx(-0.0502, abs=5e-4) and V.loc["indicateur", "se_manq_terrain"] == pytest.approx(0.0143, abs=5e-4)
196 + assert abs(V.loc["indicateur", "coef_manq_terrain"] / V.loc["indicateur", "se_manq_terrain"]) > 3
added analysis/python/tests/test_ch15.py +136 −0
@@ -0,0 +1,136 @@
1 +"""Chapitre 15 — Études de marché et segmentation : vérification des identités et des chiffres cités dans le texte et les solutions.
2 +Tolérances : indicateurs APCIQ ± 0,05 ; statistiques de test ± 0,005 (F) et ± 0,05 (W) ; silhouettes ± 0,002 ; RMSE ± 0,0005 ; montants ± 1 $."""
3 +import sys, pathlib, math
4 +import numpy as np, pandas as pd, pytest
5 +from scipy import stats
6 +sys.path.insert(0, str(pathlib.Path(__file__).resolve().parents[1] / "cases"))
7 +import ch15_segmentation as c15
8 +from traite import cas1, cas3
9 +
10 +@pytest.fixture(scope="module")
11 +def r():
12 + return c15.tout()
13 +
14 +# ------------------------------------------------------------------ A. indicateurs de marché
15 +def test_mois_inventaire_et_absorption(r):
16 + ind = r["indicateurs"].set_index("territoire")
17 + attendu = {"Province de Québec": 6.1, "RMR de Montréal": 6.6, "RMR de Québec": 3.1, "RMR de Gatineau": 6.5,
18 + "RMR de Sherbrooke": 4.3, "RMR de Trois-Rivières": 3.6, "RMR de Saguenay": 4.0}
19 + for t, mi in attendu.items():
20 + assert ind.loc[t, "mois_inventaire"] == pytest.approx(mi, abs=0.05)
21 + assert ind.loc[t, "taux_absorption_mensuel"] * ind.loc[t, "mois_inventaire"] == pytest.approx(1.0, abs=1e-12) # a = 1/MI
22 + assert ind.loc["RMR de Québec", "taux_absorption_mensuel"] * 100 == pytest.approx(32.7, abs=0.05)
23 + assert ind.loc["RMR de Montréal", "taux_absorption_mensuel"] * 100 == pytest.approx(15.2, abs=0.05)
24 + assert ind.loc["RMR de Montréal", "ratio_stock_flux"] == pytest.approx(5.0, abs=0.05)
25 + assert ind.loc["RMR de Gatineau", "ratio_stock_flux"] == pytest.approx(5.8, abs=0.05)
26 + assert ind.loc["RMR de Québec", "ratio_stock_flux"] == pytest.approx(3.2, abs=0.05)
27 + assert c15.mois_inventaire(20350, 9265) == pytest.approx(20350 / (9265 / 3))
28 +
29 +def test_parts_secteurs_gatineau(r):
30 + p = r["parts_gatineau"]
31 + assert p["ventes_apciq_t3_2026"].sum() == 1044 and p["ventes_synthetiques"].sum() == 600
32 + assert p.loc["Hull", "part_apciq"] * 100 == pytest.approx(16.7, abs=0.05) and p.loc["Hull", "part_synthetique"] * 100 == pytest.approx(31.0, abs=0.05)
33 + assert p.loc["Gatineau", "part_apciq"] * 100 == pytest.approx(33.1, abs=0.05)
34 +
35 +# ------------------------------------------------------------------ B. Chow
36 +def test_chow_identites_et_valeurs(r):
37 + ch = r["chow"]["secteur"]
38 + assert ch["G"] == 5 and ch["k0"] == 12 and ch["n"] == 600 and ch["q"] == 44 and ch["df2"] == 540 and ch["q_classique"] == 48
39 + assert ch["ssr_libre"] <= ch["ssr_contraint"] <= ch["ssr_classique"] # emboîtement des modèles
40 + assert ch["ssr_contraint"] == pytest.approx(4.0245, abs=5e-4) and ch["ssr_libre"] == pytest.approx(3.7360, abs=5e-4) and ch["ssr_classique"] == pytest.approx(6.4024, abs=5e-4)
41 + F = ((ch["ssr_contraint"] - ch["ssr_libre"]) / 44) / (ch["ssr_libre"] / 540)
42 + assert ch["F"] == pytest.approx(F) and ch["F"] == pytest.approx(0.948, abs=0.005) and ch["p"] == pytest.approx(0.570, abs=0.005)
43 + assert ch["F_classique"] == pytest.approx(8.029, abs=0.005) and ch["p_classique"] < 1e-30
44 + assert ch["wald_hc1"] == pytest.approx(54.48, abs=0.05) and ch["p_wald"] == pytest.approx(0.134, abs=0.002) and ch["q_wald"] == 44
45 + assert ch["p"] == pytest.approx(stats.f.sf(ch["F"], 44, 540))
46 + # la SCR du modèle contraint (constantes de secteur + pentes communes) est celle du modèle de référence traite.cas1.hedonique
47 + ref = r["reference"]; assert float(ref["residus"] @ ref["residus"]) == pytest.approx(ch["ssr_contraint"], rel=1e-9)
48 + assert ref["r2"] == pytest.approx(0.853, abs=0.0005) and ref["sigma"] == pytest.approx(0.0830, abs=5e-5)
49 + assert math.sqrt(ch["ssr_classique"] / (600 - 12)) == pytest.approx(0.104, abs=0.001)
50 + # coefficients par secteur cités au tableau
51 + ps = ch["par_segment"]
52 + assert ps["Aylmer"]["coef"]["ln_S"] == pytest.approx(0.513, abs=0.001) and ps["Masson-Angers"]["coef"]["ln_S"] == pytest.approx(0.558, abs=0.001)
53 + assert ps["Hull"]["coef"]["qualite"] == pytest.approx(0.064, abs=0.001) and ps["Buckingham"]["n"] == 42 and ps["Masson-Angers"]["n"] == 47
54 + assert ps["Aylmer"]["coef"]["mois"] == pytest.approx(0.0044, abs=5e-5) and ps["Masson-Angers"]["coef"]["mois"] == pytest.approx(0.0030, abs=5e-5)
55 +
56 +def test_chow_autres_partitions(r):
57 + c = r["chow"]
58 + assert c["secteur4"]["p"] == pytest.approx(0.426, abs=0.005) and c["secteur4"]["q"] == 33
59 + assert c["type"]["k0"] == 11 and c["type"]["regresseurs_retires"] == ["etages"] and c["type"]["q"] == 10
60 + assert c["type"]["p"] == pytest.approx(0.079, abs=0.005) and c["type"]["p_wald"] == pytest.approx(0.0096, abs=0.0005)
61 + assert c["periode"]["p"] == pytest.approx(0.034, abs=0.002) and c["periode"]["p_classique"] == pytest.approx(0.024, abs=0.002) and c["periode"]["p_wald"] == pytest.approx(0.037, abs=0.002)
62 + assert c["periode"]["F"] < stats.f.ppf(0.99, 11, 576) # non rejeté au seuil de 1 %
63 + assert c["kmeans"]["p"] == pytest.approx(0.588, abs=0.005) and c["kmeans"]["p_classique"] < 1e-30
64 + assert c["kmeans5"]["p"] == pytest.approx(0.427, abs=0.005) and c["kmeans5"]["p_wald"] == pytest.approx(0.057, abs=0.002)
65 + assert c["ward"]["F"] == pytest.approx(c["kmeans"]["F"]) # mêmes partitions à k = 2
66 +
67 +# ------------------------------------------------------------------ B. classification
68 +def test_classification(r):
69 + t = r["classification"]["tableau"].set_index("k")
70 + assert r["k_opt"] == 2
71 + assert t.loc[2, "silhouette_kmeans"] == pytest.approx(0.417, abs=0.002) and t.loc[3, "silhouette_kmeans"] == pytest.approx(0.225, abs=0.002)
72 + assert t.loc[5, "silhouette_kmeans"] == pytest.approx(0.188, abs=0.002) and t.loc[5, "silhouette_ward"] == pytest.approx(0.156, abs=0.002)
73 + assert t.loc[2, "ari_kmeans_ward"] == pytest.approx(1.0) and t.loc[5, "ari_kmeans_ward"] == pytest.approx(0.476, abs=0.002)
74 + assert all(np.diff(t["inertie_kmeans"].values) < 0) # l'inertie décroît avec k
75 + assert r["ari_secteur"] == pytest.approx(0.213, abs=0.002) and r["ari_secteur5"] == pytest.approx(0.124, abs=0.002)
76 + # silhouette : implémentation directe = sklearn
77 + Z = r["Z"].values; lab = r["labels_kmeans"]
78 + s_man = c15.silhouette_manuelle(Z, lab)
79 + assert np.allclose(s_man, r["silhouette_i"], atol=1e-10) and -1 <= s_man.min() and s_man.max() <= 1
80 + # segment Est = Buckingham + Masson-Angers
81 + cr = r["croise"]; est = cr.loc[cr.sum(axis=1).idxmin()]
82 + assert est["Buckingham"] == 42 and est["Masson-Angers"] == 47 and est.sum() == 89
83 + p5 = r["profil5"]; assert p5["n"].sum() == 600 and sorted(p5["n"].tolist()) == [89, 96, 98, 120, 197]
84 + assert p5.loc[p5["n"] == 89, "prix_med"].iloc[0] == 371000 and p5.loc[p5["n"] == 96, "qualite_moy"].iloc[0] == pytest.approx(1.81, abs=0.005)
85 + j = r["jeux"]
86 + assert j["coordonnees"]["ari_secteur_k5"] == pytest.approx(0.986, abs=0.002) and j["coordonnees"]["silhouette"] == pytest.approx(0.735, abs=0.002)
87 + assert j["attributs"]["ari_secteur"] == pytest.approx(0.001, abs=0.002) and j["attributs"]["silhouette"] == pytest.approx(0.221, abs=0.002)
88 +
89 +def test_validation_croisee(r):
90 + cv = r["cv"].set_index("specification")
91 + ref = cv.loc["Pentes communes + indicatrices (secteur)"]
92 + assert ref["rmse_log"] == pytest.approx(0.0839, abs=5e-4) and ref["mape_mediane_pct"] == pytest.approx(5.81, abs=0.01)
93 + assert cv.loc["Équation unique sans segment", "rmse_log"] == pytest.approx(0.1053, abs=5e-4)
94 + assert cv.loc["Une équation par segment (secteur)", "mape_mediane_pct"] == pytest.approx(6.26, abs=0.01)
95 + assert cv.loc["Une équation par segment (k-moyennes, k = 5)", "mape_mediane_pct"] == pytest.approx(6.71, abs=0.01)
96 + assert (cv["rmse_log"] >= ref["rmse_log"]).all() # la référence est la meilleure spécification
97 +
98 +# ------------------------------------------------------------------ C. marché locatif
99 +def test_locatif(r):
100 + L = r["locatif"]; rev = L["revenus"]
101 + assert rev["rbp_contrat"] == 165540 and rev["rbp_marche"] == 186036 and rev["ecart_marche"] == 20496
102 + assert L["ecart_global_pct"] == pytest.approx(11.0, abs=0.05)
103 + pt = L["par_type"]
104 + assert pt.loc["3½", "ecart_contrat_marche_pct"] == pytest.approx(15.8, abs=0.05) and pt.loc["5½", "ecart_contrat_marche_pct"] == pytest.approx(4.9, abs=0.05)
105 + assert pt.loc["4½", "ratio_marche_schl"] == pytest.approx(0.994, abs=0.0005) and pt.loc["5½", "ratio_marche_schl"] == pytest.approx(0.944, abs=0.0005)
106 + u = c15.SCHL_ZONE9["univers_2025"]; assert u["studio"] + u["1 ch."] + u["2 ch."] + u["3 ch. +"] == u["tous"] == 44554
107 + assert L["parts_univers"]["2 ch."] * 100 == pytest.approx(50.5, abs=0.05) and L["parts_univers"]["1 ch."] * 100 == pytest.approx(27.8, abs=0.05)
108 + assert L["ecart_inoccupe_occupe_1ch_zone9_pct"] == pytest.approx(22.9, abs=0.05) and L["ecart_inoccupe_occupe_rmr_pct"] == pytest.approx(27.1, abs=0.05)
109 + assert L["rattrapage_an1"]["hypothese_cas3"] == pytest.approx(3074, abs=1) and L["rattrapage_an1"]["schl_ile_6_19_2025"] == pytest.approx(1742, abs=1)
110 + assert L["horizon_renouvellement"]["schl_ile_6_19_2025"] == pytest.approx(11.8, abs=0.05)
111 + assert L["neuf_prime_pct"] == pytest.approx(55.2, abs=0.05) and L["neuf_part_univers_pct"] == pytest.approx(1.41, abs=0.005)
112 + assert 44554 / c15.SCHL_RMR["univers_2025"] * 100 == pytest.approx(6.7, abs=0.05) and 1750 / 1626 - 1 == pytest.approx(0.076, abs=0.0005)
113 +
114 +# ------------------------------------------------------------------ solutions des exercices
115 +def test_solutions_exercices():
116 + # 15.1
117 + assert c15.mois_inventaire(825, 569) == pytest.approx(4.35, abs=0.005) and c15.taux_absorption(825, 569) * 100 == pytest.approx(23.0, abs=0.05)
118 + assert c15.mois_inventaire(825 * 1.2, 569 * 0.9) == pytest.approx(5.80, abs=0.005)
119 + # 15.3
120 + vm = 9265 / 3; d = 40 / c15.JOURS_PAR_MOIS
121 + assert vm * d == pytest.approx(4059, abs=1) and (20350 - vm * d) / 20350 * 100 == pytest.approx(80.1, abs=0.05)
122 + # 15.4
123 + a = (1 + math.sqrt(2)) / 2; b = (3 + math.sqrt(10) + 4) / 3
124 + assert (b - a) / max(a, b) == pytest.approx(0.644, abs=0.001)
125 + Z = np.array([[0, 0], [0, 1], [1, 0], [4, 0], [4, 1], [5, 0]], float); lab = np.array([0, 0, 0, 1, 1, 1])
126 + assert c15.silhouette_manuelle(Z, lab)[2] == pytest.approx(0.6437, abs=1e-4)
127 + # 15.5
128 + F = ((4.0245 - 3.7360) / 44) / (3.7360 / 540); F0 = ((6.4024 - 3.7360) / 48) / (3.7360 / 540)
129 + assert F == pytest.approx(0.948, abs=0.001) and F0 == pytest.approx(8.03, abs=0.005)
130 + assert stats.f.ppf(0.95, 44, 540) == pytest.approx(1.40, abs=0.005) and stats.f.ppf(0.95, 48, 540) == pytest.approx(1.38, abs=0.005)
131 + # 15.6
132 + assert 1 - 0.85 ** 5 == pytest.approx(0.556, abs=0.001) and 1 - 0.906 ** 5 == pytest.approx(0.390, abs=0.001) and 1 - 0.915 ** 5 == pytest.approx(0.359, abs=0.001)
133 + assert 1 - 0.85 ** 10 == pytest.approx(0.803, abs=0.001) and 1 - 0.915 ** 10 == pytest.approx(0.589, abs=0.001)
134 + assert 1 - 875.25 / 1039 == pytest.approx(0.158, abs=0.001) and 1295 / 1054 - 1 == pytest.approx(0.229, abs=0.001)
135 + # probabilité d'au moins un rejet sur six tests indépendants au seuil de 5 %
136 + assert 1 - 0.95 ** 6 == pytest.approx(0.265, abs=0.001)
modified analysis/python/tests/test_ch22.py +11 −0
@@ -61,6 +61,17 @@ class TestGeometrie:
61 61 assert S["aire_etages_mefq"] * c21.PI2_PAR_M2 == pytest.approx(2046, abs=0.6)
62 62 assert S["habitable"] * c21.PI2_PAR_M2 == pytest.approx(1777, abs=0.6)
63 63
64 + def test_couts_unitaires_conventions(self, data):
65 + S, cn = data["S"], data["cn"]
66 + assert cn["total"] / S["aire_etages_mefq"] == pytest.approx(2092, abs=0.6)
67 + assert cn["total"] / S["habitable"] == pytest.approx(2408, abs=0.6)
68 + assert cn["total"] / S["nette_hors_sol"] == pytest.approx(2937, abs=0.6)
69 + assert S["habitable"] / S["aire_etages_mefq"] == pytest.approx(0.869, abs=5e-4)
70 + assert S["nette_hors_sol"] / S["habitable"] == pytest.approx(0.820, abs=5e-4)
71 + assert S["nette_hors_sol"] / S["aire_etages_mefq"] == pytest.approx(0.712, abs=5e-4)
72 + # écart de 40 % entre les coûts unitaires extrêmes
73 + assert (cn["total"] / S["nette_hors_sol"]) / (cn["total"] / S["aire_etages_mefq"]) - 1 == pytest.approx(0.404, abs=5e-3)
74 +
64 75 def test_pieces_couvrent_les_niveaux(self, data):
65 76 pcs = data["S"]["pieces"]
66 77 assert pcs.loc[pcs["niveau"] == "rez", "brute"].sum() == pytest.approx(95.04, abs=1e-6)
added analysis/python/tests/test_ch31.py +210 −0
@@ -0,0 +1,210 @@
1 +"""Tests du chapitre 31 — financement, levier, RCD et RPV (cas 3, données synthétiques ; repères BdC réels).
2 +Identités : taux périodique canadien ; paiement × facteur VA = capital ; solde prospectif = solde du tableau ; solde final nul ;
3 +Σ intérêts + capital = Σ paiements ; constante = 12·PMT/M_0 ; prêt max par RCD ⇒ RCD = seuil ; R_E = (R − Mf)/(1 − M) ;
4 +levier positif ssi y > k ; renouvellement exact : solde au terme reporté, TRI exact ≥ approximation ; barème APH ; mise de fonds ;
5 +valeurs citées dans manuscript/chapters/ch31.tex (tolérances déclarées)."""
6 +import math, pytest
7 +import numpy as np
8 +from traite import finance as fin
9 +from cases import ch31_financement as c
10 +
11 +@pytest.fixture(scope="module")
12 +def R():
13 + return c.base()
14 +
15 +@pytest.fixture(scope="module")
16 +def FI(R):
17 + return c.financement(R)
18 +
19 +# ------------------------------------------------------------------ hypothèque canadienne
20 +
21 +def test_taux_canadien_et_effectif():
22 + im = fin.taux_mensuel_canadien(0.0412)
23 + assert im == pytest.approx((1 + 0.0412 / 2) ** (1 / 6) - 1, rel=1e-14)
24 + assert im == pytest.approx(0.0034042, abs=5e-8)
25 + assert fin.taux_effectif_annuel(0.0412, 2) == pytest.approx(0.041624, abs=1e-6)
26 + assert (1 + im) ** 12 - 1 == pytest.approx(fin.taux_effectif_annuel(0.0412, 2), rel=1e-12)
27 + assert 0.0412 / 12 == pytest.approx(0.0034333, abs=5e-8)
28 +
29 +def test_identites_hypotheque(FI):
30 + h = FI["hypotheque"]; n = h.amortissement_mois; i = h.taux_periodique
31 + assert h.paiement * fin.facteur_va_annuite(i, n) == pytest.approx(h.capital, rel=1e-12)
32 + tab = h.tableau()
33 + assert tab[-1]["solde"] == pytest.approx(0.0, abs=1e-6)
34 + for k in (1, 12, 60, 120, 240, 479):
35 + assert tab[k - 1]["solde"] == pytest.approx(h.solde(k), rel=1e-9)
36 + assert sum(r["interet"] + r["capital"] for r in tab) == pytest.approx(h.paiement * n, rel=1e-9)
37 + assert sum(r["capital"] for r in tab) == pytest.approx(h.capital, rel=1e-9)
38 + assert h.constante_hypothecaire == pytest.approx(12 * h.paiement / h.capital, rel=1e-14)
39 + # capital remboursé croît géométriquement
40 + assert tab[1]["capital"] / tab[0]["capital"] == pytest.approx(1 + i, rel=1e-9)
41 +
42 +def test_valeurs_citees_financement(FI):
43 + assert FI["valeur"] == pytest.approx(2_208_566, abs=1)
44 + assert FI["rne_1"] == pytest.approx(89_956, abs=1)
45 + assert FI["pret_rpv"] == pytest.approx(1_656_424, abs=1)
46 + assert FI["pret_rcd"] == pytest.approx(1_475_959, abs=1)
47 + assert FI["pret"] == FI["pret_rcd"] and FI["contrainte"] == "RCD"
48 + assert FI["rpv"] == pytest.approx(0.6683, abs=5e-5)
49 + assert FI["mise_de_fonds"] == pytest.approx(732_606, abs=1)
50 + assert FI["paiement"] == pytest.approx(6246.93, abs=0.01)
51 + assert FI["service_dette"] == pytest.approx(74_963, abs=1)
52 + assert FI["constante"] == pytest.approx(0.050789, abs=1e-6)
53 + assert FI["interets_an1"] == pytest.approx(60_016, abs=1) and FI["capital_an1"] == pytest.approx(14_947, abs=1)
54 + assert FI["solde_60"] == pytest.approx(1_394_739, abs=1) and FI["solde_120"] == pytest.approx(1_295_148, abs=1)
55 + assert FI["rcd"] == pytest.approx(1.20, abs=1e-9)
56 + assert FI["flux_apres_dette_1"] == pytest.approx(14_993, abs=1)
57 + assert FI["cash_on_cash_1"] == pytest.approx(0.0205, abs=5e-5)
58 + assert FI["point_mort"] == pytest.approx(0.927, abs=5e-4) and FI["point_mort_hors_reserve"] == pytest.approx(0.878, abs=5e-4)
59 + assert FI["tri_fp"] == pytest.approx(0.1122, abs=5e-5)
60 + assert FI["k_annuel"] == pytest.approx(0.0406, abs=5e-5) and FI["k_annuel"] < FI["taux_effectif"]
61 + assert FI["rpv_max_par_rcd_sur_Vcap"] == pytest.approx(0.724, abs=5e-4)
62 +
63 +def test_conventions(FI):
64 + d = c.conventions(FI["pret"]).set_index("convention")
65 + ca, us = d.loc["Canada (semestrielle)"], d.loc["États-Unis (mensuelle)"]
66 + assert ca["paiement"] < us["paiement"] and ca["taux_effectif"] < us["taux_effectif"]
67 + assert us["paiement"] == pytest.approx(6279.21, abs=0.01)
68 + assert us["service_dette"] - ca["service_dette"] == pytest.approx(387, abs=1)
69 + assert us["constante"] == pytest.approx(0.051052, abs=1e-6)
70 + assert us["solde_60"] == pytest.approx(1_395_377, abs=1)
71 +
72 +def test_tableau_extraits(FI):
73 + t = c.tableau_extraits(FI["hypotheque"]).set_index("periode")
74 + assert t.loc[1, "interet"] == pytest.approx(5024.50, abs=0.01) and t.loc[1, "capital"] == pytest.approx(1222.42, abs=0.01)
75 + assert t.loc[60, "solde"] == pytest.approx(1_394_739.19, abs=0.01)
76 + assert t.loc[480, "solde"] == pytest.approx(0.0, abs=1e-6)
77 + assert t.loc[1, "part_interet"] == pytest.approx(0.804, abs=1e-3)
78 + dec = c.decomposition_annuelle(FI["hypotheque"])
79 + assert dec["interet"].sum() == pytest.approx(1_522_566, abs=1)
80 + assert (dec["interet"] + dec["capital"] - FI["service_dette"]).abs().max() < 1e-6
81 +
82 +# ------------------------------------------------------------------ ratios, prêt max, règles
83 +
84 +def test_pret_max_par_rcd_donne_rcd_seuil(R):
85 + noi = R["base"]["rne"]
86 + for i, n, q in ((0.0412, 480, 1.2), (0.0612, 300, 1.1), (0.0312, 360, 1.3)):
87 + p = fin.pret_max_par_rcd(noi, q, i, n)
88 + assert fin.rcd(noi, fin.Hypotheque(p, i, n).service_dette_annuel) == pytest.approx(q, rel=1e-12)
89 + g = c.grille_pret_max(R)
90 + assert (g.groupby(["amortissement_mois", "rcd_min"])["pret_max"].diff().dropna() < 0).all() # décroissant en taux
91 + row = g[(g["taux"] == 0.0412) & (g["amortissement_mois"] == 480) & (g["rcd_min"] == 1.2)].iloc[0]
92 + assert row["pret_max"] == pytest.approx(1_475_959, abs=1)
93 + r300 = g[(g["taux"] == 0.0412) & (g["amortissement_mois"] == 300) & (g["rcd_min"] == 1.2)].iloc[0]["pret_max"]
94 + assert row["pret_max"] / r300 == pytest.approx(1.258, abs=2e-3)
95 +
96 +def test_seuil_rpv_max_par_rcd(FI):
97 + # M_max/V = R/(RCD·f) ; seuil de liaison R = RCD · f · LTVmax = 4,571 %
98 + assert FI["rpv_max_par_rcd_sur_V"] == pytest.approx(FI["rpv"], rel=1e-9)
99 + assert 1.2 * FI["constante"] * 0.75 == pytest.approx(0.04571, abs=5e-5)
100 +
101 +def test_regles_residentielles():
102 + assert fin.taux_admissible_minimal(0.0458) == pytest.approx(0.0658) and fin.taux_admissible_minimal(0.03) == pytest.approx(0.0525)
103 + assert c.mise_de_fonds_minimale(400_000) == 20_000 and c.mise_de_fonds_minimale(650_000) == 40_000
104 + assert c.mise_de_fonds_minimale(1_000_000) == 75_000 and c.mise_de_fonds_minimale(1_500_000) == 125_000
105 + assert c.mise_de_fonds_minimale(1_600_000) == 320_000
106 + rr = c.regles_residentielles()
107 + assert rr["capacite"]["pret_au_contrat"] == pytest.approx(594_162, abs=1) and rr["capacite"]["pret_au_mqr"] == pytest.approx(491_344, abs=1)
108 + assert rr["capacite"]["reduction_mqr"] == pytest.approx(0.173, abs=5e-4)
109 +
110 +def test_prime_aph_bareme():
111 + assert c.prime_aph(0.60, 300) == pytest.approx(0.0260) and c.prime_aph(0.65, 300) == pytest.approx(0.0260)
112 + assert c.prime_aph(0.72, 480) == pytest.approx(0.0335 + 0.0075) and c.prime_aph(0.78, 480) == pytest.approx(0.0435 + 0.0075)
113 + assert c.prime_aph(0.817, 540) == pytest.approx(0.0535 + 0.0100) and c.prime_aph(0.845, 600) == pytest.approx(0.0535 + 0.0125)
114 + assert c.prime_aph(0.95, 600) == pytest.approx(0.0615 + 0.0125)
115 + # cohérence avec les valeurs du ch. 47 (RPV 78,2 / 81,7 / 84,5 % sur 40/45/50 ans → 5,10 / 6,35 / 6,60 %)
116 + assert c.prime_aph(0.782, 480) == pytest.approx(0.0510) and c.prime_aph(0.817, 540) == pytest.approx(0.0635) and c.prime_aph(0.845, 600) == pytest.approx(0.0660)
117 +
118 +def test_aph_select(R):
119 + a = c.aph_select(R).set_index("amortissement_mois")
120 + assert a.loc[480, "pret"].iloc[-1] if hasattr(a.loc[480, "pret"], "iloc") else True
121 + aph = c.aph_select(R).iloc[1:]
122 + assert list(aph["pret"].round(0)) == [1_610_138, 1_682_397, 1_741_328]
123 + assert (aph["pret_rcd"] < aph["pret_rpv"]).all() # RCD liant
124 + assert (aph["rcd_1"] - 1.10).abs().max() < 1e-9
125 + assert aph["service_dette"].iloc[0] == pytest.approx(81_778, abs=1)
126 + assert list(aph["prime_pct"].round(4)) == [0.0410, 0.0535, 0.0560]
127 + assert aph["prime"].iloc[0] == pytest.approx(66_016, abs=1) and aph["mise_de_fonds"].iloc[-1] == pytest.approx(467_238, abs=1)
128 +
129 +def test_stress_rcd(R):
130 + s = c.stress_rcd(R).set_index("marge")
131 + assert s.loc[0.0, "rcd_1"] == pytest.approx(1.2, abs=1e-9)
132 + assert s.loc[0.02, "service_dette"] == pytest.approx(97_988, abs=1) and s.loc[0.02, "rcd_1"] == pytest.approx(0.918, abs=5e-4)
133 + assert s.loc[0.02, "pret_max_qualifie"] == pytest.approx(1_129_149, abs=1)
134 + assert (s["rcd_1"].diff().dropna() < 0).all()
135 +
136 +# ------------------------------------------------------------------ levier
137 +
138 +def test_levier_courant_formule(R):
139 + l = c.levier_courant(R)
140 + assert l["RE_formule"] == pytest.approx(l["RE_direct"], rel=1e-12) and l["RE_alt"] == pytest.approx(l["RE_direct"], rel=1e-12)
141 + assert l["RE_direct"] == pytest.approx(0.0205, abs=5e-5) and l["R"] < l["f"] and l["levier_courant"] == "négatif"
142 + assert l["RE_cap"] == pytest.approx(0.0267, abs=5e-5)
143 +
144 +def test_levier_par_rpv_monotone(R, FI):
145 + t = c.levier_par_rpv(R)
146 + assert t["tri_fp"].iloc[0] == pytest.approx(0.0725, abs=1e-6) # sans levier : TRI = r
147 + assert (t["tri_fp"].diff().dropna() > 0).all() # y > k ⇒ TRI FP croît avec M
148 + assert (t["cash_on_cash_1"].diff().dropna() < 0).all() # R < f ⇒ rendement courant décroît
149 + ret = t[t["retenu"]].iloc[0]
150 + assert ret["tri_fp"] == pytest.approx(FI["tri_fp"], rel=1e-9) and ret["rcd_1"] == pytest.approx(1.2, abs=1e-9)
151 + t85 = t[np.isclose(t["rpv"], 0.85)].iloc[0]
152 + assert t85["tri_fp"] == pytest.approx(0.1524, abs=5e-5) and t85["cash_on_cash_1"] < 0 and t85["rcd_1"] < 1
153 + t75 = t[np.isclose(t["rpv"], 0.75)].iloc[0]
154 + assert t75["rcd_1"] == pytest.approx(1.07, abs=5e-3) and t75["flux_min"] == pytest.approx(-69_617, abs=1)
155 +
156 +def test_proposition_levier_positif(R, FI):
157 + """y_E > y ⇔ y > k : vérifié sur le cas 3 (y = 7,25 %, k ≈ 4,06 %) et par un contre-exemple où la dette coûte plus que l'actif."""
158 + assert FI["tri_sans_levier"] > FI["k_annuel"] and FI["tri_fp"] > FI["tri_sans_levier"]
159 + d = R["dcf"]; V = d["valeur"]
160 + cher = fin.Hypotheque(0.5 * V, 0.09, 480) # dette à 9 % > y = 7,25 %
161 + res = fin.rendement_fonds_propres_levier(d["flux"], V, cher, d["vt_nette"])
162 + assert res["tri"] < 0.0725
163 + neutre_k = 0.0725
164 + # dette dont le TRI annuel = y exactement (taux effectif y, paiements annuels) ⇒ y_E = y
165 + sd = 0.5 * V * fin.facteur_recouvrement(neutre_k, 10)
166 + flux = [-(0.5 * V)] + [cf - sd for cf in d["flux"]]; flux[-1] += d["vt_nette"]
167 + assert fin.tri(flux) == pytest.approx(0.0725, abs=1e-9)
168 +
169 +def test_risque_levier(R):
170 + rk = c.risque_levier(R)
171 + assert rk["sigma_fp"] == pytest.approx(0.08 / (1 - rk["M"])) and rk["multiplicateur"] == pytest.approx(3.01, abs=5e-3)
172 +
173 +# ------------------------------------------------------------------ renouvellement
174 +
175 +def test_renouvellement_exact(R, FI):
176 + ren = c.renouvellement(R); t = ren["tableau"].set_index("delta")
177 + assert ren["solde_terme"] == pytest.approx(FI["solde_60"], rel=1e-12) and ren["n_residuel"] == 420
178 + assert t.loc[0.0, "paiement_2"] == pytest.approx(FI["paiement"], rel=1e-9) # taux inchangé ⇒ même paiement
179 + assert t.loc[0.0, "tri_fp"] == pytest.approx(FI["tri_fp"], rel=1e-9) and t.loc[0.0, "solde_T"] == pytest.approx(FI["solde_120"], rel=1e-9)
180 + assert list(t["service_dette_2"].round(0)) == [74_963, 85_177, 95_920, 107_120] # identiques au ch. 47
181 + assert list(t["rcd_an_terme_plus_1"].round(3)) == [1.746, 1.537, 1.365, 1.222]
182 + assert list(t["tri_fp"].round(4)) == [0.1122, 0.1078, 0.1032, 0.0986]
183 + assert list(t["tri_fp_approx_ch30"].round(4)) == [0.1122, 0.1003, 0.0881, 0.0757]
184 + assert (t["tri_fp"] >= t["tri_fp_approx_ch30"] - 1e-12).all() and (t["sd_approx_ch30"] >= t["service_dette_2"] - 1e-6).all()
185 + assert t.loc[0.02, "solde_T"] == pytest.approx(1_327_042, abs=1)
186 + # le solde renouvelé suit la formule prospective avec le nouveau taux
187 + h2 = ren["detail"][0.02]["hypotheque_2"]
188 + assert h2.capital == pytest.approx(ren["solde_terme"]) and h2.solde(60) == pytest.approx(t.loc[0.02, "solde_T"], rel=1e-12)
189 + assert ren["delta_rcd_120"] == pytest.approx(0.0317, abs=1e-4) and ren["delta_rcd_100"] == pytest.approx(0.0502, abs=1e-4)
190 +
191 +def test_sensibilite_service_dette(R, FI):
192 + s = c.sensibilite_service_dette(R, taux=np.array([0.0412, 0.055, 0.0612]))
193 + b = s[(s["amortissement_mois"] == 480) & np.isclose(s["taux"], 0.0412)].iloc[0]
194 + assert b["service_dette"] == pytest.approx(FI["service_dette"], rel=1e-9) and b["rcd_1"] == pytest.approx(1.2, abs=1e-9)
195 + assert s[(s["amortissement_mois"] == 480) & np.isclose(s["taux"], 0.055)].iloc[0]["rcd_1"] < 1.0
196 + assert s[(s["amortissement_mois"] == 300) & np.isclose(s["taux"], 0.0412)].iloc[0]["rcd_1"] == pytest.approx(0.954, abs=5e-4)
197 +
198 +# ------------------------------------------------------------------ exercices
199 +
200 +def test_exercices_coherence():
201 + e = c.exercices()
202 + ca, us = e["ex1"]["Canada (semestrielle)"], e["ex1"]["États-Unis (mensuelle)"]
203 + assert ca["paiement"] == pytest.approx(4742.40, abs=0.01) and us["paiement"] == pytest.approx(4763.26, abs=0.01)
204 + assert e["ex2"]["solde_60"] == pytest.approx(944_971, abs=1) and e["ex2"]["part_remboursee_60"] == pytest.approx(0.0550, abs=5e-5)
205 + assert e["ex2"]["interets_totaux"] == pytest.approx(1_031_577, abs=1) and e["ex2"]["interets_totaux_25"] == pytest.approx(597_637, abs=1)
206 + assert e["ex3"]["pret"] == pytest.approx(1_572_079, abs=1) and e["ex3"]["pret_aph"] == pytest.approx(2_147_906, abs=1) and e["ex3"]["prime_aph"] == pytest.approx(0.0610)
207 + assert e["ex4"][0.7][0.06] == pytest.approx(0.03333, abs=1e-5) and e["ex4"][0.7][0.045] == pytest.approx(0.06833, abs=1e-5)
208 + assert e["ex5"]["mdf"] == 47_000 and e["ex5"]["paiement_mqr"] == pytest.approx(4540.60, abs=0.01) and e["ex5"]["lti_revenu_requis"] == pytest.approx(149_556, abs=1)
209 + assert e["ex6"]["service_dette_2"] == pytest.approx(90_487, abs=1) and e["ex6"]["tri_fp"] == pytest.approx(0.1055, abs=5e-5)
210 + assert e["ex7"][0]["tri_fp"] == pytest.approx(0.0949, abs=5e-5) and e["ex7"][1]["tri_fp"] == pytest.approx(0.1524, abs=5e-5)
modified assets/figures/ch14/ch14-fig-avant-apres.pdf +0 −0

Binary file not shown.

modified assets/figures/ch14/ch14-fig-geocodage.pdf +0 −0

Binary file not shown.

modified assets/figures/ch14/ch14-fig-levier-cook.pdf +0 −0

Binary file not shown.

modified assets/figures/ch14/ch14-fig-scores-robustes.pdf +0 −0

Binary file not shown.

modified assets/figures/ch14/ch14-fig-sensibilite-seuils.pdf +0 −0

Binary file not shown.

added assets/figures/ch20/ch20-fig-bootstrap-densite.pdf +0 −0

Binary file not shown.

added assets/figures/ch20/ch20-fig-fourchettes-methodes.pdf +0 −0

Binary file not shown.

added assets/figures/ch20/ch20-fig-residus-elargissement.pdf +0 −0

Binary file not shown.

added assets/figures/ch20/ch20-fig-sensibilite-poids.pdf +0 −0

Binary file not shown.

added assets/figures/ch34/ch34-fig-calcul-tal.pdf +0 −0

Binary file not shown.

added assets/figures/ch34/ch34-fig-ipc-base.pdf +0 −0

Binary file not shown.

added assets/figures/ch34/ch34-fig-rotation-esperee.pdf +0 −0

Binary file not shown.

added assets/figures/ch34/ch34-fig-trajectoires-rbp.pdf +0 −0

Binary file not shown.

added assets/figures/ch34/ch34-fig-valeur-scenarios.pdf +0 −0

Binary file not shown.

added data/processed/ipc_quebec_mensuel.csv +578 −0
Large diff (578 lines) is not loaded by default to keep the page fast.Load diff
modified manuscript/chapters/ch13.tex +1 −1
Large diff (8 lines) is not loaded by default to keep the page fast.Load diff
modified manuscript/chapters/ch14.tex +25 −25
Large diff (115 lines) is not loaded by default to keep the page fast.Load diff
modified manuscript/chapters/ch15.tex +33 −33
Large diff (146 lines) is not loaded by default to keep the page fast.Load diff
modified manuscript/chapters/ch22.tex +21 −21
Large diff (107 lines) is not loaded by default to keep the page fast.Load diff
modified manuscript/chapters/ch23.tex +569 −1
Large diff (571 lines) is not loaded by default to keep the page fast.Load diff
modified manuscript/chapters/ch24.tex +536 −1
Large diff (538 lines) is not loaded by default to keep the page fast.Load diff
modified manuscript/chapters/ch26.tex +570 −1
Large diff (572 lines) is not loaded by default to keep the page fast.Load diff
modified manuscript/chapters/ch27.tex +477 −1
Large diff (479 lines) is not loaded by default to keep the page fast.Load diff
modified manuscript/chapters/ch28.tex +522 −1
Large diff (524 lines) is not loaded by default to keep the page fast.Load diff
modified manuscript/chapters/ch31.tex +540 −1
Large diff (542 lines) is not loaded by default to keep the page fast.Load diff
modified manuscript/chapters/ch33.tex +170 −1
Large diff (172 lines) is not loaded by default to keep the page fast.Load diff
added references/lib-ch13.bib +11 −0
Large diff (11 lines) is not loaded by default to keep the page fast.Load diff
added references/lib-ch15.bib +68 −0
Large diff (68 lines) is not loaded by default to keep the page fast.Load diff
added references/lib-ch23.bib +48 −0
Large diff (48 lines) is not loaded by default to keep the page fast.Load diff
modified references/library.bib +110 −0
Large diff (119 lines) is not loaded by default to keep the page fast.Load diff
modified solutions/ch14.tex +29 −0
Large diff (30 lines) is not loaded by default to keep the page fast.Load diff
modified solutions/ch15.tex +37 −0
Large diff (38 lines) is not loaded by default to keep the page fast.Load diff
modified solutions/ch23.tex +61 −0
Large diff (62 lines) is not loaded by default to keep the page fast.Load diff
modified solutions/ch24.tex +66 −0
Large diff (67 lines) is not loaded by default to keep the page fast.Load diff
modified solutions/ch27.tex +33 −0
Large diff (34 lines) is not loaded by default to keep the page fast.Load diff
modified solutions/ch31.tex +75 −0
Large diff (76 lines) is not loaded by default to keep the page fast.Load diff
modified sources/registry.csv +5 −0
Large diff (11 lines) is not loaded by default to keep the page fast.Load diff
modified verification/calculations.csv +45 −0
Large diff (48 lines) is not loaded by default to keep the page fast.Load diff
added verification/chapters/ch13-calculations.csv +15 −0
Large diff (15 lines) is not loaded by default to keep the page fast.Load diff
added verification/chapters/ch13-claims.csv +26 −0
Large diff (26 lines) is not loaded by default to keep the page fast.Load diff
added verification/chapters/ch13-figures.csv +6 −0
Large diff (6 lines) is not loaded by default to keep the page fast.Load diff
added verification/chapters/ch13-registry.csv +6 −0
Large diff (6 lines) is not loaded by default to keep the page fast.Load diff
added verification/chapters/ch13.md +42 −0
Large diff (42 lines) is not loaded by default to keep the page fast.Load diff
added verification/chapters/ch14-calculations.csv +19 −0
Large diff (19 lines) is not loaded by default to keep the page fast.Load diff
added verification/chapters/ch14-figures.csv +7 −0
Large diff (7 lines) is not loaded by default to keep the page fast.Load diff
added verification/chapters/ch15-calculations.csv +14 −0
Large diff (14 lines) is not loaded by default to keep the page fast.Load diff
added verification/chapters/ch15-figures.csv +7 −0
Large diff (7 lines) is not loaded by default to keep the page fast.Load diff
added verification/chapters/ch22-claims.csv +19 −0
Large diff (19 lines) is not loaded by default to keep the page fast.Load diff
modified verification/chapters/ch22-figures.csv +2 −2
Large diff (7 lines) is not loaded by default to keep the page fast.Load diff
modified verification/claims.csv +43 −0
Large diff (46 lines) is not loaded by default to keep the page fast.Load diff
modified verification/figures.csv +17 −0
Large diff (21 lines) is not loaded by default to keep the page fast.Load diff
Menu

Theme