{ "cells": [ { "cell_type": "markdown", "id": "d87adfde", "metadata": {}, "source": [ "# Feature Engineering for Credit & Risk\n", "\n", "Companion notebook to [Data Science Module 2](https://rodoslay.com/data-science/02-feature-engineering-credit-risk) at **rodoslay.com**.\n", "\n", "What this notebook does:\n", "\n", "1. Generates a **synthetic loan-month panel** (5,000 loans × 36 months) with engineered-in default behavior\n", "2. Sets up a leak-proof **observation design** (snapshot date, performance window)\n", "3. Builds the feature families from the post: **ratios with hygiene, lags, rolling windows, trend slopes**\n", "4. Implements **WOE encoding with smoothing**, the right way (out-of-fold) and the wrong way (full-data)\n", "5. **Measures the leakage**: how much a pure-noise categorical appears to help when its WOE is fitted on data that includes your test set\n", "\n", "Dependencies: `numpy`, `pandas`, `matplotlib`, `scikit-learn`." ] }, { "cell_type": "code", "execution_count": 1, "id": "581166b4", "metadata": { "execution": { "iopub.execute_input": "2026-07-14T00:49:47.141376Z", "iopub.status.busy": "2026-07-14T00:49:47.141169Z", "iopub.status.idle": "2026-07-14T00:49:48.349248Z", "shell.execute_reply": "2026-07-14T00:49:48.348312Z" } }, "outputs": [], "source": [ "import numpy as np\n", "import pandas as pd\n", "import matplotlib.pyplot as plt\n", "from sklearn.linear_model import LogisticRegression\n", "from sklearn.metrics import roc_auc_score\n", "\n", "rng = np.random.default_rng(42)\n", "plt.rcParams[\"figure.figsize\"] = (9, 4.5)\n", "plt.rcParams[\"axes.grid\"] = True\n", "plt.rcParams[\"grid.alpha\"] = 0.3" ] }, { "cell_type": "markdown", "id": "574a73cd", "metadata": {}, "source": [ "## 1. A synthetic loan-month panel\n", "\n", "One row per loan per month: balance, credit limit, payment made, days past due. Each loan has a\n", "**latent riskiness** `z` that drives its behavior — riskier loans drift toward higher utilization,\n", "pay less of what's due, and go delinquent more often. Default is a monthly hazard that depends on the\n", "latent risk *and* the behavior, so behavioral features genuinely carry signal. A pure-noise\n", "`region` code (50 levels) is included deliberately — it's the bait for the leakage experiment." ] }, { "cell_type": "code", "execution_count": 2, "id": "6022908c", "metadata": { "execution": { "iopub.execute_input": "2026-07-14T00:49:48.352488Z", "iopub.status.busy": "2026-07-14T00:49:48.352120Z", "iopub.status.idle": "2026-07-14T00:49:48.421070Z", "shell.execute_reply": "2026-07-14T00:49:48.420136Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "panel rows: 164,767\n", "loans defaulting at some point in 36m: 16.6%\n" ] } ], "source": [ "N_LOANS, N_MONTHS = 5000, 36\n", "\n", "z = rng.normal(size=N_LOANS) # latent riskiness\n", "limit = np.round(rng.lognormal(mean=8.6, sigma=0.7, size=N_LOANS), -2) + 100\n", "region = rng.integers(0, 50, size=N_LOANS) # pure noise, 50 levels\n", "\n", "util = np.zeros((N_LOANS, N_MONTHS + 1))\n", "util[:, 0] = np.clip(0.35 + 0.12 * z + rng.normal(0, 0.15, N_LOANS), 0, 1.4)\n", "dpd = np.zeros((N_LOANS, N_MONTHS + 1), dtype=int)\n", "pay_ratio = np.zeros((N_LOANS, N_MONTHS + 1))\n", "alive = np.ones(N_LOANS, dtype=bool)\n", "default_month = np.full(N_LOANS, -1)\n", "\n", "def sigmoid(x):\n", " return 1.0 / (1.0 + np.exp(-x))\n", "\n", "for t in range(1, N_MONTHS + 1):\n", " drift = 0.010 * z + rng.normal(0, 0.045, N_LOANS) # risky loans drift up\n", " util[:, t] = np.clip(util[:, t - 1] + drift, 0, 1.5)\n", " pay_ratio[:, t] = np.clip(1.6 - 1.1 * util[:, t] - 0.25 * z\n", " + rng.normal(0, 0.3, N_LOANS), 0, 4)\n", " p_dpd = sigmoid(-3.4 + 1.0 * z + 1.6 * (util[:, t] - 0.5))\n", " new_dpd = rng.random(N_LOANS) < p_dpd\n", " dpd[:, t] = np.where(new_dpd, dpd[:, t - 1] + 30, 0) # cure or roll 30 days deeper\n", "\n", " hazard = sigmoid(-6.0 + 0.9 * z + 1.9 * (util[:, t] - 0.5) + 0.028 * dpd[:, t])\n", " defaults_now = alive & (rng.random(N_LOANS) < hazard)\n", " default_month[defaults_now] = t\n", " alive &= ~defaults_now\n", "\n", "panel = pd.DataFrame({\n", " \"loan_id\": np.repeat(np.arange(N_LOANS), N_MONTHS),\n", " \"month\": np.tile(np.arange(1, N_MONTHS + 1), N_LOANS),\n", " \"limit\": np.repeat(limit, N_MONTHS),\n", " \"region\": np.repeat(region, N_MONTHS),\n", " \"utilization_raw\": util[:, 1:].ravel(),\n", " \"pay_ratio\": pay_ratio[:, 1:].ravel(),\n", " \"dpd\": dpd[:, 1:].ravel(),\n", "})\n", "panel[\"balance\"] = panel[\"utilization_raw\"] * panel[\"limit\"]\n", "# rows after a loan's default month don't exist in a real servicing system\n", "dm = np.repeat(default_month, N_MONTHS)\n", "panel = panel[(dm < 0) | (panel[\"month\"] <= dm)].reset_index(drop=True)\n", "\n", "print(f\"panel rows: {len(panel):,}\")\n", "print(f\"loans defaulting at some point in 36m: {(default_month > 0).mean():.1%}\")" ] }, { "cell_type": "markdown", "id": "21ce541a", "metadata": {}, "source": [ "## 2. Observation design — before any feature exists\n", "\n", "Three decisions, fixed here once so every later cell obeys them:\n", "\n", "- **Snapshot month: t = 24.** We stand at the end of month 24 and score every loan still alive.\n", "- **Feature window: months 1–23.** Rolling windows end at **t − 1 = 23**, because month 24 is\n", " \"partially observed\" at the moment of scoring. This is the boundary rule from the post.\n", "- **Performance window: months 25–36.** Label = 1 if the loan defaults within the next 12 months." ] }, { "cell_type": "code", "execution_count": 3, "id": "58e3df54", "metadata": { "execution": { "iopub.execute_input": "2026-07-14T00:49:48.423800Z", "iopub.status.busy": "2026-07-14T00:49:48.423591Z", "iopub.status.idle": "2026-07-14T00:49:48.429562Z", "shell.execute_reply": "2026-07-14T00:49:48.428680Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "loans alive at snapshot: 4,411 12m default rate: 5.51%\n" ] } ], "source": [ "SNAP = 24\n", "in_sample = (default_month < 0) | (default_month > SNAP) # alive at snapshot\n", "label = (default_month > SNAP) & (default_month <= SNAP + 12) # default in months 25..36\n", "\n", "ids = np.arange(N_LOANS)[in_sample]\n", "y = label[in_sample].astype(int)\n", "print(f\"loans alive at snapshot: {len(ids):,} 12m default rate: {y.mean():.2%}\")" ] }, { "cell_type": "markdown", "id": "07d34eed", "metadata": {}, "source": [ "## 3. Ratios — with hygiene\n", "\n", "Utilization is balance/limit. The hygiene rules from the post: cap the **ratio** (not the inputs),\n", "and flag degenerate denominators instead of letting them produce absurd values." ] }, { "cell_type": "code", "execution_count": 4, "id": "95dae30b", "metadata": { "execution": { "iopub.execute_input": "2026-07-14T00:49:48.431919Z", "iopub.status.busy": "2026-07-14T00:49:48.431739Z", "iopub.status.idle": "2026-07-14T00:49:48.448044Z", "shell.execute_reply": "2026-07-14T00:49:48.446625Z" } }, "outputs": [], "source": [ "hist = panel[(panel[\"month\"] <= SNAP - 1) & panel[\"loan_id\"].isin(ids)]\n", "\n", "def ratio_capped(num, den, cap):\n", " degenerate = den < 1e-9\n", " r = np.where(degenerate, 0.0, num / np.maximum(den, 1e-9))\n", " return np.minimum(r, cap), degenerate.astype(int)\n", "\n", "snap_row = hist[hist[\"month\"] == SNAP - 1].set_index(\"loan_id\").reindex(ids)\n", "util_now, util_degen = ratio_capped(snap_row[\"balance\"].values,\n", " snap_row[\"limit\"].values, cap=2.0)\n", "features = pd.DataFrame(index=ids)\n", "features[\"util_now\"] = util_now\n", "features[\"util_degenerate\"] = util_degen\n", "features[\"pay_ratio_now\"] = snap_row[\"pay_ratio\"].values\n", "features[\"dpd_now\"] = snap_row[\"dpd\"].values" ] }, { "cell_type": "markdown", "id": "8171840b", "metadata": {}, "source": [ "## 4. Lags, rolling windows, trend\n", "\n", "- **Lags**: the same quantity k months earlier — lets a pointwise model see change.\n", "- **Rolling aggregates**: mean for level, **max for \"how bad did it ever get\"** (usually the stronger\n", " risk feature), counts of events. All windows end at month 23.\n", "- **Trend**: the slope of utilization over the last 6 observable months. Two borrowers at 60%\n", " utilization are different animals if one came from 30% and the other from 90%." ] }, { "cell_type": "code", "execution_count": 5, "id": "28e2dbe5", "metadata": { "execution": { "iopub.execute_input": "2026-07-14T00:49:48.450147Z", "iopub.status.busy": "2026-07-14T00:49:48.449953Z", "iopub.status.idle": "2026-07-14T00:49:48.510216Z", "shell.execute_reply": "2026-07-14T00:49:48.509220Z" } }, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
countmeanstdmin25%50%75%max
util_now4411.00.3530.3020.00.0810.2960.5601.500
util_degenerate4411.00.0000.0000.00.0000.0000.0000.000
pay_ratio_now4411.01.2650.5900.00.8361.3021.6993.074
dpd_now4411.01.3816.8760.00.0000.0000.00060.000
util_lag34411.00.3460.2880.00.0860.2980.5471.455
util_lag64411.00.3420.2730.00.1010.3050.5361.384
util_lag124411.00.3330.2440.00.1270.3140.5051.271
util_mean_3m4411.00.3500.2960.00.0820.2920.5551.459
dpd_max_3m4411.03.66610.8900.00.0000.0000.00090.000
payr_min_3m4411.01.0140.5380.00.6141.0401.4212.795
util_mean_6m4411.00.3470.2880.00.0820.2990.5511.396
dpd_max_6m4411.06.38614.1360.00.0000.0000.000120.000
\n", "
" ], "text/plain": [ " count mean std min 25% 50% 75% max\n", "util_now 4411.0 0.353 0.302 0.0 0.081 0.296 0.560 1.500\n", "util_degenerate 4411.0 0.000 0.000 0.0 0.000 0.000 0.000 0.000\n", "pay_ratio_now 4411.0 1.265 0.590 0.0 0.836 1.302 1.699 3.074\n", "dpd_now 4411.0 1.381 6.876 0.0 0.000 0.000 0.000 60.000\n", "util_lag3 4411.0 0.346 0.288 0.0 0.086 0.298 0.547 1.455\n", "util_lag6 4411.0 0.342 0.273 0.0 0.101 0.305 0.536 1.384\n", "util_lag12 4411.0 0.333 0.244 0.0 0.127 0.314 0.505 1.271\n", "util_mean_3m 4411.0 0.350 0.296 0.0 0.082 0.292 0.555 1.459\n", "dpd_max_3m 4411.0 3.666 10.890 0.0 0.000 0.000 0.000 90.000\n", "payr_min_3m 4411.0 1.014 0.538 0.0 0.614 1.040 1.421 2.795\n", "util_mean_6m 4411.0 0.347 0.288 0.0 0.082 0.299 0.551 1.396\n", "dpd_max_6m 4411.0 6.386 14.136 0.0 0.000 0.000 0.000 120.000" ] }, "execution_count": 5, "metadata": {}, "output_type": "execute_result" } ], "source": [ "g = hist.groupby(\"loan_id\")\n", "\n", "def month_value(col, m):\n", " s = hist[hist[\"month\"] == m].set_index(\"loan_id\")[col]\n", " return s.reindex(ids).values\n", "\n", "for k in (3, 6, 12):\n", " features[f\"util_lag{k}\"] = month_value(\"utilization_raw\", SNAP - 1 - k)\n", "\n", "for w in (3, 6, 12):\n", " win = hist[hist[\"month\"] >= SNAP - w] # months (24-w) .. 23\n", " gw = win.groupby(\"loan_id\")\n", " features[f\"util_mean_{w}m\"] = gw[\"utilization_raw\"].mean().reindex(ids).values\n", " features[f\"dpd_max_{w}m\"] = gw[\"dpd\"].max().reindex(ids).values\n", " features[f\"payr_min_{w}m\"] = gw[\"pay_ratio\"].min().reindex(ids).values\n", "\n", "# 6-month utilization trend: least-squares slope over months 18..23\n", "win6 = hist[hist[\"month\"] >= SNAP - 6]\n", "x = win6[\"month\"] - (SNAP - 3.5) # centered time index\n", "slopes = (win6.assign(xy=x * win6[\"utilization_raw\"], xx=x * x)\n", " .groupby(\"loan_id\")[[\"xy\", \"xx\"]].sum())\n", "features[\"util_trend_6m\"] = (slopes[\"xy\"] / slopes[\"xx\"]).reindex(ids).values\n", "\n", "features[\"region\"] = region[ids]\n", "features = features.fillna(0.0)\n", "features.describe().T.round(3).head(12)" ] }, { "cell_type": "markdown", "id": "011bf129", "metadata": {}, "source": [ "## 5. WOE encoding — with smoothing\n", "\n", "$$\\mathrm{WOE}_i = \\ln\\frac{\\%\\ \\text{goods in bin } i}{\\%\\ \\text{bads in bin } i}$$\n", "\n", "Two implementation details that separate production code from tutorial code:\n", "\n", "- **Smoothing**: a bin with 3 loans and 0 bads should shrink toward zero evidence, not claim\n", " WOE = +∞. We add a pseudo-count proportional to the portfolio's good/bad mix.\n", "- **The mapping is *fitted*** on one set of loans and *applied* to another — it's a supervised\n", " transformation, and that's exactly why it can leak (next section)." ] }, { "cell_type": "code", "execution_count": 6, "id": "97b1ca5a", "metadata": { "execution": { "iopub.execute_input": "2026-07-14T00:49:48.512582Z", "iopub.status.busy": "2026-07-14T00:49:48.512398Z", "iopub.status.idle": "2026-07-14T00:49:48.617857Z", "shell.execute_reply": "2026-07-14T00:49:48.616458Z" } }, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAwwAAAGcCAYAAAB0hicRAAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjksIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvJkbTWQAAAAlwSFlzAAAPYQAAD2EBqD+naQAASRFJREFUeJzt3Xl8EwX+//H3lEJKC6UcLW25L6EcgpxSEAoiLaIIKpd+LeWrgC58V771WNlVuVz7Q0RARZFVAQVclEVAXUEQCnIqICIKCMgtpZWjhQItNPP7g22+Ce1AC0ljwuv5eOTx6ExmMp+ZfJLmnTlimKZpCgAAAAAKEeDtAgAAAAD8cREYAAAAAFgiMAAAAACwRGAAAAAAYInAAAAAAMASgQEAAACAJQIDAAAAAEsEBgAAAACWCAwAAAAALBEYAD8wa9YsGYahzZs3e7uUa0pKSlLt2rVdxhmGoTFjxjiG89fnwIEDJVqbt5YrXd4GI0aMuOZ0nqrx8OHDCgoK0rp169z6uPBttWvXVlJSkmM4NTVVhmEoNTX1mvPGxcUpLi7OY7W524cffqhGjRqpdOnSCgsLu+b0V26bklSU99Gieu6559SuXTv3FAa/RWAAnHz88ccyDEOffvppgfuaN28uwzC0atWqAvfVrFlTsbGxLuMuXryo119/XW3atFH58uVVrlw5tWnTRq+//rouXrxY4DFq164twzAKvSUkJLhvJUvAb7/9pjFjxmjbtm3eLqVQL7/8shYtWuTtMv5Qxo0bp3bt2qlDhw7eLsUt5s2bpylTpni7DFyHc+fOacyYMUUKJe6ya9cuJSUlqV69evrHP/6hGTNmlNiyvW3kyJH64YcftGTJEm+Xgj+wQG8XAPyRdOzYUZK0du1a9enTxzE+KytLO3bsUGBgoNatW6cuXbo47jt8+LAOHz6sAQMGOMZlZ2erZ8+eWr16te655x4lJSUpICBAS5cu1ZNPPqmFCxfqiy++UEhIiMvyW7RooaeeeqpAXdHR0e5eVY/67bffNHbsWNWuXVstWrRwue8f//iH7Hb7Ved/5JFHNGDAANlsNo/U9/LLL+vBBx9U7969S3S57uCJGjMyMjR79mzNnj3bbY/pbfPmzdOOHTs0cuRIb5fiVzp16qTz58+rTJkyHlvGuXPnNHbsWEkqsT0Uqampstvtmjp1qurXr1+keXbv3q2AAN//3jUyMlL33XefXn31VfXq1cvb5eAPisAAOImOjladOnW0du1al/EbNmyQaZrq27dvgfvyh/PDhiQlJydr9erVeuONN1wOM3niiSc0bdo0jRgxQk8//bTefvttl8eqVq2a/uu//svdq/WHUrp06WtOU6pUKZUqVaoEqvljLLc4PFHjnDlzFBgYqHvvvdetjwv/ExAQoKCgIG+XcV2ys7MLfEmTLz09XZKueSiSaZq6cOGCypYt+4f+YqG4+vXrp759++rXX39V3bp1vV0O/oB8PxoDbtaxY0d9//33On/+vGPcunXr1KRJE/Xo0UMbN250+YZ83bp1MgzDcSjHkSNH9N5776lr166FHpM+fPhwdenSRe+++66OHDni1trPnTunYcOGqXLlygoNDVViYqJOnTrluH/QoEGqUqVKoYdEde/eXQ0bNrzq41sds+t8rHJqaqratGkjSRo8eLDjsKpZs2ZJKvzY2ytdeZz+mDFjLA/Xcq7n1VdfVWxsrCpXrqyyZcuqVatWWrBggctjG4ah7OxszZ49u8BjWJ0f8NZbb6lJkyay2WyKjo7W8OHDdfr06QLboGnTpvr555/VpUsXBQcHq1q1anrllVeuuq5Xmjt3rho2bKigoCC1atVKa9asueq2kS4/L/fcc4/Wrl2rtm3bKigoSHXr1tUHH3xQpGUuWrRI7dq1U7ly5Qpdp+3bt6tz584KDg5W/fr1Hdt09erVateuncqWLauGDRtqxYoVBR77+++/V48ePRQaGqpy5crpzjvv1MaNGwtdp3Xr1ik5OVnh4eEKCQlRnz59lJGRUeAxr/V8xMXF6YsvvtDBgwcdz7Fzz6Wnp+vRRx9V1apVFRQUpObNmxfYu3LgwAEZhqFXX31VM2bMUL169WSz2dSmTRt99913BWratWuXHnzwQVWqVElBQUFq3bp1iR/i0bRpU5e9n/nsdruqVaumBx980DGuKK+Vwlidw5C/jcqWLau2bdvqm2++KTBvbm6uXnzxRbVq1UoVKlRQSEiI7rjjDpfDPA8cOKDw8HBJ0tixYx3Pn/Ox+StXrtQdd9yhkJAQhYWF6b777tPOnTtdlpX/nvHzzz/roYceUsWKFV2+1HFWu3ZtjR49WpIUHh7usrz819ayZcvUunVrlS1bVu+8847jvivfD0+fPq2RI0eqRo0astlsql+/viZMmODyP6O4vbVo0SI1bdpUQUFBatq0aaGHzF5p1apVlofXzps3T4ZhaMOGDY5x3bp1kyQtXrz4mo+NmxOBAbhCx44ddfHiRW3atMkxbt26dYqNjVVsbKwyMzO1Y8cOl/saNWqkypUrS5K+/PJL5eXlKTEx0XIZiYmJunTpkpYuXeoy/uLFi/r9998L3JzDy9WMGDFCO3fu1JgxY5SYmKi5c+eqd+/eMk1T0uXDWU6cOKFly5a5zJeWlqaVK1e6Ze9GTEyMxo0bJ0kaOnSoPvzwQ3344Yfq1KnTdT/m/fff73ic/Fv+oSYRERGO6aZOnarbbrtN48aN08svv6zAwED17dtXX3zxhWOaDz/8UDabTXfccYfjsYYNG2a57DFjxmj48OGKjo7WpEmT9MADD+idd95R9+7dCwSvU6dOKSEhQc2bN9ekSZPUqFEj/eUvf9GXX35ZpPVcvXq1Ro4cqf/6r//SuHHjdOLECSUkJLj0m5W9e/fqwQcf1F133aVJkyapYsWKSkpK0k8//XTV+S5evKjvvvtOLVu2LPT+U6dO6Z577lG7du30yiuvyGazacCAAZo/f74GDBigu+++W//v//0/ZWdn68EHH9SZM2cc8/7000+644479MMPP+jZZ5/VCy+8oP379ysuLs7l9ZXvf/7nf/TDDz9o9OjReuKJJ/TZZ58VCN1FeT7+9re/qUWLFqpSpYrjOc4/n+H8+fOKi4vThx9+qIcfflgTJ05UhQoVlJSUpKlTpxaoad68eZo4caKGDRuml156SQcOHND999/v8tz/9NNPuv3227Vz504999xzmjRpkkJCQtS7d+8ifbg7depUoa/7K2/nzp276uP0799fa9asUVpamsv4tWvX6rfffnM5bLIor5Wieu+99zRs2DBFRkbqlVdeUYcOHdSrVy8dPnzYZbqsrCy9++67iouL04QJEzRmzBhlZGQoPj7ecb5TeHi4Y89rnz59HM/f/fffL0lasWKF4uPjlZ6erjFjxig5OVnr169Xhw4dCr0QQN++fXXu3Dm9/PLLGjJkSKH1T5kyxXEI6ttvv+2yPOnyoUcDBw7UXXfdpalTpxY4zDLfuXPn1LlzZ82ZM0eJiYl6/fXX1aFDB40aNUrJyckFpi9Kb3311Vd64IEHZBiGUlJS1Lt3bw0ePPiaF7iIi4tTjRo1NHfu3AL3zZ07V/Xq1VP79u0d4ypUqKB69epx0QNYMwG4+Omnn0xJ5vjx403TNM2LFy+aISEh5uzZs03TNM2qVaua06ZNM03TNLOyssxSpUqZQ4YMccw/cuRIU5L5/fffWy5j69atpiQzOTnZMa5WrVqmpEJvKSkpV6155syZpiSzVatWZm5urmP8K6+8YkoyFy9ebJqmaebl5ZnVq1c3+/fv7zL/a6+9ZhqGYf76669XXU6tWrXMQYMGFRjfuXNns3Pnzo7h7777zpRkzpw5s8C0gwYNMmvVquUyTpI5evToAuuzf//+QuvIyMgwa9asaTZr1sw8e/asY/y5c+dcpsvNzTWbNm1qdu3a1WV8SEhIoetx5XLT09PNMmXKmN27dzfz8vIc07355pumJPP999932QaSzA8++MAxLicnx4yMjDQfeOCBQtfDWf5zvXnzZse4gwcPmkFBQWafPn0sazTN/+udNWvWOMalp6ebNpvNfOqpp6663L1795qSzDfeeKPAffnrNG/ePMe4Xbt2mZLMgIAAc+PGjY7xy5YtK/Cc9+7d2yxTpoy5b98+x7jffvvNLF++vNmpU6cC69StWzfTbrc7xv/v//6vWapUKfP06dOOdSrq89GzZ88CfWaapjllyhRTkjlnzhzHuNzcXLN9+/ZmuXLlzKysLNM0TXP//v2mJLNy5crmyZMnHdMuXrzYlGR+9tlnjnF33nmn2axZM/PChQuOcXa73YyNjTUbNGhQoIYrXe2173xzfo0UZvfu3YU+l3/605/McuXKubw+ivpaufI1v2rVKlOSuWrVKsd8ERERZosWLcycnBzHdDNmzDAlubwvXLp0yWUa0zTNU6dOmVWrVjX/+7//2zEuIyPDcn1btGhhRkREmCdOnHCM++GHH8yAgAAzMTHRMW706NGmJHPgwIEFHqMw+dNnZGQUWH9J5tKlSwvMc+W2GT9+vBkSEmL+8ssvLtM999xzZqlSpcxDhw6Zplm83mrRooUZFRXleA2Ypml+9dVXpqRrvo+OGjXKtNlsLvOmp6ebgYGBhW7b7t27mzExMQU3DmCaJnsYgCvExMSocuXKjnMTfvjhB2VnZzuughQbG+v4FmbDhg3Ky8tz2dWd/w1r+fLlLZeRf19WVpbL+Hbt2mn58uUFbgMHDixS7UOHDnU5R+CJJ55QYGCg/v3vf0u6fPzxww8/rCVLlrh8Ezx37lzFxsaqTp06RVqON+Xl5WngwIE6c+aMPv30U5djksuWLev4+9SpU8rMzNQdd9yhrVu3XteyVqxYodzcXI0cOdLl5MYhQ4YoNDS0wLex5cqVc9lLU6ZMGbVt21a//vprkZbXvn17tWrVyjFcs2ZN3XfffVq2bJny8vKuOm/jxo11xx13OIbDw8PVsGHDay77xIkTkqSKFSsWen+5cuVcvplu2LChwsLCFBMT43Ipxvy/85eXl5enr776Sr1793Y5JjoqKkoPPfSQ1q5dW6D/hw4dKsMwHMN33HGH8vLydPDgQUnFfz4K8+9//1uRkZEur6nSpUvrz3/+s86ePavVq1e7TN+/f3+XbZO/jfPX8+TJk1q5cqX69eunM2fOOPYGnDhxQvHx8dqzZ4+OHj161Zrmzp1b6Ov+ytvV9lpK0i233KIWLVpo/vz5jnF5eXlasGCB7r33XpfXh7teK5s3b1Z6eroef/xxlxOhk5KSVKFCBZdpS5Uq5ZjGbrfr5MmTunTpklq3bl2k5R47dkzbtm1TUlKSKlWq5Bh/66236q677nK8zzl7/PHHi7U+halTp47i4+OvOd0nn3yiO+64QxUrVnTZM9StWzfl5eUVOLzwWr2Vv76DBg1y2ZZ33XWXGjdufM16EhMTlZOT43Ko2fz583Xp0qVC9ybn1w0UhpOegSsYhqHY2FitWbNGdrtd69atU0REhOPKGbGxsXrzzTclyREcnANDfhhw/kB+JatQUaVKFcexpNejQYMGLsPlypVTVFSUy676xMRETZgwQZ9++qkSExO1e/dubdmyRdOnT7/u5Zak559/XitXrtQXX3yhevXqudz3+eef66WXXtK2bduUk5PjGO/8IbQ48j+oXnluR5kyZVS3bl3H/fmqV69eYFkVK1bU9u3bi7S8K58/6fKHwHPnzikjI0ORkZGW89asWbPAuIoVK7qcw3I15n8OW7tSYetUoUIF1ahRo8A4SY7lZWRk6Ny5c4WeFxMTEyO73a7Dhw+rSZMmluuQ/2Eq/zGL+3wU5uDBg2rQoEGBq9vExMS4LKOoNe3du1emaeqFF17QCy+8UOgy09PTVa1aNcua3Hkp2/79++uvf/2rjh49qmrVqik1NVXp6enq37+/y3Tueq3kb68re7d06dKFnjw7e/ZsTZo0Sbt27XI59KYoX1ZYPf/S5edv2bJlBU5sdn7c8+fPKzMz02W+q72milObJO3Zs0fbt293nINxpfwTq/MVtd8Le19o2LDhNUNWo0aN1KZNG82dO1ePPvqopMvh9Pbbby/0SlCmaV73eyX8H4EBKETHjh312Wef6ccff3Scv5AvNjZWzzzzjI4ePaq1a9cqOjra5R9j/geP7du3Wx7rmv8BsijfErlb48aN1apVK8dxtnPmzFGZMmXUr1+/a85r9c8kLy+vRK4utGjRIk2YMEHjx48v8NsU33zzjXr16qVOnTrprbfeUlRUlEqXLq2ZM2dq3rx5Hq9NkuU2sPow/kdYdv65N1bBwupxPbGu3tx+Vq5VU/7JrE8//bTlt9DXukxnRkbGNfcgSZe/ALjyxPQr9e/fX6NGjdInn3yikSNH6uOPP1aFChVcXi/eeq3MmTNHSUlJ6t27t5555hlFRESoVKlSSklJ0b59+zyyTOc9KfPnz9fgwYNd7i9Kbzk/xtXY7XbdddddevbZZwu9/5ZbbnEZLol+T0xM1JNPPqkjR44oJydHGzdudHzhdaVTp06pSpUqbls2/AuBASiE8+8xrFu3zuVa7q1atZLNZlNqaqo2bdqku+++22XeHj16qFSpUvrwww8tDyH44IMPFBgY6PYfZNuzZ4/LVVLOnj2rY8eOFagxMTFRycnJOnbsmObNm6eePXtaHpLirGLFigWuDiRd/ibMOTR54luqX375RYMGDVLv3r3117/+tcD9//rXvxQUFKRly5a5XO5w5syZBaYtan21atWSdPmkR+f1y83N1f79+29ob1Bh9uzZU2DcL7/8ouDgYMtvLW9UzZo1VbZsWe3fv9+tjxseHq7g4GDt3r27wH27du1SQEBAgb0U11Kc58PqOa5Vq5a2b98uu93uspdh165dLssoqvw6Spcufd390KZNmyLtHRk9evQ1f8m3Tp06atu2rebPn68RI0Zo4cKF6t27t8trojivlWvJ31579uxR165dHeMvXryo/fv3q3nz5o5xCxYsUN26dbVw4UKX5yf/CkX5rvbcSbLsqSpVqlheNlWS4uPjtXz58iKs1fWpV6+ezp4967b3Bedte6XCtkFhBgwYoOTkZH300Uc6f/68SpcuXWBvU74rny/AGecwAIVo3bq1goKCNHfuXB09etRlD4PNZlPLli01bdo0ZWdnF7hUX40aNTR48GCtWLGiwO8sSNL06dO1cuVKPfroo6pevbpb654xY4bLbv63335bly5dUo8ePVymGzhwoAzD0JNPPqlff/21yFdHqlevnjZu3Kjc3FzHuM8//7zA1VDy/2kXFi6ux9mzZ9WnTx9Vq1bNcTnUK5UqVUqGYbh8U3vgwIFCf9E5JCSkSLV169ZNZcqU0euvv+7yrd97772nzMxM9ezZ87rWx8qGDRtcDjM4fPiwFi9erO7du3tsD07p0qXVunXra151pbhKlSql7t27a/HixS6HxB0/flzz5s1Tx44dFRoaWqzHLM7zERISUuDwE0m6++67lZaW5nKc/6VLl/TGG2+oXLly6ty5c7FqioiIUFxcnN555x0dO3aswP2FXRb2Su46hyFf//79tXHjRr3//vv6/fffC3xALM5r5Vpat26t8PBwTZ8+3eV9YdasWQVeY/k97Pzcbdq0yeXynpIUHBwsqeD7R1RUlFq0aKHZs2e73Ldjxw599dVXBb4YuVJUVJS6devmcnOnfv36acOGDQWuQiddXpdLly4V6/Gc19e5l5cvX66ff/65SI9RpUoV9ejRQ3PmzNHcuXOVkJBQ6F6EzMxM7du3z+V/HeCMPQxAIcqUKaM2bdrom2++kc1mczkRVbp8WNKkSZMkqdBre0+ePFm7du3Sn/70Jy1dutSxJ2HZsmVavHixOnfu7Jjf2dGjRzVnzpwC48uVK1fgV4kLk5ubqzvvvFP9+vXT7t279dZbb6ljx44Ffr0zPDxcCQkJ+uSTTxQWFlbkD76PPfaYFixYoISEBPXr10/79u3TnDlzCpxLUK9ePYWFhWn69OkqX768QkJC1K5du+s+qXrs2LH6+eef9fzzzxe4Tnj+5QF79uyp1157TQkJCXrooYeUnp6uadOmqX79+gXOIWjVqpVWrFih1157zfFjfc4n8OYLDw/XqFGjNHbsWCUkJKhXr16O7dqmTRu3/8he06ZNFR8frz//+c+y2Wx66623HOvvSffdd5/+9re/KSsrq9gf4q/mpZde0vLly9WxY0f96U9/UmBgoN555x3l5OQU+/cppOI9H61atdL8+fOVnJysNm3aqFy5crr33ns1dOhQvfPOO0pKStKWLVtUu3ZtLViwQOvWrdOUKVOuerECK9OmTVPHjh3VrFkzDRkyRHXr1tXx48e1YcMGHTlyRD/88MNV53fnOQzS5Q+uTz/9tJ5++mlVqlSpwAfj4rxWrqV06dJ66aWXNGzYMHXt2lX9+/fX/v37NXPmzALnMNxzzz1auHCh+vTpo549e2r//v2aPn26GjdurLNnzzqmK1u2rBo3bqz58+frlltuUaVKldS0aVM1bdpUEydOVI8ePdS+fXs9+uijOn/+vN544w1VqFDhmntfPO2ZZ57RkiVLdM899ygpKUmtWrVSdna2fvzxRy1YsEAHDhwo9iE/KSkp6tmzpzp27Kj//u//1smTJ/XGG2+oSZMmLtvsahITEx2/wTF+/PhCp1mxYoVM09R9991XrPpwE/HGpZkAXzBq1ChTkhkbG1vgvoULF5qSzPLly5uXLl0qdP6cnBxz8uTJZqtWrcyQkBAzODjYbNmypTllyhSXS5/mu9qlFQu7PKSz/MtSrl692hw6dKhZsWJFs1y5cubDDz/scvlBZx9//LEpyRw6dOi1N4aTSZMmmdWqVTNtNpvZoUMHc/PmzQUuq2qaly8R2LhxYzMwMNDlcpvXc1nVQYMGWW4b58savvfee2aDBg1Mm81mNmrUyJw5c6bjconOdu3aZXbq1MksW7asy2NYXc71zTffNBs1amSWLl3arFq1qvnEE0+Yp06dcpmmc+fOZpMmTQpsr8LWtzCSzOHDh5tz5sxxrMNtt93muHyl1bYxzcu907NnzwKPWdjzUpjjx4+bgYGB5ocfflikdbJaXv46ONu6dasZHx9vlitXzgwODja7dOlirl+/vtB1+u6771zGX3kJz3xFeT7Onj1rPvTQQ2ZYWFiB19Dx48fNwYMHm1WqVDHLlCljNmvWrMAlgPMvfTlx4sRC1/PKy1Lu27fPTExMNCMjI83SpUub1apVM++55x5zwYIFBeYvCR06dDAlmY899lih9xf1tXKty6rme+utt8w6deqYNpvNbN26tblmzZoC/We3282XX37ZrFWrlqO/P//880JfI+vXrzdbtWpllilTpsD2XrFihdmhQwezbNmyZmhoqHnvvfeaP//8s8v8VpdJtXK1y6oW1uuFbRvTNM0zZ86Yo0aNMuvXr2+WKVPGrFKlihkbG2u++uqrjvf94vbWv/71LzMmJsa02Wxm48aNzYULFxbpfTRfTk6OWbFiRbNChQrm+fPnC12X/v37mx07diz0PsA0TdMwTS+eTQbAaxYvXqzevXtrzZo1LpfjxM3p0Ucf1S+//FLoL/QC8F2XLl1SdHS07r33Xr333nsF7k9LS1OdOnX0z3/+kz0MsMQ5DMBN6h//+Ifq1q1b6CFVuPmMHj1a3333Hb/0CviZRYsWKSMjw/IcmClTpqhZs2aEBVwVexiAm8w///lPbd++XSkpKZo6dar+/Oc/e7skAICbbdq0Sdu3b9f48eNVpUqV6/4BS0AiMAA3HcMwVK5cOfXv31/Tp09XYCDXPgAAf5OUlKQ5c+aoRYsWmjVrlpo2bertkuDDfCYwpKSkaOHChdq1a5fKli2r2NhYTZgwodBffMw3a9asAj/SYrPZdOHCBU+XCwAAAPgFnzmHYfXq1Ro+fLg2btyo5cuX6+LFi+revbuys7OvOl9oaKiOHTvmuBXlx3EAAAAAXOYzxyIsXbrUZXjWrFmKiIjQli1b1KlTJ8v5DMNQZGSkp8sDAAAA/JLPBIYr5f/qYaVKla463dmzZ1WrVi3Z7Xa1bNlSL7/8spo0aWI5fU5OjnJychzDdrtdJ0+eVOXKlS1/rh4AAADwJaZp6syZM4qOjlZAwNUPOvKZcxic2e129erVS6dPn9batWstp9uwYYP27NmjW2+9VZmZmXr11Ve1Zs0a/fTTT6pevXqh84wZM8bjv6oKAAAA/BEcPnzY8nNxPp8MDE888YS+/PJLrV279por6OzixYuKiYnRwIEDLX8e/co9DJmZmapZs6YOHjyo0NDQG679Zme32/X777+rSpUq10yzuDnQE3BGP8AZ/QBn9IN7ZWVlqVatWjp9+rQqVKhw1Wl97pCkESNG6PPPP9eaNWuKFRYkqXTp0rrtttu0d+9ey2lsNptsNluB8WFhYQQGN7Db7crNzVVYWBgvdkiiJ+CKfoAz+gHO6Af3yt+GRTnk3me2tmmaGjFihD799FOtXLlSderUKfZj5OXl6ccff1RUVJQHKgQAAAD8j8/sYRg+fLjmzZunxYsXq3z58kpLS5MkVahQQWXLlpUkJSYmqlq1akpJSZEkjRs3Trfffrvq16+v06dPa+LEiTp48KAee+wxr60HAAAA4Et8JjC8/fbbkqS4uDiX8TNnzlRSUpIk6dChQy67qE6dOqUhQ4YoLS1NFStWVKtWrbR+/Xo1bty4pMoGAAAAfJrPBIainJudmprqMjx58mRNnjzZQxUBAAAA/s9nzmEAAAAAUPIIDAAAAAAsERgAAAAAWCIwAAAAALBEYAAAAABgicAAAAAAwBKBAQAAAIAlAgMAAAAASwQGAAAAAJZ85peeb2a1e3XwdgluE2AYiqleVzuP/Cp7EX6921ccWLLO2yUAAAB4BHsYAAAAAFgiMAAAAACwRGAAAAAAYInAAAAAAMASgQEAAACAJQIDAAAAAEsEBgAAAACWCAwAAAAALBEYAAAAAFgiMAAAAACwRGAAAAAAYInAAAAAAMASgQEAAACAJQIDAAAAAEsEBgAAAACWCAwAAAAALBEYAAAAAFgiMAAAAACwRGAAAAAAYInAAAAAAMASgQEAAACAJQIDAAAAAEs+ExhSUlLUpk0blS9fXhEREerdu7d27959zfk++eQTNWrUSEFBQWrWrJn+/e9/l0C1AAAAgH/wmcCwevVqDR8+XBs3btTy5ct18eJFde/eXdnZ2ZbzrF+/XgMHDtSjjz6q77//Xr1791bv3r21Y8eOEqwcAAAA8F2B3i6gqJYuXeoyPGvWLEVERGjLli3q1KlTofNMnTpVCQkJeuaZZyRJ48eP1/Lly/Xmm29q+vTpHq8ZAAAA8HU+ExiulJmZKUmqVKmS5TQbNmxQcnKyy7j4+HgtWrTIcp6cnBzl5OQ4hrOysiRJdrtddrv9Biq+fgGG4ZXlekKAYcgwDL9aJ0le6w1/YLfbZZom2xCS6Ae4oh/gjH5wr+JsR58MDHa7XSNHjlSHDh3UtGlTy+nS0tJUtWpVl3FVq1ZVWlqa5TwpKSkaO3ZsgfEZGRm6cOHC9Rd9A2Kq1/XKcj0hwDBUvUqkDEl20/R2OW6Tnp7u7RJ8lt1uV2ZmpkzTVECAzxwlCQ+hH+CMfoAz+sG9zpw5U+RpfTIwDB8+XDt27NDatWvd/tijRo1y2SuRlZWlGjVqKDw8XKGhoW5fXlHsPPKrV5brCQGGIVPSriO/+lVgiIiI8HYJPstut8swDIWHh/MPAPQDXNAPcEY/uFdQUFCRp/W5wDBixAh9/vnnWrNmjapXr37VaSMjI3X8+HGXccePH1dkZKTlPDabTTabrcD4gIAArzWnP32wlnR5d+J/bv6CN64bYxiGV19j+GOhH+CMfoAz+sF9irMNfWZrm6apESNG6NNPP9XKlStVp06da87Tvn17ff311y7jli9frvbt23uqTAAAAMCv+MwehuHDh2vevHlavHixypcv7zgPoUKFCipbtqwkKTExUdWqVVNKSook6cknn1Tnzp01adIk9ezZU//85z+1efNmzZgxw2vrAQAAAPgSn9nD8PbbbyszM1NxcXGKiopy3ObPn++Y5tChQzp27JhjODY2VvPmzdOMGTPUvHlzLViwQIsWLbrqidIAAAAA/o/P7GEwi3C8e2pqaoFxffv2Vd++fT1QEQAAAOD/fGYPAwAAAICSR2AAAAAAYInAAAAAAMASgQEAAACAJQIDAAAAAEsEBgAAAACWCAwAAAAALBEYAAAAAFgiMAAAAACwRGAAAAAAYInAAAAAAMASgQEAAACAJQIDAAAAAEsEBgAAAACWCAwAAAAALBEYAAAAAFgiMAAAAACwRGAAAAAAYInAAAAAAMASgQEAAACAJQIDAAAAAEsEBgAAAACWCAwAAAAALBEYAAAAAFgK9HYBAIqndq8O3i7BrQIMQzHV62rnkV9lN01vl+MWB5as83YJAAC4DXsYAAAAAFgiMAAAAACwRGAAAAAAYInAAAAAAMASgQEAAACAJQIDAAAAAEs+FRjWrFmje++9V9HR0TIMQ4sWLbrq9KmpqTIMo8AtLS2tZAoGAAAAfJxPBYbs7Gw1b95c06ZNK9Z8u3fv1rFjxxy3iIgID1UIAAAA+Bef+uG2Hj16qEePHsWeLyIiQmFhYe4vCAAAAPBzPhUYrleLFi2Uk5Ojpk2basyYMerQwfqXcnNycpSTk+MYzsrKkiTZ7XbZ7XaP11qYAMPwynI9IeA/h4X50zpJKtHe8Ldt54894a33Cn9gt9tlmibbEJLoB7iiH9yrONvRrwNDVFSUpk+frtatWysnJ0fvvvuu4uLitGnTJrVs2bLQeVJSUjR27NgC4zMyMnThwgVPl1yomOp1vbJcTwgwDFWvEilDkt00vV2O26Snp5fYsvypHyT/7ImS7Ad/Y7fblZmZKdM0FRDgU0fNwgPoBzijH9zrzJkzRZ7WrwNDw4YN1bBhQ8dwbGys9u3bp8mTJ+vDDz8sdJ5Ro0YpOTnZMZyVlaUaNWooPDxcoaGhHq+5MDuP/OqV5XpCgGHIlLTryK9+8+FQUomeF+NP/SD5Z09wntT1s9vtMgxD4eHhfCAA/QAX9IN7BQUFFXlavw4MhWnbtq3Wrl1reb/NZpPNZiswPiAgwGvN6S8fovKZpin7f27+oiR7w5+2Wz5/6wn+kd0YwzC8+p6LPxb6Ac7oB/cpzja86bb2tm3bFBUV5e0yAAAAAJ/gU3sYzp49q7179zqG9+/fr23btqlSpUqqWbOmRo0apaNHj+qDDz6QJE2ZMkV16tRRkyZNdOHCBb377rtauXKlvvrqK2+tAgAAAOBTfCowbN68WV26dHEM559rMGjQIM2aNUvHjh3ToUOHHPfn5ubqqaee0tGjRxUcHKxbb71VK1ascHkMAAAAANZ8KjDExcXJvMoxzrNmzXIZfvbZZ/Xss896uCoAAADAf9105zAAAAAAKDoCAwAAAABLBAYAAAAAlggMAAAAACwRGAAAAABYIjAAAAAAsERgAAAAAGCJwAAAAADAEoEBAAAAgCUCAwAAAABLBAYAAAAAlggMAAAAACwRGAAAAABYIjAAAAAAsERgAAAAAGCJwAAAAADAEoEBAAAAgCUCAwAAAABLBAYAAAAAlggMAAAAACwRGAAAAABYIjAAAAAAsERgAAAAAGCJwAAAAADAEoEBAAAAgCUCAwAAAABLBAYAAAAAlgK9XQAA4PrV7tXB2yW4VYBhKKZ6Xe088qvspuntctziwJJ13i4BAG4IexgAAAAAWCIwAAAAALBEYAAAAABgicAAAAAAwJJPBYY1a9bo3nvvVXR0tAzD0KJFi645T2pqqlq2bCmbzab69etr1qxZHq8TAAAA8Bc+FRiys7PVvHlzTZs2rUjT79+/Xz179lSXLl20bds2jRw5Uo899piWLVvm4UoBAAAA/+BTl1Xt0aOHevToUeTpp0+frjp16mjSpEmSpJiYGK1du1aTJ09WfHy8p8oEAAAA/IZPBYbi2rBhg7p16+YyLj4+XiNHjrScJycnRzk5OY7hrKwsSZLdbpfdbvdIndcSYBheWa4nBBiGDMPwq3WSVKK94W/bzh97gn64fvQDnNntdpmmyTaEJPrB3YqzHf06MKSlpalq1aou46pWraqsrCydP39eZcuWLTBPSkqKxo4dW2B8RkaGLly44LFaryamel2vLNcTAgxD1atEypD85keZJCk9Pb3EluVP/SD5Z0/QD9ePfoAzu92uzMxMmaapgACfOooaHkA/uNeZM2eKPK1fB4brMWrUKCUnJzuGs7KyVKNGDYWHhys0NNQrNe088qtXlusJAYYhU9IuP/oVV0mKiIgosWX5Uz9I/tkT9MP1ox/gzG63yzAMhYeH8wER9IObBQUFFXlavw4MkZGROn78uMu448ePKzQ0tNC9C5Jks9lks9kKjA8ICPBac/rLP818pmnK/p+bvyjJ3vCn7ZbP33qCfrgx9AOcGYbh1f/B+GOhH9ynONvQr7d2+/bt9fXXX7uMW758udq3b++ligAAAADf4lOB4ezZs9q2bZu2bdsm6fJlU7dt26ZDhw5Junw4UWJiomP6xx9/XL/++queffZZ7dq1S2+99ZY+/vhj/e///q83ygcAAAB8jk8Fhs2bN+u2227TbbfdJklKTk7WbbfdphdffFGSdOzYMUd4kKQ6deroiy++0PLly9W8eXNNmjRJ7777LpdUBQAAAIrIp85hiIuLk3mVY1oL+xXnuLg4ff/99x6sCgAAAPBfPrWHAQAAAEDJIjAAAAAAsERgAAAAAGCJwAAAAADAEoEBAAAAgCUCAwAAAABLBAYAAAAAlggMAAAAACwRGAAAAABYIjAAAAAAsERgAAAAAGCJwAAAAADAEoEBAAAAgCUCAwAAAABLBAYAAAAAlggMAAAAACwRGAAAAABYIjAAAAAAsERgAAAAAGCJwAAAAADAEoEBAAAAgCUCAwAAAABLBAYAAAAAlggMAAAAACwVOTB8/PHHys3NdQwfOXJEdrvdMXzu3Dm98sor7q0OAAAAgFcVOTAMHDhQp0+fdgw3btxYBw4ccAyfOXNGo0aNcmdtAAAAALysyIHBNM2rDgMAAADwP5zDAAAAAMASgQEAAACApcDiTLxs2TJVqFBBkmS32/X1119rx44dkuRyfgMAAAAA/1CswDBo0CCX4WHDhrkMG4Zx4xUBAAAA+MMocmBwvoQqAAAAgJuDz53DMG3aNNWuXVtBQUFq166dvv32W8tpZ82aJcMwXG5BQUElWC0AAADg24odGD755BPdf//9atq0qZo2bar7779fCxYs8ERtBcyfP1/JyckaPXq0tm7dqubNmys+Pl7p6emW84SGhurYsWOO28GDB0ukVgAAAMAfFDkw2O129e/fX/3799fPP/+s+vXrq379+vrpp5/Uv39/DRgwwOO/zfDaa69pyJAhGjx4sBo3bqzp06crODhY77//vuU8hmEoMjLScatatapHawQAAAD8SZHPYZg6dapWrFihJUuW6J577nG5b8mSJRo8eLCmTp2qkSNHurtGSVJubq62bNni8mvSAQEB6tatmzZs2GA539mzZ1WrVi3Z7Xa1bNlSL7/8spo0aWI5fU5OjnJychzDWVlZki4HJm+dxxHgRyeTB/zn0DB/WiepZM/x8bdt5489QT9cP/oBzux2u0zTZBtCEv3gbsXZjkUODDNnztTEiRMLhAVJ6tWrl1555RWPBobff/9deXl5BfYQVK1aVbt27Sp0noYNG+r999/XrbfeqszMTL366quKjY3VTz/9pOrVqxc6T0pKisaOHVtgfEZGhi5cuHDjK3IdYqrX9cpyPSHAMFS9SqQMSXY/+rXwqx0W527+1A+Sf/YE/XD96Ac4s9vtyszMlGmaCgjwudMu4Wb0g3udOXOmyNMWOTDs2bNH3bp1s7y/W7duGjFiRJEXXBLat2+v9u3bO4ZjY2MVExOjd955R+PHjy90nlGjRik5OdkxnJWVpRo1aig8PFyhoaEer7kwO4/86pXlekKAYciUtOvIr37zYUCSIiIiSmxZ/tQPkn/2BP1w/egHOLPb7TIMQ+Hh4XxABP3gZsW5EFCRA0PZsmV1+vRp1axZs9D7s7KyPHoFoipVqqhUqVI6fvy4y/jjx48rMjKySI9RunRp3Xbbbdq7d6/lNDabTTabrcD4gIAArzWnv/zTzGeapuz/ufmLkuwNf9pu+fytJ+iHG0M/wJlhGF79H4w/FvrBfYqzDYs8Zfv27fX2229b3j9t2jSXb/PdrUyZMmrVqpW+/vprx7j8X5su6nLz8vL0448/KioqylNlAgAAAH6lyHsY/va3vykuLk4nTpzQ008/rUaNGsk0Te3cuVOTJk3S4sWLtWrVKk/WquTkZA0aNEitW7dW27ZtNWXKFGVnZ2vw4MGSpMTERFWrVk0pKSmSpHHjxun2229X/fr1dfr0aU2cOFEHDx7UY4895tE6AQAAAH9R5MAQGxur+fPna+jQofrXv/7lcl/FihX10UcfqUOHDm4v0Fn//v2VkZGhF198UWlpaWrRooWWLl3qOBH60KFDLrtXTp06pSFDhigtLU0VK1ZUq1attH79ejVu3NijdQIAAAD+osiBYf/+/erTp4/i4+P11Vdf6ZdffpEk3XLLLerevbuCg4M9VqSzESNGWJ5cnZqa6jI8efJkTZ48uQSqAgAAAPxTkQNDvXr1VKtWLXXp0kVdunTRQw89ZHlpUgAAAAD+ociBYeXKlUpNTVVqaqo++ugj5ebmqm7duuratasjRPArygAAAIB/KXJgiIuLU1xcnCTpwoULWr9+vSNAzJ49WxcvXlSjRo30008/eapWAAAAACWsyIHBWVBQkLp27aqOHTuqS5cu+vLLL/XOO+9Y/uIyAAAAAN9UrMCQm5urjRs3atWqVUpNTdWmTZtUo0YNderUSW+++aY6d+7sqToBAAAAeEGRA0PXrl21adMm1alTR507d9awYcM0b948fgQNAAAA8GNFDgzffPONoqKi1LVrV8XFxalz586qXLmyJ2sDAAAA4GUB157kstOnT2vGjBkKDg7WhAkTFB0drWbNmmnEiBFasGCBMjIyPFknAAAAAC8o8h6GkJAQJSQkKCEhQZJ05swZrV27VqtWrdIrr7yihx9+WA0aNNCOHTs8ViwAAACAklXkPQxXCgkJUaVKlVSpUiVVrFhRgYGB2rlzpztrAwAAAOBlRd7DYLfbtXnzZqWmpmrVqlVat26dsrOzVa1aNXXp0kXTpk1Tly5dPFkrAAAAgBJW5MAQFham7OxsRUZGqkuXLpo8ebLi4uJUr149T9YHAAAAwIuKHBgmTpyoLl266JZbbvFkPQAAAAD+QIocGIYNG+bJOgAAAAD8AV33Sc8AAAAA/B+BAQAAAIAlAgMAAAAASwQGAAAAAJYIDAAAAAAsERgAAAAAWCIwAAAAALBEYAAAAABgicAAAAAAwBKBAQAAAIAlAgMAAAAASwQGAAAAAJYIDAAAAAAsERgAAAAAWCIwAAAAALBEYAAAAABgicAAAAAAwFKgtwsormnTpmnixIlKS0tT8+bN9cYbb6ht27aW03/yySd64YUXdODAATVo0EATJkzQ3XffXYIVAwBQMmr36uDtEtwqwDAUU72udh75VXbT9HY5bnNgyTpvlwAUi0/tYZg/f76Sk5M1evRobd26Vc2bN1d8fLzS09MLnX79+vUaOHCgHn30UX3//ffq3bu3evfurR07dpRw5QAAAIBv8qnA8Nprr2nIkCEaPHiwGjdurOnTpys4OFjvv/9+odNPnTpVCQkJeuaZZxQTE6Px48erZcuWevPNN0u4cgAAAMA3+cwhSbm5udqyZYtGjRrlGBcQEKBu3bppw4YNhc6zYcMGJScnu4yLj4/XokWLLJeTk5OjnJwcx3BWVpYk6fTp07Lb7TewBjfgUp53lusJhqFLFy9eXic/2r18+vTpkluYP/WD5Jc9QT/cAPrhxtAPPqFEe8KP2O12ZWVlqUyZMgoI8KnvvP+Q8j/jFoXPBIbff/9deXl5qlq1qsv4qlWrateuXYXOk5aWVuj0aWlplstJSUnR2LFjC4yvVavWdVSNwuzXRm+X4HYVK1b0dgk+zd96gn64MfQDnPlbP0j0BHyPzwSGkjJq1CiXvRJZWVmqUaOGDh48qNDQUC9W5h/sdrt+//13ValShW8HIImegCv6Ac7ohxuzbHBzb5fgVqYRIEXeIqX9IsP00lEfHhA/8wevLDcrK6vIX4j7TGCoUqWKSpUqpePHj7uMP378uCIjIwudJzIysljTS5LNZpPNZiswPiwsjMDgBna7Xbm5uQoLC+PNH5LoCbiiH+CMfrgxIT7zKa9oTENSmUApUDL85wg1hYWFeWW5xXlN+cyrr0yZMmrVqpW+/vprxzi73a6vv/5a7du3L3Se9u3bu0wvScuXL7ecHgAAAIArn8qeycnJGjRokFq3bq22bdtqypQpys7O1uDBgyVJiYmJqlatmlJSUiRJTz75pDp37qxJkyapZ8+e+uc//6nNmzdrxowZ3lwNAAAAwGf4VGDo37+/MjIy9OKLLyotLU0tWrTQ0qVLHSc2Hzp0yGX3SmxsrObNm6fnn39ef/3rX9WgQQMtWrRITZs29dYqAAAAAD7FpwKDJI0YMUIjRowo9L7U1NQC4/r27au+fft6uCoAAADAP/nMOQwAAAAASh6BAQAAAIAlAgMAAAAASwQGAAAAAJYIDAAAAAAsERgAAAAAWCIwAAAAALBEYAAAAABgicAAAAAAwBKBAQAAAIAlAgMAAAAASwQGAAAAAJYIDAAAAAAsERgAAAAAWCIwAAAAALBEYAAAAABgicAAAAAAwBKBAQAAAIAlAgMAAAAASwQGAAAAAJYIDAAAAAAsERgAAAAAWCIwAAAAALBEYAAAAABgicAAAAAAwBKBAQAAAIAlAgMAAAAASwQGAAAAAJYIDAAAAAAsERgAAAAAWCIwAAAAALBEYAAAAABgyWcCw8mTJ/Xwww8rNDRUYWFhevTRR3X27NmrzhMXFyfDMFxujz/+eAlVDAAAAPi+QG8XUFQPP/ywjh07puXLl+vixYsaPHiwhg4dqnnz5l11viFDhmjcuHGO4eDgYE+XCgAAAPgNnwgMO3fu1NKlS/Xdd9+pdevWkqQ33nhDd999t1599VVFR0dbzhscHKzIyMiSKhUAAADwKz4RGDZs2KCwsDBHWJCkbt26KSAgQJs2bVKfPn0s5507d67mzJmjyMhI3XvvvXrhhReuupchJydHOTk5juGsrCxJkt1ul91ud8Pa3NzsdrtM02RbwoGegDP6Ac7ohxtjGj5z5HmRmEaAZBiSn62Xt/q7OMv1icCQlpamiIgIl3GBgYGqVKmS0tLSLOd76KGHVKtWLUVHR2v79u36y1/+ot27d2vhwoWW86SkpGjs2LEFxmdkZOjChQvXvxKQdLk5MzMzZZqmAgL86wWP60NPwBn9AGf0ww2KbuztCtzLMKRKNS7/bZrercWN0tPTvbLcM2fOFHlarwaG5557ThMmTLjqNDt37rzuxx86dKjj72bNmikqKkp33nmn9u3bp3r16hU6z6hRo5ScnOwYzsrKUo0aNRQeHq7Q0NDrrgWX2e12GYah8PBw3vwhiZ6AK/oBzuiHG/Tbz96uwL3y9yz8tlMy/Wev05VfipeUoKCgIk/r1cDw1FNPKSkp6arT1K1bV5GRkQXS16VLl3Ty5MlinZ/Qrl07SdLevXstA4PNZpPNZiswPiAggDcrNzEMg+0JF/QEnNEPcEY/XD/Djz5UO5imZNr9at281dvFWa5XA0N4eLjCw8OvOV379u11+vRpbdmyRa1atZIkrVy5Una73RECimLbtm2SpKioqOuqFwAAALjZ+ERcj4mJUUJCgoYMGaJvv/1W69at04gRIzRgwADHFZKOHj2qRo0a6dtvv5Uk7du3T+PHj9eWLVt04MABLVmyRImJierUqZNuvfVWb64OAAAA4DN8IjBIl6921KhRI9155526++671bFjR82YMcNx/8WLF7V7926dO3dOklSmTBmtWLFC3bt3V6NGjfTUU0/pgQce0GeffeatVQAAAAB8jk9cJUmSKlWqdNUfaatdu7ZMpzPma9SoodWrV5dEaQAAAIDf8pk9DAAAAABKHoEBAAAAgCUCAwAAAABLBAYAAAAAlggMAAAAACwRGAAAAABYIjAAAAAAsERgAAAAAGCJwAAAAADAEoEBAAAAgCUCAwAAAABLBAYAAAAAlggMAAAAACwRGAAAAABYIjAAAAAAsERgAAAAAGCJwAAAAADAEoEBAAAAgCUCAwAAAABLBAYAAAAAlggMAAAAACwRGAAAAABYIjAAAAAAsERgAAAAAGCJwAAAAADAEoEBAAAAgCUCAwAAAABLBAYAAAAAlggMAAAAACwRGAAAAABYIjAAAAAAsERgAAAAAGDJZwLD3//+d8XGxio4OFhhYWFFmsc0Tb344ouKiopS2bJl1a1bN+3Zs8ezhQIAAAB+xGcCQ25urvr27asnnniiyPO88sorev311zV9+nRt2rRJISEhio+P14ULFzxYKQAAAOA/Ar1dQFGNHTtWkjRr1qwiTW+apqZMmaLnn39e9913nyTpgw8+UNWqVbVo0SINGDDAU6UCAAAAfsNnAkNx7d+/X2lpaerWrZtjXIUKFdSuXTtt2LDBMjDk5OQoJyfHMZyVlSVJstvtstvtni36JmC322WaJtsSDvQEnNEPcEY/3BjT8JkDSYrENAIkw5D8bL281d/FWa7fBoa0tDRJUtWqVV3GV61a1XFfYVJSUhx7M5xlZGRwKJMb2O12ZWZmyjRNBQT41wse14eegDP6Ac7ohxsU3djbFbiXYUiValz+2zS9W4sbpaene2W5Z86cKfK0Xg0Mzz33nCZMmHDVaXbu3KlGjRqVUEXSqFGjlJyc7BjOyspSjRo1FB4ertDQ0BKrw1/Z7XYZhqHw8HDe/CGJnoAr+gHO6Icb9NvP3q7AvfL3LPy2UzL9Z69TRESEV5YbFBRU5Gm9GhieeuopJSUlXXWaunXrXtdjR0ZGSpKOHz+uqKgox/jjx4+rRYsWlvPZbDbZbLYC4wMCAnizchPDMNiecEFPwBn9AGf0w/Uz/OhDtYNpSqbdr9bNW71dnOV6NTCEh4crPDzcI49dp04dRUZG6uuvv3YEhKysLG3atKlYV1oCAAAAbmY+E9cPHTqkbdu26dChQ8rLy9O2bdu0bds2nT171jFNo0aN9Omnn0q6/I3EyJEj9dJLL2nJkiX68ccflZiYqOjoaPXu3dtLawEAAAD4Fp856fnFF1/U7NmzHcO33XabJGnVqlWKi4uTJO3evVuZmZmOaZ599lllZ2dr6NChOn36tDp27KilS5cW65gtAAAA4GbmM4Fh1qxZ1/wNBvOKM+YNw9C4ceM0btw4D1YGAAAA+C+fOSQJAAAAQMkjMAAAAACwRGAAAAAAYInAAAAAAMASgQEAAACAJQIDAAAAAEsEBgAAAACWCAwAAAAALBEYAAAAAFgiMAAAAACwRGAAAAAAYInAAAAAAMASgQEAAACAJQIDAAAAAEsEBgAAAACWCAwAAAAALBEYAAAAAFgiMAAAAACwRGAAAAAAYInAAAAAAMASgQEAAACAJQIDAAAAAEsEBgAAAACWCAwAAAAALBEYAAAAAFgK9HYBAAAAcL97Pjno7RLcym63Kz09XREREQoI4DvvksTWBgAAAGCJwAAAAADAEoEBAAAAgCUCAwAAAABLBAYAAAAAlnwmMPz9739XbGysgoODFRYWVqR5kpKSZBiGyy0hIcGzhQIAAAB+xGcuq5qbm6u+ffuqffv2eu+994o8X0JCgmbOnOkYttlsnigPAAAA8Es+ExjGjh0rSZo1a1ax5rPZbIqMjPRARQAAAID/85nAcL1SU1MVERGhihUrqmvXrnrppZdUuXJly+lzcnKUk5PjGM7KypJ0+cdC7Ha7x+v1d3a7XaZpsi3hQE/AGf0AZ/QDnNEP7lWc7ejXgSEhIUH333+/6tSpo3379umvf/2revTooQ0bNqhUqVKFzpOSkuLYm+EsIyNDFy5c8HTJfs9utyszM1OmafIrjZBET8AV/QBn9AOc0Q/udebMmSJP69XA8Nxzz2nChAlXnWbnzp1q1KjRdT3+gAEDHH83a9ZMt956q+rVq6fU1FTdeeedhc4zatQoJScnO4azsrJUo0YNhYeHKzQ09LrqwP+x2+0yDEPh4eG82CGJnoAr+gHO6Ac4ox/cKygoqMjTejUwPPXUU0pKSrrqNHXr1nXb8urWrasqVapo7969loHBZrMVemJ0QEAAzekmhmGwPeGCnoAz+gHO6Ac4ox/cpzjb0KuBITw8XOHh4SW2vCNHjujEiROKiooqsWUCAAAAvsxnzmE4dOiQTp48qUOHDikvL0/btm2TJNWvX1/lypWTJDVq1EgpKSnq06ePzp49q7Fjx+qBBx5QZGSk9u3bp2effVb169dXfHx8kZdrmqak/zv5GTfGbrfrzJkzCgoK4tsBSKIn4Ip+gDP6Ac7oB/fK/2yb/1n3qkwfMWjQIFNSgduqVasc00gyZ86caZqmaZ47d87s3r27GR4ebpYuXdqsVauWOWTIEDMtLa1Yyz18+HChy+XGjRs3bty4cePGzddvhw8fvubnYeM/H7RhwW6367ffflP58uVlGIa3y/F5+SeRHz58mJPIIYmegCv6Ac7oBzijH9zLNE2dOXNG0dHR19xj4zOHJHlLQECAqlev7u0y/E5oaCgvdrigJ+CMfoAz+gHO6Af3qVChQpGm4wAwAAAAAJYIDAAAAAAsERhQomw2m0aPHl3ob13g5kRPwBn9AGf0A5zRD97DSc8AAAAALLGHAQAAAIAlAgMAAAAASwQGAAAAAJYIDAAAAAAsERhQbHFxcRo5cuRVp6ldu7amTJlSIvXAu+gHOKMf4Ix+wJXoCd9EYIBHfPfddxo6dOgNPcbChQvVvXt3Va5cWYZhaNu2be4pDiXuRvvh4sWL+stf/qJmzZopJCRE0dHRSkxM1G+//ebGKlFS3PH+MGbMGDVq1EghISGqWLGiunXrpk2bNrmpQpQkd/SDs8cff1yGYfCB04e5oyeSkpJkGIbLLSEhwU0V3nwIDPCI8PBwBQcH39BjZGdnq2PHjpowYYKbqoK33Gg/nDt3Tlu3btULL7ygrVu3auHChdq9e7d69erlxipRUtzx/nDLLbfozTff1I8//qi1a9eqdu3a6t69uzIyMtxUJUqKO/oh36effqqNGzcqOjraLY8H73BXTyQkJOjYsWOO20cffeSG6m5OBAZcl0uXLmnEiBGqUKGCqlSpohdeeEHOP+lx5e5EwzD07rvvqk+fPgoODlaDBg20ZMmSqy7jkUce0Ysvvqhu3bp5ajXgJp7uhwoVKmj58uXq16+fGjZsqNtvv11vvvmmtmzZokOHDnly1XAdSuL94aGHHlK3bt1Ut25dNWnSRK+99pqysrK0fft2T60WrlNJ9IMkHT16VP/zP/+juXPnqnTp0p5YFbhJSfWEzWZTZGSk41axYkVPrM5NgcCA6zJ79mwFBgbq22+/1dSpU/Xaa6/p3Xffveo8Y8eOVb9+/bR9+3bdfffdevjhh3Xy5MkSqhie5I1+yMzMlGEYCgsLu8Hq4W4l3Q+5ubmaMWOGKlSooObNm7tjFeBGJdEPdrtdjzzyiJ555hk1adLE3asANyup94jU1FRFRESoYcOGeuKJJ3TixAl3rsbNxQSKqXPnzmZMTIxpt9sd4/7yl7+YMTExjuFatWqZkydPdgxLMp9//nnH8NmzZ01J5pdffnnN5e3fv9+UZH7//fduqR/uVdL9YJqmef78ebNly5bmQw89dOMrALcqyX747LPPzJCQENMwDDM6Otr89ttv3bcicIuS6oeXX37ZvOuuuxzLufIx8cdRUj3x0UcfmYsXLza3b99ufvrpp2ZMTIzZpk0b89KlS+5doZsEexhwXW6//XYZhuEYbt++vfbs2aO8vDzLeW699VbH3yEhIQoNDVV6erpH60TJKMl+uHjxovr16yfTNPX222/fWOHwiJLqhy5dumjbtm1av369EhIS1K9fP95T/oA83Q9btmzR1KlTNWvWLJfl4I+rJN4jBgwYoF69eqlZs2bq3bu3Pv/8c3333XdKTU11yzrcbAgMKDFXHlNqGIbsdruXqoG3XU8/5IeFgwcPavny5QoNDfVkiShB19MPISEhql+/vm6//Xa99957CgwM1HvvvefJMlFCitMP33zzjdLT01WzZk0FBgYqMDBQBw8e1FNPPaXatWuXQLUoCTf6GaJu3bqqUqWK9u7d6+7SbgqB3i4AvunKyxdu3LhRDRo0UKlSpbxUEbypJPohPyzs2bNHq1atUuXKld322HAvb70/2O125eTkeHQZKD5P98MjjzxS4OIY8fHxeuSRRzR48GC3LAPu5Y33iCNHjujEiROKiory2DL8GYEB1+XQoUNKTk7WsGHDtHXrVr3xxhuaNGmSW5dx8uRJHTp0yHGt/d27d0uS42oH+OPwdD9cvHhRDz74oLZu3arPP/9ceXl5SktLkyRVqlRJZcqUcduycOM83Q/Z2dn6+9//rl69eikqKkq///67pk2bpqNHj6pv375uWw7cw9P9ULly5QJfIJQuXVqRkZFq2LCh25YD9/F0T5w9e1Zjx47VAw88oMjISO3bt0/PPvus6tevr/j4eLct52ZCYMB1SUxM1Pnz59W2bVuVKlVKTz75pFt/eEeSlixZ4vLt0IABAyRJo0eP1pgxY9y6LNwYT/fD0aNHHZfQa9Gihct9q1atUlxcnNuWhRvn6X4oVaqUdu3apdmzZ+v3339X5cqV1aZNG33zzTdcIecPqCT+X8C3lMR7xPbt2zV79mydPn1a0dHR6t69u8aPHy+bzea25dxMDNN0uvAtAAAAADjhpGcAAAAAlggMAAAAACwRGAAAAABYIjAAAAAAsERgAAAAAGCJwAAAAADAEoEBAAAAgCUCAwAAAABLBAYAAAAAlggMAAAAACwRGAAAAABY+v91+C2cCGrtDAAAAABJRU5ErkJggg==", "text/plain": [ "
" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "def fit_woe_numeric(x, y, n_bins=5, smooth=20):\n", " \"\"\"Quantile-bin a numeric variable and return (bin_edges, woe_per_bin).\"\"\"\n", " edges = np.unique(np.quantile(x, np.linspace(0, 1, n_bins + 1)))\n", " if len(edges) < 3: # (near-)constant variable -> single bin\n", " edges = np.array([-np.inf, np.inf])\n", " else:\n", " edges[0], edges[-1] = -np.inf, np.inf\n", " bins = np.digitize(x, edges[1:-1])\n", " n_good, n_bad = (y == 0).sum(), (y == 1).sum()\n", " prior = n_good / n_bad\n", " woe = np.zeros(len(edges) - 1)\n", " for b in range(len(edges) - 1):\n", " m = bins == b\n", " good_b, bad_b = (y[m] == 0).sum(), (y[m] == 1).sum()\n", " # smoothed odds, shrunk toward the portfolio prior\n", " woe[b] = np.log(((good_b + smooth) / (bad_b + smooth / prior)) / prior)\n", " return edges, woe\n", "\n", "def apply_woe_numeric(x, edges, woe):\n", " return woe[np.digitize(x, edges[1:-1])]\n", "\n", "def fit_woe_categorical(x, y, smooth=20):\n", " n_good, n_bad = (y == 0).sum(), (y == 1).sum()\n", " prior = n_good / n_bad\n", " mapping = {}\n", " for cat in np.unique(x):\n", " m = x == cat\n", " good_c, bad_c = (y[m] == 0).sum(), (y[m] == 1).sum()\n", " mapping[cat] = np.log(((good_c + smooth) / (bad_c + smooth / prior)) / prior)\n", " return mapping\n", "\n", "def apply_woe_categorical(x, mapping):\n", " return np.array([mapping.get(v, 0.0) for v in x]) # unseen category -> neutral evidence\n", "\n", "# Look at the utilization WOE pattern — it should be monotone decreasing\n", "edges, woe = fit_woe_numeric(features[\"util_now\"].values, y, n_bins=5)\n", "labels = [f\"bin {i+1}\" for i in range(len(woe))]\n", "plt.bar(labels, woe, color=[\"#1f4d3a\" if w >= 0 else \"#b3541e\" for w in woe])\n", "plt.axhline(0, color=\"black\", lw=0.8)\n", "plt.title(\"WOE by utilization bin (monotone = validator-friendly)\")\n", "plt.ylabel(\"WOE\"); plt.show()" ] }, { "cell_type": "markdown", "id": "03a9d317", "metadata": {}, "source": [ "## 6. The leakage experiment\n", "\n", "`region` is **pure noise** — it was generated independently of everything. A correct pipeline should\n", "find it worthless. Now compare two pipelines on a fixed train/test split:\n", "\n", "- **LEAKY**: fit every WOE mapping on **all loans (train + test)**, then train on train, evaluate on test.\n", "- **CLEAN**: fit every WOE mapping on **train only**.\n", "\n", "The leaky pipeline lets each test loan's own outcome sneak into the WOE value of its region — the\n", "encoding partially memorizes the test set." ] }, { "cell_type": "code", "execution_count": 7, "id": "e06455f8", "metadata": { "execution": { "iopub.execute_input": "2026-07-14T00:49:48.620161Z", "iopub.status.busy": "2026-07-14T00:49:48.619938Z", "iopub.status.idle": "2026-07-14T00:49:48.667446Z", "shell.execute_reply": "2026-07-14T00:49:48.667052Z" } }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "test AUC, LEAKY pipeline (WOE fitted on train+test): 0.8387\n", "test AUC, CLEAN pipeline (WOE fitted on train only): 0.7877\n", "leakage bonus: +0.0510 AUC — from a pure-noise variable\n" ] } ], "source": [ "num_cols = [c for c in features.columns if c not in (\"region\",)]\n", "split = rng.random(len(ids)) < 0.7\n", "Xn = features[num_cols].values\n", "reg = features[\"region\"].values\n", "\n", "def run_pipeline(fit_mask):\n", " \"\"\"Fit WOE mappings on rows where fit_mask is True; train on train; return test AUC.\"\"\"\n", " cols = []\n", " for j in range(Xn.shape[1]):\n", " e, w = fit_woe_numeric(Xn[fit_mask, j], y[fit_mask], n_bins=5)\n", " cols.append(apply_woe_numeric(Xn[:, j], e, w))\n", " mapping = fit_woe_categorical(reg[fit_mask], y[fit_mask], smooth=5)\n", " cols.append(apply_woe_categorical(reg, mapping))\n", " X = np.column_stack(cols)\n", " model = LogisticRegression(max_iter=2000).fit(X[split], y[split])\n", " return roc_auc_score(y[~split], model.predict_proba(X[~split])[:, 1])\n", "\n", "auc_leaky = run_pipeline(fit_mask=np.ones(len(ids), dtype=bool)) # WOE fitted on ALL rows\n", "auc_clean = run_pipeline(fit_mask=split) # WOE fitted on train only\n", "\n", "print(f\"test AUC, LEAKY pipeline (WOE fitted on train+test): {auc_leaky:.4f}\")\n", "print(f\"test AUC, CLEAN pipeline (WOE fitted on train only): {auc_clean:.4f}\")\n", "print(f\"leakage bonus: {auc_leaky - auc_clean:+.4f} AUC — from a pure-noise variable\")" ] }, { "cell_type": "markdown", "id": "f94354d1", "metadata": {}, "source": [ "The gap is the **leakage bonus**: performance that exists only because the encoding saw the test\n", "outcomes. In production there is no test set to leak from — the bonus evaporates and the model\n", "underdelivers vs. its validation report. With higher-cardinality categoricals, smaller samples, or\n", "less smoothing, the gap grows; try `smooth=1` and `n_bins=20` above and watch it widen.\n", "\n", "The standard fix in scorecard shops is **out-of-fold encoding**: split the training data into K folds,\n", "encode each fold with a mapping fitted on the other K−1, and fit the final mapping on all training data\n", "only for scoring genuinely new loans." ] }, { "cell_type": "markdown", "id": "f2759b0b", "metadata": {}, "source": [ "## 7. Which features carry the signal?\n", "\n", "A quick univariate scan — each feature's solo AUC on the test set (clean encodings). Note how the\n", "**max/min window features and the trend** compete with the level features, and that `region` sits at\n", "~0.50, exactly where a noise variable belongs." ] }, { "cell_type": "code", "execution_count": 8, "id": "107a1ab1", "metadata": { "execution": { "iopub.execute_input": "2026-07-14T00:49:48.670565Z", "iopub.status.busy": "2026-07-14T00:49:48.670392Z", "iopub.status.idle": "2026-07-14T00:49:48.711848Z", "shell.execute_reply": "2026-07-14T00:49:48.711295Z" } }, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
featuresolo test AUC
0payr_min_6m0.8050
1payr_min_3m0.7956
2payr_min_12m0.7931
3pay_ratio_now0.7873
4util_mean_3m0.7855
5util_now0.7808
6util_lag30.7799
7util_mean_12m0.7737
8util_mean_6m0.7735
9util_lag60.7727
10util_lag120.7645
11dpd_max_12m0.6796
12util_trend_6m0.6675
13region (noise)0.5428
14util_degenerate0.5000
15dpd_now0.5000
16dpd_max_3m0.5000
17dpd_max_6m0.5000
\n", "
" ], "text/plain": [ " feature solo test AUC\n", "0 payr_min_6m 0.8050\n", "1 payr_min_3m 0.7956\n", "2 payr_min_12m 0.7931\n", "3 pay_ratio_now 0.7873\n", "4 util_mean_3m 0.7855\n", "5 util_now 0.7808\n", "6 util_lag3 0.7799\n", "7 util_mean_12m 0.7737\n", "8 util_mean_6m 0.7735\n", "9 util_lag6 0.7727\n", "10 util_lag12 0.7645\n", "11 dpd_max_12m 0.6796\n", "12 util_trend_6m 0.6675\n", "13 region (noise) 0.5428\n", "14 util_degenerate 0.5000\n", "15 dpd_now 0.5000\n", "16 dpd_max_3m 0.5000\n", "17 dpd_max_6m 0.5000" ] }, "execution_count": 8, "metadata": {}, "output_type": "execute_result" } ], "source": [ "rows = []\n", "for j, c in enumerate(num_cols):\n", " e, w = fit_woe_numeric(Xn[split, j], y[split], n_bins=5)\n", " v = apply_woe_numeric(Xn[~split, j], e, w)\n", " auc = roc_auc_score(y[~split], -v) # negative: high WOE = safe\n", " rows.append((c, max(auc, 1 - auc)))\n", "mapping = fit_woe_categorical(reg[split], y[split])\n", "v = apply_woe_categorical(reg[~split], mapping)\n", "auc = roc_auc_score(y[~split], -v)\n", "rows.append((\"region (noise)\", max(auc, 1 - auc)))\n", "\n", "pd.DataFrame(rows, columns=[\"feature\", \"solo test AUC\"]) \\\n", " .sort_values(\"solo test AUC\", ascending=False).reset_index(drop=True).round(4)" ] }, { "cell_type": "markdown", "id": "5930421e", "metadata": {}, "source": [ "## Exercises\n", "\n", "1. **Boundary leakage, quantified.** Rebuild the rolling features with windows ending at month 24\n", " (the snapshot month itself) instead of 23, and — worse — at month 26. How much AUC do you \"gain\"?\n", " That gain is fictional; explain to an imaginary validator why.\n", "2. **Window-length tournament.** Generate 1/3/6/12/24-month versions of `util_mean` and `dpd_max`;\n", " which lengths survive a forward-selection loop?\n", "3. **Monotone binning.** Modify `fit_woe_numeric` to merge adjacent bins until WOE is monotone.\n", " How much in-sample IV do you lose, and does *test* AUC actually drop?\n", "4. **Out-of-fold encoding.** Implement K-fold WOE encoding inside the training set and compare its\n", " cross-validated AUC to the naive within-fold version — the same experiment as Section 6, one level down.\n", "5. **Feature dictionary.** Write the markdown table for every feature built here: name, formula,\n", " window, availability lag, missing rule. (Yes, really. This is the habit that matters most.)" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.11.15" } }, "nbformat": 4, "nbformat_minor": 5 }