Credit-risk tutorial¶

The full workflow: mine candidate constraints from data, review and accept them, price a missing-value flip, build a rating ladder, and generate diverse alternatives.

In [1]:

import numpy as np

rng = np.random.default_rng(42)
n = 6000
names = [
    "income_monthly", "utilization", "n_active_loans", "n_loans_total",
    "max_dpd_30d", "max_dpd_12m", "months_since_last_delinq", "age",
]
income = rng.lognormal(8.3, 0.5, n).round(-1)
utilization = rng.beta(2, 3, n).round(3)
n_total = rng.poisson(4, n).astype(float) + 1
n_active = np.minimum(np.floor(n_total * rng.beta(3, 2, n)), n_total)
dpd_12m = np.floor(rng.exponential(6, n)) * (rng.random(n) < 0.4)
dpd_30d = np.floor(dpd_12m * rng.beta(2, 4, n))
months_delinq = rng.exponential(14, n).round(0)
months_delinq[dpd_12m == 0] = np.nan          # no delinquency -> no record
age = rng.integers(21, 75, n).astype(float)

X = np.column_stack([income, utilization, n_active, n_total,
                     dpd_30d, dpd_12m, months_delinq, age])
risk = (
    -0.9 * np.log(income / 4000)
    + 2.2 * utilization
    + 0.35 * dpd_30d + 0.15 * dpd_12m
    - 0.02 * np.nan_to_num(months_delinq, nan=36.0)
    - 0.015 * age
)
y = (risk + rng.logistic(scale=0.8, size=n) > np.median(risk)).astype(int)
X.shape, y.mean().round(3)

import numpy as np rng = np.random.default_rng(42) n = 6000 names = [ "income_monthly", "utilization", "n_active_loans", "n_loans_total", "max_dpd_30d", "max_dpd_12m", "months_since_last_delinq", "age", ] income = rng.lognormal(8.3, 0.5, n).round(-1) utilization = rng.beta(2, 3, n).round(3) n_total = rng.poisson(4, n).astype(float) + 1 n_active = np.minimum(np.floor(n_total * rng.beta(3, 2, n)), n_total) dpd_12m = np.floor(rng.exponential(6, n)) * (rng.random(n) < 0.4) dpd_30d = np.floor(dpd_12m * rng.beta(2, 4, n)) months_delinq = rng.exponential(14, n).round(0) months_delinq[dpd_12m == 0] = np.nan # no delinquency -> no record age = rng.integers(21, 75, n).astype(float) X = np.column_stack([income, utilization, n_active, n_total, dpd_30d, dpd_12m, months_delinq, age]) risk = ( -0.9 * np.log(income / 4000) + 2.2 * utilization + 0.35 * dpd_30d + 0.15 * dpd_12m - 0.02 * np.nan_to_num(months_delinq, nan=36.0) - 0.015 * age ) y = (risk + rng.logistic(scale=0.8, size=n) > np.median(risk)).astype(int) X.shape, y.mean().round(3)

Out[1]:

((6000, 8), np.float64(0.534))

In [2]:

import xgboost as xgb

model = xgb.XGBClassifier(n_estimators=120, max_depth=4, learning_rate=0.2,
                          random_state=0)
model.fit(X, y)
model.get_booster().feature_names = list(names)   # domain names for explanations
proba = model.predict_proba(X)[:, 1]
cutoff = 0.30                     # the credit policy's PD cutoff
applicant = X[int(np.argmax(proba))]      # a clearly declined applicant
float(model.predict_proba(applicant.reshape(1, -1))[0, 1])

import xgboost as xgb model = xgb.XGBClassifier(n_estimators=120, max_depth=4, learning_rate=0.2, random_state=0) model.fit(X, y) model.get_booster().feature_names = list(names) # domain names for explanations proba = model.predict_proba(X)[:, 1] cutoff = 0.30 # the credit policy's PD cutoff applicant = X[int(np.argmax(proba))] # a clearly declined applicant float(model.predict_proba(applicant.reshape(1, -1))[0, 1])

Out[2]:

0.9998914003372192

1. Mine candidate constraints¶

Invariants are mined from the sample and returned for review — nothing is auto-applied. Note the DPD hierarchy and loan-count order arriving as ready-to-paste code.

In [3]:

import treecf

