Take a model from a table to an endpoint

import pandas as pd

df = pd.read_parquet("raw.parquet")

print(len(df), "rows", df["target"].mean().round(4), "base rate")

# Missingness by column, worst first. Anything above 5% is a decision you
# have to make and record, not a default to fill with zero.
missing = df.isna().mean().sort_values(ascending=False)
print(missing[missing > 0.05])

# Constant and near-constant columns: they carry no signal and will survive
# every feature selection step you write, wasting your time later.
for c in df.columns:
    top = df[c].value_counts(normalize=True, dropna=False).iloc[0]
    if top > 0.99:
        print("near-constant:", c, round(top, 4))

# Duplicated rows, and duplicated entities with different targets, which is
# a labelling problem rather than a deduplication one.
print("exact duplicate rows:", df.duplicated().sum())
conflict = df.groupby("customer_id")["target"].nunique()
print("ids with conflicting labels:", (conflict > 1).sum())

from sklearn.model_selection import GroupShuffleSplit

def time_split(df, ts="applied_at", group="customer_id",
               val_start="2024-07-01", test_start="2024-10-01"):
    "Train on the past, validate on the next quarter, test on the one after."
    train = df[df[ts] < val_start]
    val   = df[(df[ts] >= val_start) & (df[ts] < test_start)]
    test  = df[df[ts] >= test_start]

    # A customer whose first loan is in train and second in test leaks: the
    # model has seen that customer's behaviour and the outcome attached to it.
    leaked = set(train[group]) & set(val[group]) | set(train[group]) & set(test[group])
    train = train[~train[group].isin(leaked)]
    return train, val, test


def test_split_does_not_leak():
    train, val, test = time_split(df)
    assert set(train["customer_id"]).isdisjoint(val["customer_id"])
    assert set(train["customer_id"]).isdisjoint(test["customer_id"])
    assert train["applied_at"].max() < val["applied_at"].min()
    assert val["applied_at"].max() < test["applied_at"].min()

from sklearn.dummy import DummyClassifier
from sklearn.linear_model import LogisticRegression
from sklearn.metrics import average_precision_score, roc_auc_score

def evaluate(model, X, y):
    p = model.predict_proba(X)[:, 1]
    return {"auc": round(roc_auc_score(y, p), 4),
            "ap": round(average_precision_score(y, p), 4)}

dummy = DummyClassifier(strategy="prior").fit(X_train, y_train)
print("dummy   ", evaluate(dummy, X_val, y_val))
# dummy    {'auc': 0.5, 'ap': 0.031}   <- average precision equals the base rate

logit = LogisticRegression(max_iter=1000).fit(X_train_numeric, y_train)
print("logistic", evaluate(logit, X_val_numeric, y_val))

from sklearn.compose import ColumnTransformer
from sklearn.pipeline import Pipeline
from sklearn.impute import SimpleImputer
from sklearn.preprocessing import OneHotEncoder, StandardScaler
from sklearn.ensemble import HistGradientBoostingClassifier

numeric = ["amount", "term_months", "income", "days_since_last_application"]
categorical = ["product", "channel", "region"]

pre = ColumnTransformer([
    ("num", Pipeline([
        # The median is learned from train only, because the Pipeline is
        # fitted on train only. That is the whole reason for this structure.
        ("impute", SimpleImputer(strategy="median")),
        ("scale", StandardScaler()),
    ]), numeric),
    ("cat", Pipeline([
        ("impute", SimpleImputer(strategy="constant", fill_value="MISSING")),
        # Categories unseen at training time must not raise in production.
        ("onehot", OneHotEncoder(handle_unknown="infrequent_if_exist",
                                 min_frequency=20)),
    ]), categorical),
])

pipe = Pipeline([("pre", pre),
                 ("model", HistGradientBoostingClassifier(random_state=0))])

pipe.fit(X_train, y_train)
print(evaluate(pipe, X_val, y_val))

import numpy as np
from sklearn.metrics import precision_recall_curve

COST_FP = 12       # a manual review
COST_FN = 800      # the loss we did not prevent

p = pipe.predict_proba(X_val)[:, 1]

thresholds = np.linspace(0.01, 0.99, 99)
costs = []
for t in thresholds:
    pred = p >= t
    fp = int(((pred == 1) & (y_val == 0)).sum())
    fn = int(((pred == 0) & (y_val == 1)).sum())
    costs.append(fp * COST_FP + fn * COST_FN)

best = thresholds[int(np.argmin(costs))]
print("threshold", round(best, 3), "expected cost", min(costs))
print("cost at the default 0.5:", costs[49])

import itertools, json, random
from sklearn.base import clone

space = {
    "model__learning_rate": [0.01, 0.03, 0.1, 0.3],
    "model__max_leaf_nodes": [15, 31, 63, 127],
    "model__min_samples_leaf": [10, 20, 50, 100],
    "model__l2_regularization": [0.0, 0.1, 1.0],
}

rng = random.Random(0)
combos = list(itertools.product(*space.values()))
rng.shuffle(combos)

results = []
for combo in combos[:30]:
    params = dict(zip(space.keys(), combo))
    m = clone(pipe).set_params(**params).fit(X_train, y_train)
    row = {"params": params,
           "train": evaluate(m, X_train, y_train)["ap"],
           "val": evaluate(m, X_val, y_val)["ap"]}
    results.append(row)
    print(json.dumps(row))          # or log it to MLflow / Weights and Biases

best = max(results, key=lambda r: r["val"])
print("best", best)
# A large train-val gap here is overfitting: fewer leaves, more
# min_samples_leaf, more l2 - in that order.

import joblib, pandas as pd
from fastapi import FastAPI
from pydantic import BaseModel

MODEL_VERSION = "2024-11-03-hgb-r0"
pipe = joblib.load("model.joblib")
THRESHOLD = 0.082          # chosen in the metric stage, not invented here

class Application(BaseModel):
    amount: float
    term_months: int
    income: float | None = None      # nullable: the imputer handles it
    days_since_last_application: float | None = None
    product: str
    channel: str
    region: str

app = FastAPI()

@app.get("/health")
def health():
    return {"ok": True, "model_version": MODEL_VERSION}

@app.post("/predict")
def predict(a: Application):
    row = pd.DataFrame([a.model_dump()])
    p = float(pipe.predict_proba(row)[0, 1])
    # Return the probability as well as the decision: downstream teams will
    # want to move the threshold without redeploying the model.
    return {"probability": p, "decision": p >= THRESHOLD,
            "model_version": MODEL_VERSION}

import numpy as np

def psi(reference: np.ndarray, current: np.ndarray, bins: int = 10) -> float:
    "Population stability index. Below 0.1 is stable, above 0.25 is a real shift."
    edges = np.quantile(reference, np.linspace(0, 1, bins + 1))
    edges[0], edges[-1] = -np.inf, np.inf
    r = np.histogram(reference, edges)[0] / len(reference)
    c = np.histogram(current, edges)[0] / len(current)
    eps = 1e-6                       # empty bins otherwise give infinity
    return float(np.sum((c - r) * np.log((c + eps) / (r + eps))))

daily = predictions_log[predictions_log["date"] == yesterday]
for col in numeric:
    print(col, round(psi(X_train[col].to_numpy(), daily[col].to_numpy()), 3))

print("positive rate", round(daily["decision"].mean(), 4),
      "vs training", round(y_train.mean(), 4))
print("null rate", daily[numeric].isna().mean().round(3).to_dict())