from gridiron.data import pbp, nfl
PCOLS = ["season_type", "game_id", "play_id", "posteam", "defteam", "play_type", "down", "ydstogo", "yardline_100",
"yrdln", "pass", "qb_dropback", "xpass", "week", "qtr", "time", "fixed_drive", "fixed_drive_result",
"yards_gained", "epa", "success", "sack", "interception", "fumble_lost", "touchdown", "pass_location",
"run_location", "score_differential", "desc", "rusher_player_name", "receiver_player_name"]
FTN_COLS = ["game_id", "play_id", "starting_hash", "qb_location", "n_offense_backfield", "is_motion", "is_no_huddle",
"is_play_action", "is_screen_pass", "is_rpo", "n_pass_rushers", "n_blitzers"]
PART_COLS = ["game_id", "play_id", "offense_personnel", "defense_personnel", "defenders_in_box",
"defense_coverage_type", "defense_man_zone_type", "was_pressure", "time_to_throw"]
def load_season(season):
"""Runs and dropbacks, with FTN charting and participation attached."""
d = pbp([season], columns=PCOLS)
d = d[d.play_type.isin(["pass", "run"]) & d.down.notna()].copy()
d["play_id"] = d.play_id.astype(float)
f = nfl().load_ftn_charting([season]).to_pandas()
f = f.rename(columns={"nflverse_game_id": "game_id", "nflverse_play_id": "play_id"})
f["play_id"] = f.play_id.astype(float)
p = nfl().load_participation([season]).to_pandas().rename(columns={"nflverse_game_id": "game_id"})
p["play_id"] = p.play_id.astype(float)
return (d.merge(f[FTN_COLS], on=["game_id", "play_id"], how="left")
.merge(p[PART_COLS], on=["game_id", "play_id"], how="left"))
def count_pos(s, keys):
if not isinstance(s, str):
return np.nan
c = {k: int(v) for v, k in re.findall(r"(\d+) ([A-Z]+)", s)}
return sum(c.get(k, 0) for k in keys)
def pers_code(s):
"""FTN-era personnel string ('1 C, 2 G, 1 QB, 1 RB, 2 T, 1 TE, 3 WR') -> '11'."""
if not isinstance(s, str):
return "?"
c = {k: int(v) for v, k in re.findall(r"(\d+) ([A-Z]+)", s)}
return f"{c.get('RB', 0) + c.get('FB', 0)}{c.get('TE', 0)}"
COV = {"COVER_0": "C0", "COVER_1": "C1", "COVER_2": "C2", "2_MAN": "2M", "COVER_3": "C3", "COVER_4": "C4",
"COVER_6": "C6", "COVER_9": "C9", "COMBO": "CMB", "BLOWN": "BLN", "PREVENT": "PRV"}
TWO_HIGH = {"C2", "2M", "C4", "C6"}
ONE_HIGH = {"C0", "C1", "C3", "C9"} # C0 has nobody deep; "not two-high" for the shell call
def chart_codes(d):
"""The public charting, written in the chart sheet's short codes."""
d = d.copy()
d["pers"] = d.offense_personnel.map(pers_code)
nb = d.n_offense_backfield
d["align"] = d.qb_location.fillna("?") + nb.map(lambda v: "?" if pd.isna(v) else str(int(min(v, 2))))
d["mot"] = np.where(d.is_motion.fillna(False).astype(bool), "M", "–")
d["ndb"] = d.defense_personnel.map(lambda s: count_pos(s, ("CB", "FS", "SS", "S", "DB")))
d["dpers"] = d.ndb.map(lambda v: {4: "base", 5: "nickel", 6: "dime", 7: "quarter"}.get(int(v), "?")
if pd.notna(v) else "?")
d["box"] = d.defenders_in_box.map(lambda v: "?" if pd.isna(v) else str(int(v)))
d["cv"] = d.defense_coverage_type.map(COV)
d["rush"] = d.n_pass_rushers
kind = np.where(d.qb_dropback == 1, "P", "R")
kind = np.where(d.is_play_action.fillna(False).astype(bool) & (kind == "P"), "PA", kind)
kind = np.where(d.is_screen_pass.fillna(False).astype(bool), "SCR", kind)
d["kind"] = kind
d["is_pass"] = (d.qb_dropback == 1).astype(int)
def res(row):
if row.sack == 1:
return "sk FUM" if row.fumble_lost == 1 else f"sk {int(row.yards_gained):+d}".replace("-", "−")
if row.interception == 1:
return "INT"
if row.fumble_lost == 1:
return "FUM"
if row.touchdown == 1:
return "TD"
return f"{int(row.yards_gained):+d}".replace("-", "−")
d["res"] = d.apply(res, axis=1)
d["hash"] = d.starting_hash.where(d.starting_hash.isin(["L", "M", "R"]), "?") # FTN writes 0 when unknown
return d
S24 = load_season(2024)
S25 = load_season(2025)
# ------------------------------------------------------------------ Super Bowl LX, both offenses
SB_ID = "2025_22_SEA_NE"
SB = chart_codes(S25[S25.game_id == SB_ID]).sort_values("play_id").reset_index(drop=True)
NE = SB[SB.posteam == "NE"]
SEA = SB[SB.posteam == "SEA"]
NE_DB = NE[NE.qb_dropback == 1]
NE_LAB = NE_DB[NE_DB.cv.isin(TWO_HIGH | ONE_HIGH)]
SBX = dict(
ne_snaps=len(NE), ne_db=len(NE_DB), ne_pass=float(NE.is_pass.mean()),
four=int((NE_DB.rush == 4).sum()), five=int((NE_DB.rush >= 5).sum()),
sacks=int(NE.sack.sum()), sacks4=int(((NE.sack == 1) & (NE.rush == 4)).sum()),
press=float(NE_DB.was_pressure.fillna(False).astype(float).mean()),
two=float(NE_LAB.cv.isin(TWO_HIGH).mean()), n_lab=len(NE_LAB),
man=float((NE_DB.defense_man_zone_type == "MAN_COVERAGE").sum()
/ NE_DB.defense_man_zone_type.isin(["MAN_COVERAGE", "ZONE_COVERAGE"]).sum()),
sea_snaps=len(SEA), sea_run=float(1 - SEA.is_pass.mean()), sea_uc=float((SEA.qb_location == "U").mean()),
walker_n=int((SEA.rusher_player_name == "K.Walker").sum()),
walker_y=int(SEA.loc[SEA.rusher_player_name == "K.Walker", "yards_gained"].sum()),
ne_succ=float(NE.success.mean()), sea_succ=float(SEA.success.mean()),
)
_R25 = S25[S25.season_type == "REG"]
SBX["ne_reg_pass"] = float(_R25.loc[_R25.posteam == "NE", "qb_dropback"].mean())
SBX["sea_reg_run"] = float(1 - _R25.loc[_R25.posteam == "SEA", "qb_dropback"].mean())
# Seattle's defense against dropbacks in the 2025 regular season, and the league's
REG = S25[(S25.season_type == "REG") & (S25.qb_dropback == 1)].copy()
REG["c"] = REG.defense_coverage_type.map(COV)
REG_LAB = REG[REG.c.isin(TWO_HIGH | ONE_HIGH)]
SEA_REG = REG_LAB[REG_LAB.defteam == "SEA"]
SEASON = dict(sea_two=float(SEA_REG.c.isin(TWO_HIGH).mean()), lg_two=float(REG_LAB.c.isin(TWO_HIGH).mean()),
sea_four=float((REG[(REG.defteam == "SEA") & (REG.n_pass_rushers > 0)].n_pass_rushers == 4).mean()),
lg_four=float((REG[REG.n_pass_rushers > 0].n_pass_rushers == 4).mean()),
sea_press=float(REG.loc[REG.defteam == "SEA", "was_pressure"].fillna(False).astype(float).mean()),
lg_press=float(REG.was_pressure.fillna(False).astype(float).mean()))
def sea_two_rate(down):
"""Seattle's 2025 regular-season two-high rate on dropbacks, by down (1st, 2nd, 3rd/4th)."""
dd = SEA_REG.down.clip(upper=3)
return float(SEA_REG.loc[dd == min(down, 3), "c"].isin(TWO_HIGH).mean())
# The strip-sack: third-and-6, end of the third quarter
SS_ROW = NE[(NE.qtr == 3) & (NE.sack == 1) & (NE.fumble_lost == 1)].iloc[0]
# ------------------------------------------------------------------ the pre-snap model (situation + formation)
def feats(d):
xp = np.clip(d.xpass.astype(float), 1e-3, 1 - 1e-3)
nb = d.n_offense_backfield.fillna(1)
return np.column_stack([np.ones(len(d)), np.log(xp / (1 - xp)), d.qb_location == "S", d.qb_location == "P",
nb == 0, nb >= 2]).astype(float)
def fit(X, y, lam=1.0, iters=30):
"""Ridge-penalised logistic regression by Newton's method (numpy only)."""
w = np.zeros(X.shape[1])
for _ in range(iters):
p = 1 / (1 + np.exp(-X @ w))
H = X.T @ (X * (p * (1 - p))[:, None]) + lam * np.eye(X.shape[1])
w -= np.linalg.solve(H, X.T @ (p - y) + lam * w)
return w
TR = S24[(S24.season_type == "REG") & S24.xpass.notna()]
W_INFO = fit(feats(TR), TR["pass"].values.astype(float))
def p_info(d):
return 1 / (1 + np.exp(-feats(d) @ W_INFO))
def conf_of(p, cap=0.9):
return float(np.clip(np.round(max(p, 1 - p) * 20) / 20, 0.5, cap))
SB["p_rule"] = p_info(SB)
SB["call"] = np.where(SB.p_rule >= 0.5, "P", "R")
SB["conf"] = SB.p_rule.map(conf_of)
# ------------------------------------------------------------------ a simulated reader: every Panthers snap of 2025
CAR = S25[(S25.posteam == "CAR") & S25.xpass.notna()].sort_values(["week", "play_id"]).reset_index(drop=True)
CAR["p_info"] = p_info(CAR)
GAMES = list(dict.fromkeys(CAR.game_id))
CAR["g"] = CAR.game_id.map({g: i for i, g in enumerate(GAMES)})
N_GAMES = len(GAMES)
LG = np.log(CAR.p_info / (1 - CAR.p_info)).values
Y = CAR["pass"].values.astype(int)
def simulate(k0, k1, s0, s1, cap0, cap1, seed):
"""A reader who sees what the model sees, through noise, and over- or under-states it.
k > 1 overstates, k < 1 understates; both settle over the season from (k0, s0, cap0) to (k1, s1, cap1)."""
rng = np.random.default_rng(seed)
frac = (CAR.g / (N_GAMES - 1)).values
k, sd = k0 + (k1 - k0) * frac, s0 + (s1 - s0) * frac
cap = np.where(frac < 0.45, cap0, cap1)
pr = 1 / (1 + np.exp(-(k * LG + rng.normal(0, 1, len(CAR)) * sd)))
conf = np.clip(np.round(np.maximum(pr, 1 - pr) * 20) / 20, 0.5, cap)
return np.where(pr >= 0.5, "P", "R"), conf
ME_CALL, ME_CONF = simulate(1.5, 0.95, 0.6, 0.4, 0.95, 0.90, seed=15)
A_CALL, A_CONF = simulate(2.2, 2.2, 0.3, 0.3, 0.99, 0.99, seed=7)
B_CALL, B_CONF = simulate(0.35, 0.35, 0.1, 0.1, 0.85, 0.85, seed=7)
LOGDF = pd.DataFrame({
"game": CAR.g + 1, "week": CAR.week.astype(int), "opp": CAR.defteam, "down": CAR.down.astype(int),
"togo": CAR.ydstogo.astype(int), "align": CAR.qb_location.fillna("?") + CAR.n_offense_backfield.map(
lambda v: "?" if pd.isna(v) else str(int(min(v, 2)))),
"call": ME_CALL, "conf": ME_CONF,
"clue": np.where(np.abs(CAR.p_info - CAR.xpass) > 0.10, "align", "sit"),
"result": np.where(Y == 1, "P", "R")})
HIT = (LOGDF.call == LOGDF.result).values.astype(float)
def brier_of(conf, hit):
return float(np.mean((np.asarray(conf) - np.asarray(hit)) ** 2))
BASE_P = float(Y.mean())
REF = dict(half=0.25, base=BASE_P * (1 - BASE_P), xpass=float(np.mean((CAR.xpass.values - Y) ** 2)),
info=float(np.mean((CAR.p_info.values - Y) ** 2)), me=brier_of(ME_CONF, HIT))
FIRST = (CAR.g < 8).values
ME = dict(n=len(CAR), acc=float(HIT.mean()), conf=float(ME_CONF.mean()),
b1=brier_of(ME_CONF[FIRST], HIT[FIRST]), b2=brier_of(ME_CONF[~FIRST], HIT[~FIRST]),
a1=float(HIT[FIRST].mean()), a2=float(HIT[~FIRST].mean()),
c1=float(ME_CONF[FIRST].mean()), c2=float(ME_CONF[~FIRST].mean()),
xpass_acc=float(((CAR.xpass >= 0.5).astype(int) == Y).mean()),
info_acc=float(((CAR.p_info >= 0.5).astype(int) == Y).mean()))
DD_RUN = ((CAR.ydstogo <= 3) & (CAR.down != 3)) | ((CAR.down == 1) & (CAR.ydstogo <= 6))
ME["dd_acc"] = float((np.where(DD_RUN, 0, 1) == Y).mean())
_AL = (LOGDF.clue == "align").values
ME["align_me"] = float(HIT[_AL].mean())
ME["align_x"] = float(((CAR.xpass.values[_AL] >= 0.5).astype(int) == Y[_AL]).mean())
ME["always"] = BASE_P
def top_band(conf, hit, lo):
m = np.asarray(conf) >= lo
return float(np.asarray(hit)[m].mean()), int(m.sum())
ME["top1"], ME["n_top1"] = top_band(ME_CONF[FIRST], HIT[FIRST], 0.9)
ME["top2"], ME["n_top2"] = top_band(ME_CONF[~FIRST], HIT[~FIRST], 0.9)
ME["mid1"], ME["n_mid1"] = (lambda m: (float(HIT[FIRST][m].mean()), int(m.sum())))(
(ME_CONF[FIRST] >= 0.7) & (ME_CONF[FIRST] < 0.8))
# Readers A and B for the drill
HIT_A = ((A_CALL == "P") == (Y == 1)).astype(float)
HIT_B = ((B_CALL == "P") == (Y == 1)).astype(float)
def recalibrated(conf, hit):
"""Brier score if every stated value were replaced by its own hit rate (perfect calibration, same calls)."""
t = pd.DataFrame({"c": conf, "h": hit})
return brier_of(t.groupby("c").h.transform("mean"), t.h)
AB = dict(a_acc=float(HIT_A.mean()), b_acc=float(HIT_B.mean()), a_b=brier_of(A_CONF, HIT_A),
b_b=brier_of(B_CONF, HIT_B), a_conf=float(A_CONF.mean()), b_conf=float(B_CONF.mean()),
a_fix=recalibrated(A_CONF, HIT_A), b_fix=recalibrated(B_CONF, HIT_B),
a_90=top_band(A_CONF, HIT_A, 0.9), b_top=float(B_CONF.max()))
def murphy(f, o):
"""Brier = reliability − resolution + uncertainty (Murphy 1973), grouping by the distinct forecast values."""
t = pd.DataFrame({"f": f, "o": o})
g = t.groupby("f").o.agg(["size", "mean"])
ob = t.o.mean()
rel = float((g["size"] * (g.index - g["mean"]) ** 2).sum() / len(t))
res = float((g["size"] * (g["mean"] - ob) ** 2).sum() / len(t))
return rel, res, float(ob * (1 - ob))
F_PASS = np.where(ME_CALL == "P", ME_CONF, 1 - ME_CONF) # the same calls, as a chance of a pass
MUR1 = murphy(F_PASS[FIRST], Y[FIRST])
MUR2 = murphy(F_PASS[~FIRST], Y[~FIRST])
# the season in four stretches (the same split the scoring code prints)
_PART = pd.cut(LOGDF.game, [0, 4, 8, 13, 18], labels=False)
QTR = [brier_of(ME_CONF[_PART == i], HIT[_PART == i]) for i in range(4)]
QTR_ACC = [float(HIT[_PART == i].mean()) for i in range(4)]
QTR_CONF = [float(ME_CONF[_PART == i].mean()) for i in range(4)]
def group_cost(c, h):
"""Expected Brier score per call for a group stated at c that came true a share h of the time."""
return h * (1 - c) ** 2 + (1 - h) * c ** 2
_A90 = A_CONF >= 0.9
AB["a90_c"], AB["a90_share"] = float(A_CONF[_A90].mean()), float(_A90.mean())
AB["a90_cost"] = group_cost(AB["a90_c"], AB["a_90"][0])
AB["a90_fixcost"] = group_cost(AB["a_90"][0], AB["a_90"][0])
_B_LO = B_CONF <= 0.65
AB["b_lo_c"], AB["b_lo_h"], AB["b_lo_share"] = (float(B_CONF[_B_LO].mean()), float(HIT_B[_B_LO].mean()),
float(_B_LO.mean()))
AB["b_lo_cost"] = group_cost(AB["b_lo_c"], AB["b_lo_h"])
# Seattle's 2025 two-high rate by down, as printed on the scouting card in the second-layer code below
CARD = {d_: round(sea_two_rate(d_), 2) for d_ in (1, 2, 3)}
