TypeScript 90.1%
JavaScript 3.2%
Python 3.2%
CSS 1.8%
HTML 1.8%
1// UQO Éval — modélisation hédonique : orchestration (échantillon → design → estimation → diagnostics → rapport).2import { buildSample, type Obs } from "./data";3import { buildDesign, type Column, type Design } from "./design";4import { breuschPagan, fitMetrics, histogram, jarqueBera, qqPoints, smearing, toDollars, vif, type FitMetrics } from "./diagnostics";5import { covariance, huber, lad, ols, penalized, sar, sarPredict, type Fit } from "./estimators";6import { hcat, mean, median, quantile, sd, sortedCopy, takeRows, takeVec, type Mat } from "./linalg";7import { KnnIndex, lag, moranI, type Moran } from "./spatial";8import { LIMITS, METHODS, type HedonicSpec } from "./spec";9import { fSf, stars, tPValue, tQuantile } from "./stats";1011type Txt = { fr: string; en: string };1213export interface Coefficient {14 name: string;15 label: Txt;16 kind: Column["kind"] | "rho";17 var?: string;18 level?: string;19 base?: string;20 n?: number;21 coef: number;22 se: number | null;23 t: number | null;24 p: number | null;25 stars: string;26 ciLo: number | null;27 ciHi: number | null;28 /** effet économique : % du prix et $ au prix médian, pour le pas indiqué */29 effectPct: number | null;30 effectDollar: number | null;31 effectLo: number | null; // IC 95 % de l'effet en %32 effectHi: number | null;33 step: Txt;34}3536export interface Descriptive {37 key: string;38 label: Txt;39 n: number;40 mean: number;41 sd: number;42 min: number;43 p25: number;44 median: number;45 p75: number;46 max: number;47}4849export interface HedonicResult {50 spec: HedonicSpec;51 sample: {52 n: number;53 available: number;54 trimmed: number;55 dropped: { reason: Txt; n: number }[];56 period: { from: string; to: string };57 byType: { type: string; n: number }[];58 descriptives: Descriptive[];59 medianPrice: number;60 center: { lat: number; lng: number };61 };62 model: {63 method: HedonicSpec["method"];64 methodLabel: Txt;65 depvar: HedonicSpec["depvar"];66 equation: string;67 n: number;68 k: number;69 df: number;70 seType: HedonicSpec["se"];71 nClusters: number;72 coefficients: Coefficient[];73 fit: FitMetrics;74 fstat: { f: number; df1: number; df2: number; p: number } | null;75 lambda: { value: number; grid: { lambda: number; cvErr: number }[] } | null;76 lassoSelected: string[] | null;77 postLasso: Coefficient[] | null;78 rho: Coefficient | null;79 robust: { scale: number; downweightedPct: number } | null;80 notes: Txt[];81 };82 diagnostics: {83 bp: { stat: number; df: number; p: number } | null;84 jb: { stat: number; p: number; skew: number; kurt: number };85 moran: Moran | null;86 vif: { name: string; label: Txt; vif: number }[];87 };88 cv: { folds: number; metrics: FitMetrics; gapMdape: number } | null;89 lab: { n: number; lab: { mdape: number; rmse: number; within10: number; within20: number }; ours: { mdape: number; rmse: number; within10: number; within20: number }; label: Txt } | null;90 timeIndex: { level: string; n: number; coef: number; index: number; lo: number | null; hi: number | null }[] | null;91 geoEffects: { level: string; n: number; coef: number; pct: number; se: number | null }[] | null;92 plots: {93 residFitted: { f: number; e: number }[];94 qq: { q: number; e: number }[];95 hist: { lo: number; hi: number; n: number }[];96 map: { lat: number; lng: number; e: number }[];97 actualPred: { a: number; p: number }[];98 };99 insights: Txt[];100 code: { r: string; python: string; stata: string };101 timingMs: number;102}103104const T = (fr: string, en: string): Txt => ({ fr, en });105const fmtPct = (x: number, d = 1) => `${x >= 0 ? "+" : "−"}${Math.abs(x).toFixed(d)} %`;106const fmtMoney = (x: number) => `${Math.round(x).toLocaleString("fr-CA")} $`;107const fmtMoneyEn = (x: number) => `$${Math.round(x).toLocaleString("en-CA")}`;108109function descriptives(design: Design): Descriptive[] {110 const { X, cols } = design;111 const out: Descriptive[] = [];112 const pushVec = (key: string, label: Txt, v: Float64Array) => {113 const s = sortedCopy(v);114 out.push({ key, label, n: v.length, mean: mean(v), sd: sd(v), min: s[0], p25: quantile(s, 0.25), median: quantile(s, 0.5), p75: quantile(s, 0.75), max: s[s.length - 1] });115 };116 pushVec("price", T("Prix ($)", "Price ($)"), Float64Array.from(design.obs, (o) => o.price));117 for (const c of cols) {118 if (c.kind !== "cont" || c.var === "age2") continue;119 const v = new Float64Array(X.n);120 const j = cols.indexOf(c);121 for (let i = 0; i < X.n; i++) v[i] = c.log ? Math.exp(X.a[i * X.k + j]) : X.a[i * X.k + j];122 pushVec(c.var!, c.log ? { fr: c.label.fr.replace(/^ln /, ""), en: c.label.en.replace(/^ln /, "") } : c.label, v);123 }124 return out;125}126127function coefTable(cols: Column[], beta: Float64Array, V: Float64Array | null, df: number, logModel: boolean, medP: number, offset = 0): Coefficient[] {128 const k = cols.length;129 const tcrit = V ? tQuantile(0.975, df) : NaN;130 return cols.map((c, j) => {131 const b = beta[j + offset];132 const se = V ? Math.sqrt(Math.max(0, V[(j + offset) * (k + offset) + (j + offset)])) : null;133 const t = se ? b / se : null;134 const p = t != null ? tPValue(t, df) : null;135 let step = T("par unité", "per unit");136 // effet économique d'un coefficient b (en % du prix et en $ au prix médian)137 let eff: (x: number) => { pct: number; dollar: number } | null = () => null;138 if (c.kind === "intercept") {139 step = T("—", "—");140 } else if (c.kind === "cont" && c.log) {141 step = T("pour +10 %", "for +10%");142 eff = logModel ? (x) => ({ pct: 100 * (Math.pow(1.1, x) - 1), dollar: medP * (Math.pow(1.1, x) - 1) }) : (x) => ({ pct: (100 * x * Math.log(1.1)) / medP, dollar: x * Math.log(1.1) });143 } else {144 if (c.kind !== "cont" && c.kind !== "trend") step = T("vs référence", "vs base");145 eff = logModel ? (x) => ({ pct: 100 * (Math.exp(x) - 1), dollar: medP * (Math.exp(x) - 1) }) : (x) => ({ pct: (100 * x) / medP, dollar: x });146 }147 const e0 = eff(b);148 const eLo = se != null ? eff(b - tcrit * se) : null;149 const eHi = se != null ? eff(b + tcrit * se) : null;150 const effectPct = e0?.pct ?? null;151 const effectDollar = e0?.dollar ?? null;152 return {153 name: c.name,154 label: c.label,155 kind: c.kind,156 var: c.var,157 level: c.level,158 base: c.base,159 n: c.n,160 coef: b,161 se,162 t,163 p,164 stars: p != null ? stars(p) : "",165 ciLo: se != null ? b - tcrit * se : null,166 ciHi: se != null ? b + tcrit * se : null,167 effectPct,168 effectDollar,169 effectLo: eLo?.pct ?? null,170 effectHi: eHi?.pct ?? null,171 step,172 };173 });174}175176function equationOf(spec: HedonicSpec, cols: Column[]): string {177 const lhs = spec.depvar === "lnprice" ? "ln(Prix)" : "Prix";178 const terms = cols.filter((c) => c.kind === "cont" || c.kind === "trend").map((c) => c.name);179 const cats = [...new Set(cols.filter((c) => c.kind === "dummy").map((c) => c.var))];180 const fe: string[] = [];181 if (spec.fe.time !== "none") fe.push(`δ_${spec.fe.time === "year" ? "année" : spec.fe.time === "quarter" ? "trimestre" : "mois"}`);182 if (spec.fe.geo !== "none") fe.push(`γ_zone`);183 const rho = spec.method === "sar" ? "ρ·W·" + lhs + " + " : "";184 return `${lhs} = ${rho}β₀ + ${terms.map((t, i) => `β${i + 1}·${t}`).join(" + ")}${cats.length ? " + " + cats.map((c) => `Σ β·1[${c}]`).join(" + ") : ""}${fe.length ? " + " + fe.join(" + ") : ""} + ε`;185}186187function codeSnippets(spec: HedonicSpec, cols: Column[]): HedonicResult["code"] {188 const cont = cols.filter((c) => c.kind === "cont").map((c) => (c.var === "age2" ? "I(age^2/100)" : c.log ? `log(${c.var})` : c.var!));189 const trend = cols.some((c) => c.kind === "trend") ? ["poly(x_km, 2)", "poly(y_km, 2)", "x_km:y_km"] : [];190 const cats = [...new Set(cols.filter((c) => c.kind === "dummy").map((c) => c.var!))];191 const fe: string[] = [];192 if (spec.fe.time !== "none") fe.push(spec.fe.time);193 if (spec.fe.geo !== "none") fe.push("zone");194 const y = spec.depvar === "lnprice" ? "log(price)" : "price";195 const rhs = [...cont, ...trend, ...cats].join(" + ") || "1";196 const vcovR = spec.se === "hc1" ? 'vcov = "hetero"' : spec.se === "cluster" ? "vcov = ~zone" : 'vcov = "iid"';197 const r =198 spec.method === "ols" || spec.method === "huber" || spec.method === "lad"199 ? `library(fixest)${spec.method !== "ols" ? "\nlibrary(MASS); library(quantreg)" : ""}200d <- read.csv("echantillon_uqo_eval.csv")201${spec.method === "ols" ? `m <- feols(${y} ~ ${rhs}${fe.length ? " | " + fe.join(" + ") : ""}, data = d, ${vcovR})\netable(m)` : spec.method === "huber" ? `m <- rlm(${y} ~ ${rhs}${fe.map((f) => ` + factor(${f})`).join("")}, data = d, psi = psi.huber)\nsummary(m)` : `m <- rq(${y} ~ ${rhs}${fe.map((f) => ` + factor(${f})`).join("")}, tau = 0.5, data = d)\nsummary(m, se = "boot")`}`202 : spec.method === "sar"203 ? `library(spatialreg); library(spdep)204d <- read.csv("echantillon_uqo_eval.csv")205nb <- knn2nb(knearneigh(cbind(d$lng, d$lat), k = 8)); W <- nb2listw(nb, style = "W")206m <- stsls(${y} ~ ${rhs}${fe.map((f) => ` + factor(${f})`).join("")}, data = d, listw = W) # SAR par 2SLS (Kelejian-Prucha)207summary(m)`208 : `library(glmnet)209d <- read.csv("echantillon_uqo_eval.csv")210X <- model.matrix(~ ${rhs}${fe.map((f) => ` + factor(${f})`).join("")}, d)[, -1]211cv <- cv.glmnet(X, ${y.replace("price", "d$price")}, alpha = ${spec.method === "lasso" ? 1 : 0}, nfolds = 5)212coef(cv, s = "lambda.min")`;213 const pyRhs = [...cols.filter((c) => c.kind === "cont").map((c) => (c.var === "age2" ? "I(age**2/100)" : c.log ? `np.log(${c.var})` : c.var!)), ...cats.map((c) => `C(${c})`), ...fe.map((f) => `C(${f})`)].join(" + ") || "1";214 const pyY = spec.depvar === "lnprice" ? "np.log(price)" : "price";215 const python =216 spec.method === "ridge" || spec.method === "lasso"217 ? `import pandas as pd, numpy as np218from sklearn.linear_model import ${spec.method === "lasso" ? "LassoCV" : "RidgeCV"}219d = pd.read_csv("echantillon_uqo_eval.csv")220X = pd.get_dummies(d[[${[...cols.filter((c) => c.kind === "cont" && c.var !== "age2").map((c) => `"${c.var}"`), ...cats.map((c) => `"${c}"`), ...fe.map((f) => `"${f}"`)].join(", ")}]], drop_first=True)221y = ${pyY.replace("price", "d.price")}222m = ${spec.method === "lasso" ? "LassoCV(cv=5)" : "RidgeCV(alphas=np.logspace(-4, 1, 24))"}.fit((X - X.mean()) / X.std(), y)223print(dict(zip(X.columns, m.coef_)))`224 : `import pandas as pd, numpy as np225import statsmodels.formula.api as smf226d = pd.read_csv("echantillon_uqo_eval.csv")227m = smf.${spec.method === "lad" ? "quantreg" : "ols"}("${pyY} ~ ${pyRhs}", data=d).fit(${spec.method === "lad" ? "q=0.5" : spec.se === "hc1" ? 'cov_type="HC1"' : spec.se === "cluster" ? 'cov_type="cluster", cov_kwds={"groups": d["zone"]}' : ""})228print(m.summary())`;229 const stY = spec.depvar === "lnprice" ? "lnprice" : "price";230 const stX = [...cols.filter((c) => c.kind === "cont").map((c) => (c.var === "age2" ? "c.age#c.age" : c.log ? `ln_${c.var}` : c.var!)), ...cats.map((c) => `i.${c}`)].join(" ");231 const stata =232 spec.method === "ols"233 ? `import delimited echantillon_uqo_eval.csv, clear234${spec.depvar === "lnprice" ? "gen lnprice = ln(price)\n" : ""}${cols.filter((c) => c.kind === "cont" && c.log).map((c) => `gen ln_${c.var} = ln(${c.var})`).join("\n")}235${fe.length ? `reghdfe ${stY} ${stX}, absorb(${fe.join(" ")}) vce(${spec.se === "cluster" ? "cluster zone" : "robust"})` : `regress ${stY} ${stX}, vce(${spec.se === "cluster" ? "cluster zone" : "robust"})`}`236 : spec.method === "lad"237 ? `import delimited echantillon_uqo_eval.csv, clear\nqreg ${stY} ${stX} ${fe.map((f) => `i.${f}`).join(" ")}, quantile(.5)`238 : spec.method === "huber"239 ? `import delimited echantillon_uqo_eval.csv, clear\nrreg ${stY} ${stX} ${fe.map((f) => `i.${f}`).join(" ")}`240 : spec.method === "sar"241 ? `import delimited echantillon_uqo_eval.csv, clear\nspset, modify coordsys(latlong, kilometers)\nspmatrix create idistance W, knn(8) normalize(row)\nspivregress ${stY} ${stX} ${fe.map((f) => `i.${f}`).join(" ")}, dvarlag(W)`242 : `import delimited echantillon_uqo_eval.csv, clear\n${spec.method === "lasso" ? "lasso" : "elasticnet"} linear ${stY} ${stX} ${fe.map((f) => `i.${f}`).join(" ")}, selection(cv, folds(5))${spec.method === "ridge" ? " alpha(0)" : ""}`;243 return { r, python, stata };244}245246/* ------------------------------ estimation ------------------------------ */247interface Est {248 fit: Fit;249 V: Float64Array | null;250 cols: Column[]; // colonnes correspondant à fit.beta (SAR : ρ en tête)251 offsetX: number; // 1 si ρ en tête252 lambda: HedonicResult["model"]["lambda"];253 lassoSelected: string[] | null;254 postLasso: Coefficient[] | null;255 rho: Coefficient | null;256 robust: HedonicResult["model"]["robust"];257 nb: import("./spatial").KnnIndex | null;258 sarRes: ReturnType<typeof sar> | null;259 notes: Txt[];260}261262function estimateAll(design: Design, spec: HedonicSpec, medP: number): Est {263 const { X, y, cols, clusters, nClusters } = design;264 const logModel = spec.depvar === "lnprice";265 const notes: Txt[] = [];266 const lat = Float64Array.from(design.obs, (o) => o.lat);267 const lng = Float64Array.from(design.obs, (o) => o.lng);268 if (spec.method === "ols" || spec.method === "huber") {269 const fit = spec.method === "ols" ? ols(X, y) : huber(X, y);270 const V = covariance(fit, spec.se, clusters, nClusters);271 if (spec.method === "huber") notes.push(T("Huber : écarts-types de type sandwich sur la dernière itération pondérée (approximation usuelle).", "Huber: sandwich-type standard errors on the last weighted iteration (usual approximation)."));272 return { fit, V, cols, offsetX: 0, lambda: null, lassoSelected: null, postLasso: null, rho: null, robust: spec.method === "huber" ? { scale: fit.extra.scale, downweightedPct: fit.extra.downweightedPct } : null, nb: null, sarRes: null, notes };273 }274 if (spec.method === "lad") {275 const fit = lad(X, y);276 notes.push(T("Régression médiane : coefficients par IRLS ; les écarts-types exigeraient un bootstrap (non calculé ici).", "Median regression: IRLS coefficients; standard errors would require a bootstrap (not computed here)."));277 return { fit, V: null, cols, offsetX: 0, lambda: null, lassoSelected: null, postLasso: null, rho: null, robust: null, nb: null, sarRes: null, notes };278 }279 if (spec.method === "ridge" || spec.method === "lasso") {280 const res = penalized(X, y, spec.method);281 notes.push(T(`λ choisi par validation croisée à 5 blocs sur une grille de ${res.grid.length} valeurs (échelle standardisée).`, `λ chosen by 5-fold cross-validation over a ${res.grid.length}-value grid (standardised scale).`));282 let postLasso: Coefficient[] | null = null;283 let selected: string[] | null = null;284 if (spec.method === "lasso") {285 selected = res.selected.filter((j) => j > 0).map((j) => cols[j].name);286 if (res.postFit) {287 const pcols = res.selected.map((j) => cols[j]);288 const Vp = covariance(res.postFit, spec.se, clusters, nClusters);289 postLasso = coefTable(pcols, res.postFit.beta, Vp, res.postFit.df, logModel, medP);290 }291 notes.push(T("Coefficients LASSO biaisés vers zéro ; l'inférence post-LASSO (MCO sur les variables retenues) est indicative.", "LASSO coefficients are shrunk toward zero; post-LASSO inference (OLS on kept variables) is indicative."));292 } else notes.push(T("Ridge : pas d'écarts-types (estimateur biaisé) ; lire les coefficients comme des effets stabilisés.", "Ridge: no standard errors (biased estimator); read coefficients as stabilised effects."));293 return { fit: res.fit, V: null, cols, offsetX: 0, lambda: { value: res.lambda, grid: res.grid }, lassoSelected: selected, postLasso, rho: null, robust: null, nb: null, sarRes: null, notes };294 }295 // SAR296 if (X.n > LIMITS.maxN) throw new Error("SAR : échantillon trop grand");297 const res = sar(X, y, lat, lng, 8);298 const V = covariance(res.fit, spec.se === "cluster" ? "cluster" : spec.se, clusters, nClusters);299 const k = res.fit.k;300 const se = V ? Math.sqrt(Math.max(0, V[0])) : null;301 const t = se ? res.rho / se : null;302 const p = t != null ? tPValue(t, res.fit.df) : null;303 const rho: Coefficient = { name: "rho", label: T("ρ — autocorrélation spatiale (W·y, 8 voisins)", "ρ — spatial autocorrelation (W·y, 8 neighbours)"), kind: "rho", coef: res.rho, se, t, p, stars: p != null ? stars(p) : "", ciLo: se != null ? res.rho - 1.96 * se : null, ciHi: se != null ? res.rho + 1.96 * se : null, effectPct: null, effectDollar: null, effectLo: null, effectHi: null, step: T("—", "—") };304 void k;305 notes.push(T("SAR estimé par doubles moindres carrés (instruments WX, W²X) : les β sont des effets directs ; l'effet total ≈ β/(1−ρ).", "SAR estimated by two-stage least squares (WX, W²X instruments): β are direct effects; total effect ≈ β/(1−ρ)."));306 return { fit: res.fit, V, cols, offsetX: 1, lambda: null, lassoSelected: null, postLasso: null, rho, robust: null, nb: res.index, sarRes: res, notes };307}308309/* ------------------------------ validation croisée ------------------------------ */310function crossValidate(design: Design, spec: HedonicSpec, folds: number, medP: number): FitMetrics {311 const { X, y } = design;312 const n = X.n;313 const logModel = spec.depvar === "lnprice";314 const prices = Float64Array.from(design.obs, (o) => o.price);315 const lat = Float64Array.from(design.obs, (o) => o.lat);316 const lng = Float64Array.from(design.obs, (o) => o.lng);317 const predAll = new Float64Array(n);318 const fittedAll = new Float64Array(n);319 for (let f = 0; f < folds; f++) {320 const tr: number[] = [];321 const te: number[] = [];322 for (let i = 0; i < n; i++) (i % folds === f ? te : tr).push(i);323 const Xtr = takeRows(X, tr);324 const ytr = takeVec(y, tr);325 const Xte = takeRows(X, te);326 let predTe: Float64Array;327 let smear = 1;328 if (spec.method === "sar") {329 const r = sar(Xtr, ytr, takeVec(lat, tr), takeVec(lng, tr), 8);330 smear = smearing(r.fit.resid, logModel);331 predTe = sarPredict(r, Xte, takeVec(lat, te), takeVec(lng, te), ytr, r.fit.beta);332 } else {333 let beta: Float64Array;334 let resid: Float64Array;335 if (spec.method === "ols") ({ beta, resid } = ols(Xtr, ytr));336 else if (spec.method === "huber") ({ beta, resid } = huber(Xtr, ytr, 15));337 else if (spec.method === "lad") ({ beta, resid } = lad(Xtr, ytr, 25));338 else ({ beta, resid } = penalized(Xtr, ytr, spec.method, 3).fit);339 smear = smearing(resid, logModel);340 predTe = new Float64Array(te.length);341 for (let r = 0; r < te.length; r++) {342 let s = 0;343 for (let j = 0; j < Xte.k; j++) s += Xte.a[r * Xte.k + j] * beta[j];344 predTe[r] = s;345 }346 }347 for (let r = 0; r < te.length; r++) {348 fittedAll[te[r]] = predTe[r];349 predAll[te[r]] = logModel ? Math.exp(predTe[r]) * smear : predTe[r];350 }351 }352 const resid = new Float64Array(n);353 for (let i = 0; i < n; i++) resid[i] = y[i] - fittedAll[i];354 void medP;355 // métriques hors échantillon : on passe des prédictions déjà en $ via smear=1 et fitted=ln(pred)356 const fittedForMetrics = logModel ? Float64Array.from(predAll, Math.log) : predAll;357 return fitMetrics(y, fittedForMetrics, resid, prices, logModel, X.k, 1);358}359360/* --------------------------------- rapport --------------------------------- */361export function runHedonic(spec: HedonicSpec): HedonicResult {362 const t0 = Date.now();363 const sample = buildSample(spec);364 if (sample.obs.length < LIMITS.minN) throw new Error(`échantillon insuffisant : ${sample.obs.length} observations (min ${LIMITS.minN}). Élargissez la zone, la période ou les types.`);365 const design = buildDesign(sample.obs, spec);366 const { X, y, cols, obs } = design;367 const n = X.n;368 const logModel = spec.depvar === "lnprice";369 const prices = Float64Array.from(obs, (o) => o.price);370 const medP = median(prices);371372 const est = estimateAll(design, spec, medP);373 const fit = est.fit;374 const coefficients = coefTable(cols, fit.beta, est.V, fit.df, logModel, medP, est.offsetX);375 const metrics = fitMetrics(y, fit.fitted, fit.resid, prices, logModel, fit.k);376377 // F global (MCO/Huber seulement)378 let fstat: HedonicResult["model"]["fstat"] = null;379 if ((spec.method === "ols" || spec.method === "huber") && fit.k > 1) {380 const df1 = fit.k - 1;381 const df2 = fit.df;382 const f = (metrics.r2 / df1) / ((1 - metrics.r2) / df2);383 fstat = { f, df1, df2, p: fSf(f, df1, df2) };384 }385386 // diagnostics387 let bp: HedonicResult["diagnostics"]["bp"] = null;388 try {389 bp = breuschPagan(X, fit.resid);390 } catch {}391 const jb = jarqueBera(fit.resid);392 const lat = Float64Array.from(obs, (o) => o.lat);393 const lng = Float64Array.from(obs, (o) => o.lng);394 let moran: Moran | null = null;395 try {396 const index = est.nb ?? new KnnIndex(lat, lng);397 const nb = est.sarRes ? est.sarRes.nb : index.selfNeighbours(8);398 moran = moranI(nb, fit.resid);399 } catch {}400 let vifs: HedonicResult["diagnostics"]["vif"] = [];401 try {402 const contCols = cols.map((c, j) => ({ c, j })).filter(({ c }) => c.kind === "cont" || c.kind === "trend");403 const olsFit = spec.method === "ols" ? fit : ols(X, y);404 if (olsFit.bread && contCols.length) {405 const v = vif(X, olsFit.bread, contCols.map(({ j }) => j));406 vifs = contCols.map(({ c }, i) => ({ name: c.name, label: c.label, vif: v[i] }));407 }408 } catch {}409410 // validation croisée411 let cv: HedonicResult["cv"] = null;412 if (spec.cv) {413 try {414 const m = crossValidate(design, spec, 5, medP);415 cv = { folds: 5, metrics: m, gapMdape: m.mdape - metrics.mdape };416 } catch {}417 }418419 // comparaison au laboratoire420 let lab: HedonicResult["lab"] = null;421 if (spec.compareLab) {422 const idx: number[] = [];423 for (let i = 0; i < n; i++) if (obs[i].lab && obs[i].lab! > 0) idx.push(i);424 if (idx.length >= 30) {425 const pred = toDollars(fit.fitted, logModel, metrics.smear);426 const m = (get: (i: number) => number) => {427 const ape: number[] = [];428 let se = 0;429 let w10 = 0;430 let w20 = 0;431 for (const i of idx) {432 const a = Math.abs(get(i) - prices[i]) / prices[i];433 ape.push(a);434 se += (get(i) - prices[i]) ** 2;435 if (a <= 0.1) w10++;436 if (a <= 0.2) w20++;437 }438 return { mdape: median(ape), rmse: Math.sqrt(se / idx.length), within10: w10 / idx.length, within20: w20 / idx.length };439 };440 lab = {441 n: idx.length,442 lab: m((i) => obs[i].lab!),443 ours: m((i) => pred[i]),444 label: spec.source === "sales" ? T("LightGBM du laboratoire (est_2026 ramenée à la date de vente par l'indice)", "Lab LightGBM (est_2026 deflated to the sale date with the index)") : T("Mesure UQO Éval de l'annonce (hybride 65/35)", "UQO Éval measure of the listing (65/35 hybrid)"),445 };446 }447 }448449 // indice temporel et effets géo450 const V = est.V;451 const kk = fit.k;452 const idxOf = (name: string) => cols.findIndex((c) => c.name === name);453 let timeIndex: HedonicResult["timeIndex"] = null;454 if (spec.fe.time !== "none" && design.timeLevels.length) {455 timeIndex = design.timeLevels.map((l) => {456 if (l.base) return { level: l.level, n: l.n, coef: 0, index: 100, lo: null, hi: null };457 const j = idxOf(`t=${l.level}`) + est.offsetX;458 const b = fit.beta[j];459 const se = V ? Math.sqrt(Math.max(0, V[j * kk + j])) : null;460 const toIdx = (x: number) => (logModel ? 100 * Math.exp(x) : 100 * (1 + x / medP));461 return { level: l.level, n: l.n, coef: b, index: toIdx(b), lo: se != null ? toIdx(b - 1.96 * se) : null, hi: se != null ? toIdx(b + 1.96 * se) : null };462 });463 }464 let geoEffects: HedonicResult["geoEffects"] = null;465 if (spec.fe.geo !== "none" && design.geoLevels.length) {466 geoEffects = design.geoLevels467 .map((l) => {468 if (l.base) return { level: l.level, n: l.n, coef: 0, pct: 0, se: null };469 const j = idxOf(`g=${l.level}`) + est.offsetX;470 const b = fit.beta[j];471 const se = V ? Math.sqrt(Math.max(0, V[j * kk + j])) : null;472 return { level: l.level, n: l.n, coef: b, pct: logModel ? 100 * (Math.exp(b) - 1) : (100 * b) / medP, se };473 })474 .sort((a, b) => b.pct - a.pct);475 }476477 // graphiques478 const stride = Math.max(1, Math.floor(n / 1200));479 const residFitted: { f: number; e: number }[] = [];480 const actualPred: { a: number; p: number }[] = [];481 const predD = toDollars(fit.fitted, logModel, metrics.smear);482 for (let i = 0; i < n; i += stride) {483 residFitted.push({ f: fit.fitted[i], e: fit.resid[i] / metrics.sigma });484 actualPred.push({ a: prices[i], p: predD[i] });485 }486 const strideMap = Math.max(1, Math.floor(n / 2500));487 const map: { lat: number; lng: number; e: number }[] = [];488 for (let i = 0; i < n; i += strideMap) map.push({ lat: obs[i].lat, lng: obs[i].lng, e: fit.resid[i] / metrics.sigma });489490 // insights491 // typographie française : virgule décimale492 const insights = buildInsights(spec, design, coefficients, metrics, cv, lab, bp, jb, moran, vifs, timeIndex, geoEffects, est, medP).map((s) => ({493 fr: s.fr.replace(/(\d)\.(\d)/g, "$1,$2"),494 en: s.en,495 }));496497 const byTypeMap = new Map<string, number>();498 for (const o of obs) byTypeMap.set(o.type, (byTypeMap.get(o.type) ?? 0) + 1);499 const dates = obs.map((o) => o.date).sort();500 const methodDef = METHODS.find((m) => m.key === spec.method)!;501502 return {503 spec,504 sample: {505 n,506 available: sample.available,507 trimmed: sample.trimmed,508 dropped: design.dropped,509 period: { from: dates[0], to: dates[dates.length - 1] },510 byType: [...byTypeMap.entries()].map(([type, n2]) => ({ type, n: n2 })).sort((a, b) => b.n - a.n),511 descriptives: descriptives(design),512 medianPrice: medP,513 center: design.center,514 },515 model: {516 method: spec.method,517 methodLabel: T(methodDef.fr, methodDef.en),518 depvar: spec.depvar,519 equation: equationOf(spec, cols),520 n,521 k: fit.k,522 df: fit.df,523 seType: spec.se,524 nClusters: design.nClusters,525 coefficients,526 fit: metrics,527 fstat,528 lambda: est.lambda,529 lassoSelected: est.lassoSelected,530 postLasso: est.postLasso,531 rho: est.rho,532 robust: est.robust,533 notes: est.notes,534 },535 diagnostics: { bp, jb, moran, vif: vifs },536 cv,537 lab,538 timeIndex,539 geoEffects,540 plots: { residFitted, qq: qqPoints(fit.resid, metrics.sigma), hist: histogram(Float64Array.from(fit.resid, (e) => e / metrics.sigma)), map, actualPred },541 insights,542 code: codeSnippets(spec, cols),543 timingMs: Date.now() - t0,544 };545}546547/* --------------------------------- lecture --------------------------------- */548function buildInsights(549 spec: HedonicSpec,550 design: Design,551 coefs: Coefficient[],552 m: FitMetrics,553 cv: HedonicResult["cv"],554 lab: HedonicResult["lab"],555 bp: HedonicResult["diagnostics"]["bp"],556 jb: HedonicResult["diagnostics"]["jb"],557 moran: Moran | null,558 vifs: HedonicResult["diagnostics"]["vif"],559 timeIndex: HedonicResult["timeIndex"],560 geo: HedonicResult["geoEffects"],561 est: Est,562 medP: number563): Txt[] {564 const out: Txt[] = [];565 const pc = (x: number) => (100 * x).toFixed(1);566 // 1. qualité567 out.push(568 T(569 `Le modèle explique ${pc(m.r2)} % de la variance ${spec.depvar === "lnprice" ? "du log du prix" : "du prix"} (R² ajusté ${pc(m.adjR2)} %) sur ${m.n.toLocaleString("fr-CA")} observations. Erreur médiane absolue : ${pc(m.mdape)} % ; ${pc(m.within10)} % des prédictions à ±10 %, ${pc(m.within20)} % à ±20 %.`,570 `The model explains ${pc(m.r2)}% of the variance ${spec.depvar === "lnprice" ? "of log price" : "of price"} (adjusted R² ${pc(m.adjR2)}%) over ${m.n.toLocaleString("en-CA")} observations. Median absolute error: ${pc(m.mdape)}%; ${pc(m.within10)}% of predictions within ±10%, ${pc(m.within20)}% within ±20%.`571 )572 );573 if (cv) {574 const gap = cv.gapMdape * 100;575 out.push(576 T(577 `Hors échantillon (validation croisée à ${cv.folds} blocs), l'erreur médiane passe à ${pc(cv.metrics.mdape)} % (${gap >= 0 ? "+" : "−"}${Math.abs(gap).toFixed(1)} pt) — ${gap > 3 ? "écart notable : le modèle sur-apprend (trop d'effets fixes fins pour l'échantillon ?)" : "écart faible : le modèle généralise bien"}.`,578 `Out of sample (${cv.folds}-fold cross-validation), the median error moves to ${pc(cv.metrics.mdape)}% (${gap >= 0 ? "+" : "−"}${Math.abs(gap).toFixed(1)} pt) — ${gap > 3 ? "notable gap: the model overfits (too many fine fixed effects for the sample?)" : "small gap: the model generalises well"}.`579 )580 );581 }582 // 2. effets clés583 const key = (v: string) => coefs.find((c) => c.var === v && c.kind === "cont");584 const a = key("aire");585 if (a && a.effectPct != null)586 out.push(587 a.step.fr === "pour +10 %"588 ? T(`Élasticité prix-superficie : ${a.coef.toFixed(3)} — +10 % d'aire d'étages ⇒ ${fmtPct(a.effectPct)} du prix (≈ ${fmtMoney(a.effectDollar!)} au prix médian de ${fmtMoney(medP)})${a.p != null ? `, p ${a.p < 0.001 ? "< 0,001" : "= " + a.p.toFixed(3)}` : ""}.`, `Price-area elasticity: ${a.coef.toFixed(3)} — +10% floor area ⇒ ${fmtPct(a.effectPct)} of price (≈ ${fmtMoneyEn(a.effectDollar!)} at the median price of ${fmtMoneyEn(medP)})${a.p != null ? `, p ${a.p < 0.001 ? "< 0.001" : "= " + a.p.toFixed(3)}` : ""}.`)589 : T(`Chaque m² d'aire d'étages supplémentaire vaut ${fmtPct(a.effectPct, 2)} du prix, soit ≈ ${fmtMoney(a.effectDollar!)} au prix médian.`, `Each extra m² of floor area is worth ${fmtPct(a.effectPct, 2)} of price, ≈ ${fmtMoneyEn(a.effectDollar!)} at the median price.`)590 );591 const te = key("terrain");592 if (te && te.effectPct != null) out.push(T(`Terrain : ${te.step.fr === "pour +10 %" ? `élasticité ${te.coef.toFixed(3)} (+10 % ⇒ ${fmtPct(te.effectPct)})` : `${fmtMoney(te.effectDollar!)} par m²`} — ${Math.abs(te.coef) < (te.step.fr === "pour +10 %" ? 0.3 : 50) ? "rendement décroissant marqué, typique de la rente foncière résidentielle" : "contribution forte du sol"}.`, `Lot: ${te.step.en === "for +10%" ? `elasticity ${te.coef.toFixed(3)} (+10% ⇒ ${fmtPct(te.effectPct)})` : `${fmtMoneyEn(te.effectDollar!)} per m²`} — ${Math.abs(te.coef) < (te.step.en === "for +10%" ? 0.3 : 50) ? "strong diminishing returns, typical of residential land rent" : "strong land contribution"}.`));593 const ag = key("age");594 const ag2 = coefs.find((c) => c.var === "age2");595 if (ag && ag.effectPct != null) {596 if (ag2) {597 // sommet de la parabole : d/dage = b1 + 2·b2·age/100 = 0 → age* = −50·b1/b2598 const turn = ag2.coef !== 0 ? (-50 * ag.coef) / ag2.coef : NaN;599 out.push(T(`Âge : ${fmtPct(ag.effectPct, 2)} par année au départ, courbure ${ag2.coef >= 0 ? "convexe" : "concave"}${Number.isFinite(turn) && turn > 0 && turn < 200 ? ` — la dépréciation ${ag2.coef > 0 ? "s'inverse (effet vintage) vers" : "s'accélère au-delà de"} ${turn.toFixed(0)} ans` : ""}.`, `Age: ${fmtPct(ag.effectPct, 2)} per year initially, ${ag2.coef >= 0 ? "convex" : "concave"} curvature${Number.isFinite(turn) && turn > 0 && turn < 200 ? ` — depreciation ${ag2.coef > 0 ? "reverses (vintage effect) around" : "accelerates beyond"} ${turn.toFixed(0)} years` : ""}.`));600 } else out.push(T(`Dépréciation : ${fmtPct(ag.effectPct, 2)} du prix par année d'âge (≈ ${fmtMoney(ag.effectDollar!)} au prix médian).`, `Depreciation: ${fmtPct(ag.effectPct, 2)} of price per year of age (≈ ${fmtMoneyEn(ag.effectDollar!)} at the median price).`));601 }602 const dc = key("dist_centre");603 if (dc && dc.effectPct != null) out.push(T(`Gradient de rente : ${fmtPct(dc.effectPct, 2)} du prix ${dc.step.fr} de distance au centre${dc.effectPct < 0 ? " — la centralité est valorisée" : " — la périphérie est ici valorisée (grands terrains, secteurs cossus ?)"}.`, `Rent gradient: ${fmtPct(dc.effectPct, 2)} of price ${dc.step.en} of distance to the centre${dc.effectPct < 0 ? " — centrality is valued" : " — the periphery is valued here (large lots, affluent sectors?)"}.`));604 const dummies = coefs.filter((c) => c.kind === "dummy" && c.p != null && c.p < 0.05 && c.effectPct != null).sort((x, y2) => Math.abs(y2.effectPct!) - Math.abs(x.effectPct!));605 if (dummies.length) {606 const d = dummies[0];607 out.push(T(`Attribut le plus discriminant : « ${d.label.fr} » ⇒ ${fmtPct(d.effectPct!)} vs ${d.base} (significatif à 5 %).`, `Most discriminating attribute: “${d.label.en}” ⇒ ${fmtPct(d.effectPct!)} vs ${d.base} (significant at 5%).`));608 }609 // 3. temps / géo610 if (timeIndex && timeIndex.length > 1) {611 const last = timeIndex[timeIndex.length - 1];612 const peak = timeIndex.reduce((p, c) => (c.index > p.index ? c : p), timeIndex[0]);613 out.push(T(`Indice hédonique (qualité constante) : ${last.index.toFixed(1)} en ${last.level} (base 100 = ${timeIndex[0].level}), sommet ${peak.index.toFixed(1)} en ${peak.level}. C'est l'appréciation nette des changements de composition des ventes.`, `Hedonic (constant-quality) index: ${last.index.toFixed(1)} in ${last.level} (base 100 = ${timeIndex[0].level}), peak ${peak.index.toFixed(1)} in ${peak.level}. This is appreciation net of changes in the sales mix.`));614 }615 if (geo && geo.length > 2) {616 const top = geo[0];617 const bottom = geo[geo.length - 1];618 out.push(T(`Prime de localisation : ${top.level} ${fmtPct(top.pct)} vs ${bottom.level} ${fmtPct(bottom.pct)} (référence : ${design.geoLevels.find((l) => l.base)?.level}) — l'écart de ${(top.pct - bottom.pct).toFixed(0)} pts est l'effet « quartier » à caractéristiques égales.`, `Location premium: ${top.level} ${fmtPct(top.pct)} vs ${bottom.level} ${fmtPct(bottom.pct)} (base: ${design.geoLevels.find((l) => l.base)?.level}) — the ${(top.pct - bottom.pct).toFixed(0)}-pt gap is the “neighbourhood” effect holding characteristics constant.`));619 }620 // 4. diagnostics621 if (bp) out.push(bp.p < 0.05 ? T(`Hétéroscédasticité détectée (Breusch-Pagan p ${bp.p < 0.001 ? "< 0,001" : "= " + bp.p.toFixed(3)}) : ${spec.se === "classic" ? "passez aux écarts-types robustes (HC1) — les p-valeurs classiques sont trompeuses" : "vos écarts-types robustes sont donc de mise"}.`, `Heteroskedasticity detected (Breusch-Pagan p ${bp.p < 0.001 ? "< 0.001" : "= " + bp.p.toFixed(3)}): ${spec.se === "classic" ? "switch to robust (HC1) standard errors — classical p-values are misleading" : "your robust standard errors are therefore warranted"}.`) : T(`Pas d'hétéroscédasticité détectable (Breusch-Pagan p = ${bp.p.toFixed(3)}).`, `No detectable heteroskedasticity (Breusch-Pagan p = ${bp.p.toFixed(3)}).`));622 out.push(jb.p < 0.05 ? T(`Résidus non normaux (Jarque-Bera, asymétrie ${jb.skew.toFixed(2)}, aplatissement ${jb.kurt.toFixed(1)}) — sans conséquence sur les coefficients en grand échantillon, mais ${jb.kurt > 5 ? "les queues épaisses suggèrent des ventes atypiques : essayez la régression robuste ou quantile" : "les intervalles de prédiction individuels sont approximatifs"}.`, `Non-normal residuals (Jarque-Bera, skewness ${jb.skew.toFixed(2)}, kurtosis ${jb.kurt.toFixed(1)}) — harmless for coefficients in large samples, but ${jb.kurt > 5 ? "fat tails suggest atypical sales: try robust or quantile regression" : "individual prediction intervals are approximate"}.`) : T("Résidus compatibles avec la normalité (Jarque-Bera).", "Residuals compatible with normality (Jarque-Bera)."));623 if (moran && Number.isFinite(moran.z)) out.push(moran.p < 0.05 ? T(`Autocorrélation spatiale des résidus (I de Moran = ${moran.I.toFixed(3)}, z = ${moran.z.toFixed(1)}) : les erreurs voisines se ressemblent — il manque de la localisation au modèle (${spec.method === "sar" ? "malgré le terme ρWy" : spec.fe.geo === "none" ? "ajoutez des effets fixes géo ou passez au SAR" : "affinez la grille ou passez au SAR"}).`, `Spatial autocorrelation of residuals (Moran's I = ${moran.I.toFixed(3)}, z = ${moran.z.toFixed(1)}): neighbouring errors look alike — the model lacks location (${spec.method === "sar" ? "despite the ρWy term" : spec.fe.geo === "none" ? "add geo fixed effects or switch to SAR" : "refine the grid or switch to SAR"}).`) : T(`Pas d'autocorrélation spatiale résiduelle notable (I de Moran = ${moran.I.toFixed(3)}, p = ${moran.p.toFixed(2)}) — la localisation est bien captée.`, `No notable residual spatial autocorrelation (Moran's I = ${moran.I.toFixed(3)}, p = ${moran.p.toFixed(2)}) — location is well captured.`));624 const bad = vifs.filter((v) => v.vif > 10 && v.name !== "age2" && v.name !== "age" && !v.name.startsWith("trend_"));625 const mech = vifs.some((v) => v.vif > 10 && (v.name === "age2" || v.name === "age" || v.name.startsWith("trend_")));626 if (mech && !bad.length) out.push(T("Les VIF élevés de l'âge² ou de la surface de tendance sont mécaniques (termes polynomiaux) — sans gravité pour l'interprétation des effets marginaux.", "High VIFs for age² or the trend surface are mechanical (polynomial terms) — harmless for interpreting marginal effects."));627 if (bad.length) out.push(T(`Colinéarité : VIF > 10 pour ${bad.map((b) => b.label.fr).join(", ")} — coefficients instables ; ${spec.method === "ridge" ? "le ridge l'atténue" : "envisagez le ridge ou retirez une variable redondante"}.`, `Collinearity: VIF > 10 for ${bad.map((b) => b.label.en).join(", ")} — unstable coefficients; ${spec.method === "ridge" ? "ridge mitigates it" : "consider ridge or drop a redundant variable"}.`));628 // 5. méthode629 if (est.rho) out.push(T(`ρ = ${est.rho.coef.toFixed(3)}${est.rho.p != null ? ` (p ${est.rho.p < 0.001 ? "< 0,001" : "= " + est.rho.p.toFixed(3)})` : ""} : ${Math.abs(est.rho.coef) > 0.2 ? "forte dépendance spatiale — le prix des 8 voisins « explique » une part importante du prix ; multiplicateur spatial ≈ " + (1 / (1 - est.rho.coef)).toFixed(2) : "dépendance spatiale modérée"}.`, `ρ = ${est.rho.coef.toFixed(3)}${est.rho.p != null ? ` (p ${est.rho.p < 0.001 ? "< 0.001" : "= " + est.rho.p.toFixed(3)})` : ""}: ${Math.abs(est.rho.coef) > 0.2 ? "strong spatial dependence — neighbours' prices “explain” a large share of price; spatial multiplier ≈ " + (1 / (1 - est.rho.coef)).toFixed(2) : "moderate spatial dependence"}.`));630 if (est.robust) out.push(T(`Huber : ${est.robust.downweightedPct.toFixed(1)} % des observations sous-pondérées (résidus > 1,345·σ̂).`, `Huber: ${est.robust.downweightedPct.toFixed(1)}% of observations downweighted (residuals > 1.345·σ̂).`));631 if (est.lassoSelected) out.push(T(`LASSO : ${est.lassoSelected.length} variables retenues sur ${coefs.length - 1} (λ = ${est.lambda!.value.toExponential(2)}).`, `LASSO: ${est.lassoSelected.length} variables kept out of ${coefs.length - 1} (λ = ${est.lambda!.value.toExponential(2)}).`));632 // 6. laboratoire633 if (lab) {634 const better = lab.ours.mdape < lab.lab.mdape;635 out.push(T(`Face au laboratoire (${lab.n.toLocaleString("fr-CA")} observations communes) : votre modèle ${pc(lab.ours.mdape)} % d'erreur médiane vs ${pc(lab.lab.mdape)} % pour ${lab.label.fr} — ${better ? "vous faites mieux en échantillon (attention : le laboratoire est évalué hors échantillon, comparez avec votre validation croisée)" : "le laboratoire reste devant : il exploite 690 000 ventes et des interactions non linéaires (arbres)"}.`, `Against the lab (${lab.n.toLocaleString("en-CA")} common observations): your model ${pc(lab.ours.mdape)}% median error vs ${pc(lab.lab.mdape)}% for ${lab.label.en} — ${better ? "you do better in-sample (caution: the lab is evaluated out of sample; compare with your cross-validation)" : "the lab stays ahead: it uses 690,000 sales and non-linear interactions (trees)"}.`));636 }637 return out;638}639640export type { Obs };641export { hcat, lag, type Mat };642