Python 57%
TeX 42.9%
1#!/usr/bin/env python32# -*- coding: utf-8 -*-3"""Figures + tableau d'ablation détaillé pour le papier v5."""4import numpy as np, pandas as pd, warnings, json5import matplotlib6matplotlib.use("Agg")7import matplotlib.pyplot as plt8from sklearn.linear_model import LinearRegression9from sklearn.ensemble import HistGradientBoostingRegressor10from sklearn.model_selection import KFold, cross_val_predict11from sklearn.metrics import r2_score12from sklearn.inspection import permutation_importance13warnings.filterwarnings("ignore")14plt.rcParams.update({"font.family":"serif","font.size":10,"axes.grid":True,15 "grid.alpha":0.25,"axes.axisbelow":True,"figure.dpi":150,"savefig.bbox":"tight"})16C_OLS,C_GBM,C_BL="#3B6EA5","#C4462E","#4A8C5B"1718# build df via v5_ml prefix19exec(open("v5_ml.py").read().split("kf=KFold")[0])20kf=KFold(n_splits=5,shuffle=True,random_state=42)21onehot=pd.get_dummies(df[cats],dtype=float).reset_index(drop=True)2223struct=["sqft_final","has_sqft","lot_sqft","beds","baths","parking","storeys","age","age2"]24amen_=list(AMEN.keys())25S0=struct+amen_26S1=S0+rv+["n_images_c","virtual_tour"]27S2=S1+gv28S3=S2+vv29sets=[("Structural\n+ amenity",S0),("+ Rooms\n+ listing",S1),("+ Text\n(LLM)",S2),("+ Vision\n(façade)",S3)]3031def ols_r2(cols):32 Xi=df[cols].copy()33 for c in cols: Xi[c]=Xi[c].fillna(df.groupby("ptype")[c].transform("median")).fillna(df[c].median())34 X=pd.concat([Xi.reset_index(drop=True),onehot],axis=1)35 p=cross_val_predict(LinearRegression(),X,y,cv=kf,n_jobs=-1);return r2_score(y,p),p36def gbm_r2(cols):37 X=pd.concat([df[cols].reset_index(drop=True),onehot],axis=1)38 m=HistGradientBoostingRegressor(max_iter=600,learning_rate=0.05,max_leaf_nodes=63,l2_regularization=1.0,random_state=42)39 p=cross_val_predict(m,X,y,cv=kf,n_jobs=-1);return r2_score(y,p),p4041ols=[];gbm=[];gpred=None42for nm,cols in sets:43 ro,_=ols_r2(cols);rg,pg=gbm_r2(cols);ols.append(ro);gbm.append(rg)44 if nm.startswith("+ Vision"): gpred=pg45 print(f"{nm.replace(chr(10),' '):<22} OLS={ro:.4f} GBM={rg:.4f}")4647# blend numbers (déjà calculés) : S2 texte, S3 texte+vision48blend={2:0.760,3:0.779}4950# ===== Figure 1 : ablation =====51x=np.arange(len(sets))52fig,ax=plt.subplots(figsize=(6.4,3.9))53ax.plot(x,ols,"-o",color=C_OLS,lw=2,label="OLS (linear)")54ax.plot(x,gbm,"-s",color=C_GBM,lw=2,label="Gradient boosting")55ax.plot([2,3],[blend[2],blend[3]],"-D",color=C_BL,lw=2,label="Blend (NN+GBM)")56for xi,(ro,rg) in enumerate(zip(ols,gbm)):57 ax.annotate(f"{rg:.3f}",(xi,rg),textcoords="offset points",xytext=(0,7),ha="center",fontsize=8,color=C_GBM)58 ax.annotate(f"{ro:.3f}",(xi,ro),textcoords="offset points",xytext=(0,-13),ha="center",fontsize=8,color=C_OLS)59ax.annotate(f"{blend[3]:.3f}",(3,blend[3]),textcoords="offset points",xytext=(0,7),ha="center",fontsize=8,color=C_BL)60ax.set_xticks(x);ax.set_xticklabels([s[0] for s in sets])61ax.set_ylabel(r"Out-of-sample $R^2$ (5-fold CV)");ax.set_ylim(0.55,0.80)62ax.legend(frameon=False,loc="lower right",fontsize=9)63ax.set_title("Incremental gain from each feature layer",fontsize=11)64fig.savefig("paper/figures/fig_ablation.pdf");plt.close(fig)65print("-> fig_ablation.pdf")6667# ===== Figure 2 : calibration + résidus (blend v5) =====68try: bl=np.load("data/processed/v5_blend_oos.npy")69except Exception: bl=gpred70yp=np.exp(bl);ya=np.exp(y)71fig,axs=plt.subplots(1,2,figsize=(7.6,3.6))72ax=axs[0]73hb=ax.hexbin(ya/1e3,yp/1e3,gridsize=45,bins="log",cmap="Blues",xscale="log",yscale="log",mincnt=1)74lim=[60,6000];ax.plot(lim,lim,"--",color="black",lw=1)75ax.set_xlim(lim);ax.set_ylim(lim);ax.set_xlabel("Actual price (\\$000)");ax.set_ylabel("Predicted price (\\$000)")76ax.set_title(f"Calibration (blend, $R^2$={r2_score(y,bl):.3f})",fontsize=10)77r=(yp/ya-1)*10078ax2=axs[1]79ax2.hist(np.clip(r,-60,60),bins=50,color=C_BL,alpha=0.85,edgecolor="white",lw=0.3)80ax2.axvline(0,color="black",lw=1,ls="--")81ax2.set_xlabel("Prediction error (\\%)");ax2.set_ylabel("Listings")82ax2.set_title(f"Errors (median |.|={np.median(np.abs(r)):.1f}\\%)",fontsize=10)83fig.savefig("paper/figures/fig_calibration.pdf");plt.close(fig)84print("-> fig_calibration.pdf")8586# ===== Figure 3 : importance par bloc (GBM v5) =====87Xfull=pd.concat([df[S3].reset_index(drop=True),onehot],axis=1)88gb=HistGradientBoostingRegressor(max_iter=600,learning_rate=0.05,max_leaf_nodes=63,l2_regularization=1.0,random_state=42).fit(Xfull,y)89idx=np.random.RandomState(0).choice(len(Xfull),4000,replace=False)90pi=permutation_importance(gb,Xfull.iloc[idx],y[idx],n_repeats=3,random_state=0,n_jobs=-1)91imp=pd.Series(pi.importances_mean,index=Xfull.columns).clip(lower=0)92def block(c):93 if c.startswith("city_fe") or c.startswith("ptype_fe") or c.startswith("source"): return "Location/type FE"94 if c in ("sqft_final","has_sqft","lot_sqft"): return "Size"95 if c in ("beds","baths","parking","storeys","age","age2"): return "Structural"96 if c in amen_: return "Amenities (rules)"97 if c in rv or c in ("n_images_c","virtual_tour"): return "Rooms / listing"98 if c in gv: return "Text (LLM)"99 if c in vv: return "Vision (façade)"100 return "Other"101bl_imp=imp.groupby(imp.index.map(block)).sum().sort_values()102bl_imp=bl_imp/bl_imp.sum()*100103fig,ax=plt.subplots(figsize=(6.4,3.6))104colors=["#8899AA"]*len(bl_imp)105for i,n in enumerate(bl_imp.index):106 if n=="Vision (façade)": colors[i]=C_GBM107 if n=="Text (LLM)": colors[i]=C_BL108ax.barh(bl_imp.index,bl_imp.values,color=colors,edgecolor="white")109for i,v in enumerate(bl_imp.values): ax.text(v+0.4,i,f"{v:.0f}\\%",va="center",fontsize=8)110ax.set_xlabel("Share of permutation importance (\\%)")111ax.set_title("Feature-block importance (gradient boosting, v5)",fontsize=11)112ax.set_xlim(0,max(bl_imp.values)*1.15)113fig.savefig("paper/figures/fig_importance.pdf");plt.close(fig)114print("-> fig_importance.pdf")115print("blocs:",{k:round(v,1) for k,v in bl_imp.sort_values(ascending=False).items()})116117# ===== Table ablation détaillée =====118def med(p): return np.median(np.abs(np.exp(p-y)-1))*100119rowsL=[]120labels=["Structural + amenity","\\quad + Rooms + listing","\\quad + Text (LLM)","\\quad + Vision (façade)"]121for i,(nm,cols) in enumerate(sets):122 ro,po=ols_r2(cols);rg,pg=gbm_r2(cols)123 bcell = f"{blend[i]:.3f}" if i in blend else "---"124 rowsL.append(f" {labels[i]} & {len(cols)} & {ro:.3f} & {rg:.3f} & {bcell} \\\\")125tab=r"""% Auto-generated v5126\begin{table}[htbp]127 \centering128 \begin{threeparttable}129 \caption{Feature Ablation: Out-of-Sample $R^2$ by Feature Layer}130 \label{tab:ablation}131 \begin{tabular}{@{}lcccc@{}}132 \toprule133 Feature set & \# feats & OLS & GBM & Blend \\134 \midrule135"""+"\n".join(rowsL)+r"""136 \bottomrule137 \end{tabular}138 \begin{tablenotes}[flushleft]139 \footnotesize140 \item \textit{Notes:} Each row adds a feature layer to the row above (nested).141 Entries are out-of-sample $R^2$ on $\ln(\text{Price})$ from 5-fold cross-validation142 with a common seed. ``\# feats'' counts numeric/indicator features before the city,143 type, and source fixed effects (added in every specification). The blend averages144 the robust neural network and gradient boosting and is reported where both were145 trained. GBM and the neural network handle missing living area natively; OLS uses146 within-type median imputation.147 \end{tablenotes}148 \end{threeparttable}149\end{table}150"""151open("paper/tables/ablation.tex","w").write(tab)152print("-> paper/tables/ablation.tex")153