mined = treecf.suggest_constraints(X, feature_names=names)
for s in mined[:8]:
    print(s.as_code())

import treecf mined = treecf.suggest_constraints(X, feature_names=names) for s in mined[:8]: print(s.as_code())

constraint("n_active_loans <= n_loans_total")  # support=1.0000, n=6000
constraint("max_dpd_30d <= max_dpd_12m")  # support=1.0000, n=6000
constraint("max_dpd_12m <= income_monthly")  # support=1.0000, n=6000
constraint("months_since_last_delinq <= income_monthly")  # support=1.0000, n=2044
constraint("age <= income_monthly")  # support=1.0000, n=6000
constraint("utilization <= n_loans_total")  # support=1.0000, n=6000
constraint("n_loans_total <= age")  # support=1.0000, n=6000
constraint("max_dpd_30d <= age")  # support=1.0000, n=6000

In [4]:

accepted = [s.constraint for s in mined
            if s.kind == "order" and s.support == 1.0]
accepted

accepted = [s.constraint for s in mined if s.kind == "order" and s.support == 1.0] accepted

Out[4]:

[Linear(coefficients={'n_active_loans': 1.0, 'n_loans_total': -1.0}, op='<=', rhs=0.0, missing_policy='satisfied'),
 Linear(coefficients={'max_dpd_30d': 1.0, 'max_dpd_12m': -1.0}, op='<=', rhs=0.0, missing_policy='satisfied'),
 Linear(coefficients={'max_dpd_12m': 1.0, 'income_monthly': -1.0}, op='<=', rhs=0.0, missing_policy='satisfied'),
 Linear(coefficients={'months_since_last_delinq': 1.0, 'income_monthly': -1.0}, op='<=', rhs=0.0, missing_policy='satisfied'),
 Linear(coefficients={'age': 1.0, 'income_monthly': -1.0}, op='<=', rhs=0.0, missing_policy='satisfied'),
 Linear(coefficients={'utilization': 1.0, 'n_loans_total': -1.0}, op='<=', rhs=0.0, missing_policy='satisfied'),
 Linear(coefficients={'n_loans_total': 1.0, 'age': -1.0}, op='<=', rhs=0.0, missing_policy='satisfied'),
 Linear(coefficients={'max_dpd_30d': 1.0, 'age': -1.0}, op='<=', rhs=0.0, missing_policy='satisfied')]

2. Explain with a priced missing-value flip¶

months_since_last_delinq is NaN when there is no delinquency record — and reaching that state is a legitimate recommendation with an explicit price (delta_miss).

In [5]:

from treecf import AllowMissing, Explainer, Freeze, Target

exp = Explainer(
    model,
    background=X,
    constraints=accepted + [
        Freeze("age"),
        AllowMissing("months_since_last_delinq", delta_miss=2.0),
    ],
    value_policy={"max_dpd_30d": "integer", "max_dpd_12m": "integer",
                  "n_active_loans": "integer", "n_loans_total": "integer"},
)
res = exp.explain(applicant, target=Target.probability(range=(0.0, cutoff)),
                  seed=0)
res.changes

from treecf import AllowMissing, Explainer, Freeze, Target exp = Explainer( model, background=X, constraints=accepted + [ Freeze("age"), AllowMissing("months_since_last_delinq", delta_miss=2.0), ], value_policy={"max_dpd_30d": "integer", "max_dpd_12m": "integer", "n_active_loans": "integer", "n_loans_total": "integer"}, ) res = exp.explain(applicant, target=Target.probability(range=(0.0, cutoff)), seed=0) res.changes

Out[5]:

{'income_monthly': (1640.0, 1902.739347036747),
 'utilization': (0.186, 0.1419999897480011),
 'n_loans_total': (2.0, 1.0),
 'max_dpd_30d': (15.0, 3.0),
 'max_dpd_12m': (24.0, 3.0),
 'months_since_last_delinq': (26.0, nan)}

3. The rating ladder¶

One search per band: the increasing cost of each better grade.

In [6]:

ladder = exp.explain(applicant, seed=0, target=Target.bands({
    "approve": (0.00, 0.30),
    "prime":   (0.00, 0.15),
    "super":   (0.00, 0.05),
}))
{k: (round(v.distance, 3) if hasattr(v, 'distance') else v.reason)
 for k, v in ladder.items()}

