Repository Public Research preview
1# Author: Simon-Pierre Boucher — contact@spboucher.ai2"""305_descriptive_analysis.py4==========================5Generate descriptive statistics tables and figures for the Airbnb-rent6linkage analysis.78Tables -> results/tables/ (.tex + .csv)9Figures -> figures/ (.pdf)1011Inputs:12 data/processed/merged_analysis.parquet (preferred)13 data/processed/merged_spatial.parquet (fallback)14 data/processed/airbnb_clean.parquet15 data/processed/rent_clean.parquet16"""1718import sys19from pathlib import Path2021sys.path.insert(0, str(Path(__file__).resolve().parent.parent))2223from src.plotting import use_descriptive_style2425use_descriptive_style()2627import matplotlib.pyplot as plt # noqa: E40228import matplotlib.ticker as mticker # noqa: E40229import numpy as np # noqa: E40230import pandas as pd # noqa: E4023132from src.config import ( # noqa: E40233 AIRBNB_CLEAN,34 RENT_CLEAN,35 MERGED_ANALYSIS,36 MERGED_SPATIAL,37 FIG_DIR,38 TABLE_DIR,39 require,40)4142PROCESSED_HINT = "Run scripts/04_merge_data.py first."434445# ===================================================================46# Helper: summary statistics47# ===================================================================4849def summary_stats(df: pd.DataFrame, cols: list[str], labels: dict | None = None) -> pd.DataFrame:50 """51 Compute N, mean, sd, min, p25, median, p75, max for *cols*.52 Columns absent from *df* are silently skipped (matches the published53 tables, which omit variables not present in the processed data).54 """55 rows = []56 for c in cols:57 if c not in df.columns:58 continue59 s = df[c].dropna().astype(float) # booleans -> shares60 label = (labels or {}).get(c, c)61 rows.append({62 "Variable": label,63 "N": int(len(s)),64 "Mean": s.mean(),65 "SD": s.std(),66 "Min": s.min(),67 "P25": s.quantile(0.25),68 "Median": s.median(),69 "P75": s.quantile(0.75),70 "Max": s.max(),71 })72 return pd.DataFrame(rows)737475def save_table(df: pd.DataFrame, name: str):76 """Save as .csv and .tex (a tabular fragment, wrapped by the paper)."""77 csv_path = TABLE_DIR / f"{name}.csv"78 tex_path = TABLE_DIR / f"{name}.tex"79 df.to_csv(csv_path, index=False)80 # LaTeX: format floats to 2 decimals except integers81 float_fmt = "%.2f"82 df.to_latex(tex_path, index=False, float_format=float_fmt)83 print(f" Saved {csv_path.name} + {tex_path.name}")848586# ===================================================================87# Tables88# ===================================================================8990def make_tables(airbnb: pd.DataFrame, rent: pd.DataFrame, merged: pd.DataFrame):91 print("\n--- Tables ---")9293 # 1. Summary statistics: Airbnb94 airbnb_vars = [95 "price_numeric", "rating", "bedrooms", "bathrooms",96 "guests_count", "amenities_count", "num_images", "num_reviews",97 "quality_score", "is_superhost", "is_entire_home",98 ]99 airbnb_labels = {100 "price_numeric": "Nightly price (CAD)",101 "is_superhost": "Superhost (share)",102 "is_entire_home": "Entire home (share)",103 "rating": "Rating",104 "bedrooms": "Bedrooms",105 "bathrooms": "Bathrooms",106 "guests_count": "Max guests",107 "amenities_count": "Amenities count",108 "num_images": "Number of images",109 "num_reviews": "Number of reviews",110 "quality_score": "Quality score",111 }112 tbl_airbnb = summary_stats(airbnb, airbnb_vars, airbnb_labels)113 save_table(tbl_airbnb, "summary_stats_airbnb")114115 # 2. Summary statistics: Rent116 # 2026-10-17: the original requested a non-existent column ``size_sqft``;117 # the interior size is stored as ``size_interior_sqft`` (3,699 listings).118 rent_vars = [119 "monthly_rent", "bedrooms", "bathrooms", "size_interior_sqft",120 ]121 rent_labels = {122 "monthly_rent": "Monthly rent (CAD)",123 "bedrooms": "Bedrooms",124 "bathrooms": "Bathrooms",125 "size_interior_sqft": "Interior size (sq ft)",126 }127 tbl_rent = summary_stats(rent, rent_vars, rent_labels)128 save_table(tbl_rent, "summary_stats_rent")129130 # 3. Correlation matrix131 corr_vars = []132 # Rent outcome133 for v in ["log_rent", "monthly_rent"]:134 if v in merged.columns:135 corr_vars.append(v)136 break # prefer log_rent137138 # Airbnb exposure measures (500m buffer)139 for v in ["airbnb_count_500m", "airbnb_density_500m",140 "mean_airbnb_price_500m", "share_entire_home_500m"]:141 if v in merged.columns:142 corr_vars.append(v)143144 # Property characteristics145 for v in ["bedrooms", "bathrooms", "size_sqft"]:146 if v in merged.columns:147 corr_vars.append(v)148149 if len(corr_vars) >= 3:150 corr_df = merged[corr_vars].dropna()151 corr_mat = corr_df.corr()152 corr_mat.round(3).to_csv(TABLE_DIR / "correlation_matrix.csv")153 # 2026-10-17: emit a publication-ready fragment (numbered columns,154 # readable row labels) instead of raw pandas output with unescaped155 # underscores, which could not be \input in the paper.156 corr_labels = {157 "log_rent": "Log monthly rent",158 "airbnb_count_500m": "Airbnb count (500 m)",159 "airbnb_density_500m": "Airbnb density (500 m)",160 "mean_airbnb_price_500m": "Mean Airbnb price (500 m)",161 "share_entire_home_500m": "Entire-home share (500 m)",162 "bedrooms": "Bedrooms",163 "bathrooms": "Bathrooms",164 }165 k = len(corr_vars)166 lines = [r"\begin{tabular}{l" + "r" * k + "}", r"\toprule",167 " & " + " & ".join(f"({i + 1})" for i in range(k)) + r" \\",168 r"\midrule"]169 for i, v in enumerate(corr_vars):170 row = [f"({i + 1}) {corr_labels.get(v, v)}"]171 for j, w in enumerate(corr_vars):172 row.append(f"{corr_mat.loc[v, w]:.2f}" if j <= i else "")173 lines.append(" & ".join(row) + r" \\")174 lines += [r"\bottomrule", r"\end{tabular}",175 r"\parbox{\textwidth}{\footnotesize Pearson correlations on the "176 + f"{len(corr_df):,}" + r" rental listings with at least one Airbnb "177 r"listing within 500\,m (the composition variables are undefined "178 r"otherwise). At a fixed radius the density is a rescaling of the "179 r"count, hence the correlation of 1.00.}"]180 (TABLE_DIR / "correlation_matrix.tex").write_text("\n".join(lines) + "\n",181 encoding="utf-8")182 print(" Saved correlation_matrix.csv + correlation_matrix.tex")183 else:184 print(" Skipping correlation matrix — not enough overlapping columns.")185186 # 4. Exposure summary by buffer radius (2026-10-17: the paper previously187 # described these statistics in prose only)188 rows = []189 for tag, lab in [("250m", "250 m"), ("500m", "500 m"), ("1000m", "1 km"), ("2000m", "2 km")]:190 c = merged[f"airbnb_count_{tag}"].dropna()191 rows.append({192 "Radius": lab, "N": int(len(c)), "Mean": c.mean(), "SD": c.std(),193 "Median": c.median(), "P90": c.quantile(0.90), "Max": c.max(),194 "Share zero": (c == 0).mean(),195 })196 expo = pd.DataFrame(rows)197 expo.to_csv(TABLE_DIR / "exposure_summary.csv", index=False)198 lines = [r"\begin{tabular}{lrrrrrrr}", r"\toprule",199 r"Radius & N & Mean & SD & Median & P90 & Max & Share zero \\", r"\midrule"]200 for _, r in expo.iterrows():201 lines.append(f"{r['Radius']} & {r['N']:,} & {r['Mean']:.2f} & {r['SD']:.2f} & "202 f"{r['Median']:.0f} & {r['P90']:.0f} & {r['Max']:.0f} & {r['Share zero']:.3f} \\\\")203 lines += [r"\bottomrule", r"\end{tabular}",204 r"\parbox{\textwidth}{\footnotesize Number of Airbnb listings within each "205 r"radius of a rental listing (Haversine distance), over all rental listings "206 r"with valid coordinates. ``Share zero'' is the fraction of rental listings "207 r"with no Airbnb listing within the radius.}"]208 (TABLE_DIR / "exposure_summary.tex").write_text("\n".join(lines) + "\n", encoding="utf-8")209 print(" Saved exposure_summary.csv + exposure_summary.tex")210211212# ===================================================================213# Figures214# ===================================================================215216def fig_dist_airbnb_price(airbnb: pd.DataFrame):217 """Histogram of Airbnb nightly prices."""218 prices = airbnb["price_numeric"].dropna()219 # Trim extreme outliers (above 99th percentile) for readability220 p99 = prices.quantile(0.99)221 prices_trim = prices[prices <= p99]222223 fig, ax = plt.subplots(figsize=(6, 4))224 ax.hist(prices_trim, bins=60, color="#4C72B0", edgecolor="white", linewidth=0.5)225 ax.set_xlabel("Nightly price (CAD)")226 ax.set_ylabel("Count")227 ax.set_title("Distribution of Airbnb Nightly Prices")228 ax.xaxis.set_major_formatter(mticker.StrMethodFormatter("${x:,.0f}"))229 fig.tight_layout()230 fig.savefig(FIG_DIR / "dist_airbnb_price.pdf")231 plt.close(fig)232 print(" dist_airbnb_price.pdf")233234235def fig_dist_airbnb_logprice(airbnb: pd.DataFrame):236 """Histogram of log nightly prices (added 2026-10-17): documents the237 bimodality of the scraped price field (second mode above $1,000)."""238 prices = airbnb["price_numeric"].dropna()239 fig, ax = plt.subplots(figsize=(6, 4))240 ax.hist(np.log(prices), bins=40, color="#4C72B0", edgecolor="white", linewidth=0.5)241 ax.axvline(np.log(1000), color="firebrick", linestyle="--", linewidth=1.2,242 label="$1,000 per night")243 ticks = [50, 100, 200, 500, 1000, 2000, 5000, 10000]244 ax.set_xticks(np.log(ticks))245 ax.set_xticklabels([f"${t:,}" for t in ticks], rotation=0, fontsize=8)246 ax.set_xlabel("Nightly price (CAD, log scale)")247 ax.set_ylabel("Count")248 ax.set_title("Distribution of Airbnb Nightly Prices (log scale)")249 ax.legend(fontsize=8, frameon=False)250 fig.tight_layout()251 fig.savefig(FIG_DIR / "dist_airbnb_logprice.pdf")252 plt.close(fig)253 print(" dist_airbnb_logprice.pdf")254255256def fig_dist_rent(rent: pd.DataFrame):257 """Histogram of monthly rents."""258 rents = rent["monthly_rent"].dropna()259 p99 = rents.quantile(0.99)260 rents_trim = rents[rents <= p99]261262 fig, ax = plt.subplots(figsize=(6, 4))263 ax.hist(rents_trim, bins=60, color="#DD8452", edgecolor="white", linewidth=0.5)264 ax.set_xlabel("Monthly rent (CAD)")265 ax.set_ylabel("Count")266 ax.set_title("Distribution of Monthly Rents")267 ax.xaxis.set_major_formatter(mticker.StrMethodFormatter("${x:,.0f}"))268 fig.tight_layout()269 fig.savefig(FIG_DIR / "dist_rent.pdf")270 plt.close(fig)271 print(" dist_rent.pdf")272273274def fig_airbnb_by_city(airbnb: pd.DataFrame):275 """Bar chart of Airbnb listing count by top 15 cities."""276 city_col = "city_clean" if "city_clean" in airbnb.columns else "city"277 counts = airbnb[city_col].value_counts().head(15)278279 fig, ax = plt.subplots(figsize=(8, 5))280 counts.sort_values().plot.barh(ax=ax, color="#4C72B0", edgecolor="white")281 ax.set_xlabel("Number of listings")282 ax.set_ylabel("")283 ax.set_title("Airbnb Listings by City (Top 15)")284 fig.tight_layout()285 fig.savefig(FIG_DIR / "airbnb_by_city.pdf")286 plt.close(fig)287 print(" airbnb_by_city.pdf")288289290def fig_rent_by_city(rent: pd.DataFrame):291 """Bar chart of mean rent by top 15 cities."""292 city_counts = rent["city"].value_counts()293 top_cities = city_counts.head(15).index294 sub = rent[rent["city"].isin(top_cities)]295 means = sub.groupby("city")["monthly_rent"].mean().sort_values()296297 fig, ax = plt.subplots(figsize=(8, 5))298 means.plot.barh(ax=ax, color="#DD8452", edgecolor="white")299 ax.set_xlabel("Mean monthly rent (CAD)")300 ax.set_ylabel("")301 ax.set_title("Mean Monthly Rent by City (Top 15)")302 ax.xaxis.set_major_formatter(mticker.StrMethodFormatter("${x:,.0f}"))303 fig.tight_layout()304 fig.savefig(FIG_DIR / "rent_by_city.pdf")305 plt.close(fig)306 print(" rent_by_city.pdf")307308309def fig_scatter_airbnb_rent(merged: pd.DataFrame):310 """Scatter plot of mean Airbnb count (500m) vs log rent311 at a neighbourhood level (postal code or small-area average)."""312 if "airbnb_count_500m" not in merged.columns:313 print(" Skipping scatter_airbnb_rent.pdf — no airbnb_count_500m column.")314 return315316 rent_col = "log_rent" if "log_rent" in merged.columns else "monthly_rent"317318 # Aggregate at postal-code level if available, else use city319 group_col = None320 for candidate in ["postal_code", "city"]:321 if candidate in merged.columns:322 group_col = candidate323 break324325 if group_col is None:326 print(" Skipping scatter_airbnb_rent.pdf — no grouping column found.")327 return328329 agg = merged.groupby(group_col).agg(330 mean_airbnb_count=("airbnb_count_500m", "mean"),331 mean_rent=(rent_col, "mean"),332 n=("airbnb_count_500m", "size"),333 ).reset_index()334335 # Keep neighbourhoods with at least 5 observations336 agg = agg[agg["n"] >= 5]337338 fig, ax = plt.subplots(figsize=(6, 5))339 ax.scatter(agg["mean_airbnb_count"], agg["mean_rent"],340 alpha=0.6, s=20, color="#4C72B0", edgecolors="none")341 ax.set_xlabel("Mean Airbnb count within 500 m")342 y_label = "Log monthly rent" if rent_col == "log_rent" else "Mean monthly rent (CAD)"343 ax.set_ylabel(y_label)344 ax.set_title("Airbnb Density vs. Rent")345 fig.tight_layout()346 fig.savefig(FIG_DIR / "scatter_airbnb_rent.pdf")347 plt.close(fig)348 print(" scatter_airbnb_rent.pdf")349350351def fig_map_airbnb(airbnb: pd.DataFrame):352 """Scatter map of Airbnb listings coloured by price."""353 lon_col = "lon" if "lon" in airbnb.columns else "long"354 valid = airbnb.dropna(subset=["lat", lon_col, "price_numeric"])355356 # Cap at 99th percentile for colour scale357 p99 = valid["price_numeric"].quantile(0.99)358359 fig, ax = plt.subplots(figsize=(8, 7))360 sc = ax.scatter(361 valid[lon_col], valid["lat"],362 c=valid["price_numeric"].clip(upper=p99),363 cmap="YlOrRd", s=4, alpha=0.6, edgecolors="none",364 )365 cbar = fig.colorbar(sc, ax=ax, shrink=0.7)366 cbar.set_label("Nightly price (CAD)")367 ax.set_xlabel("Longitude")368 ax.set_ylabel("Latitude")369 ax.set_title("Airbnb Listings — Quebec")370 fig.tight_layout()371 fig.savefig(FIG_DIR / "map_airbnb.pdf")372 plt.close(fig)373 print(" map_airbnb.pdf")374375376def fig_map_rent(rent: pd.DataFrame):377 """Scatter map of rental listings coloured by rent."""378 valid = rent.dropna(subset=["lat", "lon", "monthly_rent"])379 p99 = valid["monthly_rent"].quantile(0.99)380381 fig, ax = plt.subplots(figsize=(8, 7))382 sc = ax.scatter(383 valid["lon"], valid["lat"],384 c=valid["monthly_rent"].clip(upper=p99),385 cmap="YlGnBu", s=4, alpha=0.6, edgecolors="none",386 )387 cbar = fig.colorbar(sc, ax=ax, shrink=0.7)388 cbar.set_label("Monthly rent (CAD)")389 ax.set_xlabel("Longitude")390 ax.set_ylabel("Latitude")391 ax.set_title("Rental Listings — Quebec")392 fig.tight_layout()393 fig.savefig(FIG_DIR / "map_rent.pdf")394 plt.close(fig)395 print(" map_rent.pdf")396397398# ===================================================================399# Main400# ===================================================================401402def main():403 print("=" * 70)404 print("05 DESCRIPTIVE ANALYSIS")405 print("=" * 70)406407 # ------------------------------------------------------------------408 # Load data409 # ------------------------------------------------------------------410 print("\n[1] Loading data ...")411412 airbnb = pd.read_parquet(require(AIRBNB_CLEAN, PROCESSED_HINT))413 if "long" in airbnb.columns and "lon" not in airbnb.columns:414 airbnb = airbnb.rename(columns={"long": "lon"})415416 rent = pd.read_parquet(require(RENT_CLEAN, PROCESSED_HINT))417418 merged_path = MERGED_ANALYSIS419 if not merged_path.exists():420 merged_path = MERGED_SPATIAL421 merged = pd.read_parquet(require(merged_path, PROCESSED_HINT))422 if "long" in merged.columns and "lon" not in merged.columns:423 merged = merged.rename(columns={"long": "lon"})424425 print(f" Airbnb : {airbnb.shape}")426 print(f" Rent : {rent.shape}")427 print(f" Merged : {merged.shape} ({merged_path.name})")428429 # ------------------------------------------------------------------430 # Tables431 # ------------------------------------------------------------------432 make_tables(airbnb, rent, merged)433434 # ------------------------------------------------------------------435 # Figures436 # ------------------------------------------------------------------437 print("\n--- Figures ---")438 fig_dist_airbnb_price(airbnb)439 fig_dist_airbnb_logprice(airbnb)440 fig_dist_rent(rent)441 fig_airbnb_by_city(airbnb)442 fig_rent_by_city(rent)443 fig_scatter_airbnb_rent(merged)444 fig_map_airbnb(airbnb)445 fig_map_rent(rent)446447 print("\nDone.\n")448449450if __name__ == "__main__":451 main()452