ladder = exp.explain(applicant, seed=0, target=Target.bands({ "approve": (0.00, 0.30), "prime": (0.00, 0.15), "super": (0.00, 0.05), })) {k: (round(v.distance, 3) if hasattr(v, 'distance') else v.reason) for k, v in ladder.items()}

Out[6]:

{'approve': 12.805, 'prime': 12.822, 'super': 15.293}

In [7]:

from treecf.viz import plot_ladder

plot_ladder(ladder);

from treecf.viz import plot_ladder plot_ladder(ladder);

4. Alternative plans, and which levers are essential¶

Block each of the primary plan's biggest levers in turn and re-solve. Levers with a workaround yield an alternative plan (at a higher cost); levers with none are essential — approval is unreachable without them. Both answers are useful to an applicant.

In [8]:

from treecf import Infeasible

base = accepted + [Freeze("age"),
                   AllowMissing("months_since_last_delinq", delta_miss=2.0)]
policy = {"max_dpd_30d": "integer", "max_dpd_12m": "integer",
          "n_active_loans": "integer", "n_loans_total": "integer"}
idx = {f: i for i, f in enumerate(names)}
levers = sorted(res.changes,
                key=lambda f: abs(res.changes[f][1] - res.changes[f][0])
                / exp.sigma[idx[f]], reverse=True)[:3]
alternatives, essential = [res], []
for lever in levers:
    alt = Explainer(model, background=X, value_policy=policy,
                    constraints=base + [Freeze(lever)])
    cand = alt.explain(applicant, seed=0, time_budget_s=30,
                       target=Target.probability(range=(0.0, cutoff)))
    if isinstance(cand, Infeasible):
        essential.append(lever)      # no plan exists without this lever
    else:
        alternatives.append(cand)
print("essential levers:", essential)
[(round(a.distance, 2), sorted(a.changes)) for a in alternatives]

from treecf import Infeasible base = accepted + [Freeze("age"), AllowMissing("months_since_last_delinq", delta_miss=2.0)] policy = {"max_dpd_30d": "integer", "max_dpd_12m": "integer", "n_active_loans": "integer", "n_loans_total": "integer"} idx = {f: i for i, f in enumerate(names)} levers = sorted(res.changes, key=lambda f: abs(res.changes[f][1] - res.changes[f][0]) / exp.sigma[idx[f]], reverse=True)[:3] alternatives, essential = [res], [] for lever in levers: alt = Explainer(model, background=X, value_policy=policy, constraints=base + [Freeze(lever)]) cand = alt.explain(applicant, seed=0, time_budget_s=30, target=Target.probability(range=(0.0, cutoff))) if isinstance(cand, Infeasible): essential.append(lever) # no plan exists without this lever else: alternatives.append(cand) print("essential levers:", essential) [(round(a.distance, 2), sorted(a.changes)) for a in alternatives]

essential levers: ['max_dpd_12m', 'max_dpd_30d']

Out[8]:

[(12.8,
  ['income_monthly',
   'max_dpd_12m',
   'max_dpd_30d',
   'months_since_last_delinq',
   'n_loans_total',
   'utilization']),
 (13.21,
  ['income_monthly',
   'max_dpd_12m',
   'max_dpd_30d',
   'months_since_last_delinq',
   'n_active_loans',
   'utilization'])]

In [9]:

from treecf.viz import plot_counterfactuals

plot_counterfactuals(alternatives);

from treecf.viz import plot_counterfactuals plot_counterfactuals(alternatives);

Compare the plans directly: every alternative's changes on shared feature axes (one color per plan, gray dots mark the factual values), and the cost-vs-outcome trade-off — what each plan asks for and what it buys.

In [10]:

from treecf.viz import plot_alternatives, plot_tradeoff

plot_alternatives(alternatives, explainer=exp);   # deltas in sigma units

from treecf.viz import plot_alternatives, plot_tradeoff plot_alternatives(alternatives, explainer=exp); # deltas in sigma units

In [11]:

plot_tradeoff(alternatives, target=Target.probability(range=(0.0, cutoff)));

plot_tradeoff(alternatives, target=Target.probability(range=(0.0, cutoff)));

5. Mass-producing counterfactuals for a day's declines¶

Score a day's applications, take the declines, and produce (up to) two recourse plans per applicant in one call — then store the batch and look plans up later without recomputing.

In [12]:

import time

declined = np.flatnonzero(proba > cutoff)[:200]     # today's declines
app_ids = [f"APP-{i:05d}" for i in declined]

start = time.perf_counter()
batch = exp.explain_batch(
    X[declined],
    target=Target.probability(range=(0.0, cutoff)),
    n_per_example=2,           # counterfactuals per applicant
    diversity="seeds",
    ids=app_ids,
    seed=0,
)
wall = time.perf_counter() - start
feasible = sum(r.feasible for r in batch)
print(f"{len(batch)} plans for {len(declined)} applicants "
      f"in {wall:.1f}s ({1000 * wall / len(declined):.0f} ms/applicant)")
print(f"feasible plans: {feasible}")

import time declined = np.flatnonzero(proba > cutoff)[:200] # today's declines app_ids = [f"APP-{i:05d}" for i in declined] start = time.perf_counter() batch = exp.explain_batch( X[declined], target=Target.probability(range=(0.0, cutoff)), n_per_example=2, # counterfactuals per applicant diversity="seeds", ids=app_ids, seed=0, ) wall = time.perf_counter() - start feasible = sum(r.feasible for r in batch) print(f"{len(batch)} plans for {len(declined)} applicants " f"in {wall:.1f}s ({1000 * wall / len(declined):.0f} ms/applicant)") print(f"feasible plans: {feasible}")

269 plans for 200 applicants in 3.6s (18 ms/applicant)
feasible plans: 269

In [13]:

import pathlib, tempfile

store_path = pathlib.Path(tempfile.mkdtemp()) / "counterfactuals_today.json"
batch.save(store_path)                        # store once...

from treecf import BatchResult
stored = BatchResult.load(store_path)
stored.for_id(app_ids[0])                     # ...look up any time

import pathlib, tempfile store_path = pathlib.Path(tempfile.mkdtemp()) / "counterfactuals_today.json" batch.save(store_path) # store once... from treecf import BatchResult stored = BatchResult.load(store_path) stored.for_id(app_ids[0]) # ...look up any time

Out[13]:

[BatchRecord(id='APP-00000', k=0, feasible=True, x_cf=array([4.69e+03, 3.13e-01, 4.00e+00, 8.00e+00, 0.00e+00, 0.00e+00,
        5.40e+01, 2.90e+01]), changes={'months_since_last_delinq': (nan, 54.0)}, distance=0.2857142857142857, n_changed=1, score_raw=-0.9703074229268243, score_prob=0.27481923054079793, seed=0, blocked_lever=None, coalition=None)]

In [14]:

stored.to_frame().head(6)      # or analyze the whole day as a DataFrame

stored.to_frame().head(6) # or analyze the whole day as a DataFrame

Out[14]:

	id	k	feasible	distance	n_changed	score_raw	score_prob	seed	blocked_lever	coalition	changed_features	cf_income_monthly	cf_utilization	cf_n_active_loans	cf_n_loans_total	cf_max_dpd_30d	cf_max_dpd_12m	cf_months_since_last_delinq	cf_age
0	APP-00000	0	True	0.285714	1	-0.970307	0.274819	0	None	None	[months_since_last_delinq]	4690.0	0.313	4.0	8.0	0.0	0.0	54.0	29.0
1	APP-00001	0	True	0.285714	1	-1.356183	0.204861	1009	None	None	[months_since_last_delinq]	2390.0	0.222	8.0	9.0	0.0	0.0	54.0	36.0
2	APP-00002	0	True	0.285714	1	-1.235083	0.225293	2018	None	None	[months_since_last_delinq]	5860.0	0.072	8.0	10.0	0.0	0.0	61.0	39.0
3	APP-00002	1	True	1.174497	2	-1.052479	0.258749	2019	None	None	[n_active_loans, utilization]	5860.0	0.046	7.0	10.0	0.0	0.0	NaN	39.0
4	APP-00003	0	True	1.000000	1	-1.320359	0.210759	3027	None	None	[n_active_loans]	6440.0	0.142	7.0	12.0	0.0	0.0	NaN	54.0
5	APP-00004	0	True	0.547253	2	-1.202409	0.231047	4036	None	None	[income_monthly, months_since_last_delinq]	1860.0	0.579	3.0	4.0	0.0	0.0	61.0	28.0

Visualizing the batch¶

Four batch-level views: which levers the plans use (and in which direction), the effort each change costs per plan, cost/sparsity/feasibility at a glance, and how far each lever actually moves across applicants.

In [15]:

from treecf.viz_batch import plot_batch_levers, plot_batch_matrix

plot_batch_levers(batch);

from treecf.viz_batch import plot_batch_levers, plot_batch_matrix plot_batch_levers(batch);

In [16]:

plot_batch_matrix(batch, explainer=exp, max_row_labels=15);

plot_batch_matrix(batch, explainer=exp, max_row_labels=15);

In [17]:

from treecf.viz_batch import plot_batch_deltas, plot_batch_summary

plot_batch_summary(batch);

from treecf.viz_batch import plot_batch_deltas, plot_batch_summary plot_batch_summary(batch);

In [18]:

plot_batch_deltas(batch, explainer=exp);   # deltas in sigma units

plot_batch_deltas(batch, explainer=exp); # deltas in sigma units

6. Coalitions: grouped recourse¶

A single plan can mix unrelated levers — income, credit usage, and debt history in one instruction. Coalitions split recourse by what the applicant controls together: one counterfactual per named feature group, each allowed to change only its own group. An infeasible group is a finding in itself: that front alone cannot reach the target. The "(all levers)" baseline shows what the grouping costs versus the unrestricted optimum. Opt-in — plain explain is unchanged.

In [19]:

groups = {
    "debt history": ["max_dpd_30d", "max_dpd_12m", "months_since_last_delinq"],
    "credit usage": ["utilization", "n_active_loans", "n_loans_total"],
    "income":       ["income_monthly"],
}
grouped = exp.explain_coalitions(
    applicant, target=Target.probability(range=(0.0, cutoff)),
    coalitions=groups, include_full=True, seed=0,
)
{name: (round(out.distance, 2) if hasattr(out, 'distance') else 'infeasible')
 for name, out in grouped.items()}

groups = { "debt history": ["max_dpd_30d", "max_dpd_12m", "months_since_last_delinq"], "credit usage": ["utilization", "n_active_loans", "n_loans_total"], "income": ["income_monthly"], } grouped = exp.explain_coalitions( applicant, target=Target.probability(range=(0.0, cutoff)), coalitions=groups, include_full=True, seed=0, ) {name: (round(out.distance, 2) if hasattr(out, 'distance') else 'infeasible') for name, out in grouped.items()}

Out[19]:

{'(all levers)': 12.8,
 'debt history': 'infeasible',
 'credit usage': 'infeasible',
 'income': 'infeasible'}

In [20]:

plot_alternatives(grouped, explainer=exp);   # coalition names label the plans

plot_alternatives(grouped, explainer=exp); # coalition names label the plans

In [21]:

plot_tradeoff(grouped, target=Target.probability(range=(0.0, cutoff)));

plot_tradeoff(grouped, target=Target.probability(range=(0.0, cutoff)));

The same mode scales to the whole batch: one record per coalition per applicant, with the group name in the coalition column.

In [22]:

grouped_batch = exp.explain_batch(
    X[declined][:20], target=Target.probability(range=(0.0, cutoff)),
    diversity="coalitions", coalitions=groups, include_full=True,
    ids=app_ids[:20], seed=0,
)
grouped_batch.to_frame()[
    ["id", "k", "coalition", "feasible", "distance", "n_changed"]
].head(8)

grouped_batch = exp.explain_batch( X[declined][:20], target=Target.probability(range=(0.0, cutoff)), diversity="coalitions", coalitions=groups, include_full=True, ids=app_ids[:20], seed=0, ) grouped_batch.to_frame()[ ["id", "k", "coalition", "feasible", "distance", "n_changed"] ].head(8)

Out[22]:

	id	k	coalition	feasible	distance	n_changed
0	APP-00000	0	(all levers)	True	0.285714	1.0
1	APP-00000	1	debt history	True	0.285714	1.0
2	APP-00000	2	credit usage	True	0.758389	1.0
3	APP-00000	3	income	True	1.038462	1.0
4	APP-00001	0	(all levers)	True	0.285714	1.0
5	APP-00001	1	debt history	True	0.285714	1.0
6	APP-00001	2	credit usage	True	2.046980	2.0
7	APP-00001	3	income	True	2.230769	1.0