The puck is cylindrical

There is a German phrase about soccer: Der Ball ist rund translating into "the ball is round". Anything can happen, the game resists certainty, and it is always a probability function without absolutes. Hockey takes that same premise with an amplifier. The puck is cylindrical. It bounces, slides, flips, stands on edges, disappears under skates, and spends half the night being chased by players changing one at a time.

That makes hockey analytics a credit assignment problem with a long list of extenuating circumstances. The sport is shift-oriented and deployment-oriented. Goalies contort and flail their impossibly nimble limbs. Power plays and empty nets create separate little micro games within the game. A goal is scored by one player, assisted by at most two, and enabled by a dozen invisible decisions that never touch a box score: the defenseman who killed a rush thirty seconds earlier, the smooth skater who did the clean zone entry upon the prior shift leading to an offensive zone draw, the winger whose forecheck forced the turnover, the center who lost the faceoff but won the board battle after.

Plus-minus, a naive approach for on-ice credit, is so confounded by teammates, opponents, usage, and goaltending that it mostly measures whose line you're on. The modern answer is RAPM, Regularized Adjusted Plus-Minus, which stops asking "who scored" and starts asking "what changes when you're on the ice, after accounting for everyone else who is also on the ice." Answering that question properly requires reconstructing every second of a season, which is a data engineering problem before it's ever a statistics problem.

So let's try to make sense of the madness. nhl_2526 is a dataset that pulls the entire 2025-26 NHL season from the league's public API into oleander, builds an expected-goals model, chops 1,312 games into 518,282 stints, and solves one very large, very sparse ridge regression. This post walks through the whole thing, including the actual math on one actual Celebrini shift.

Pulling down a season

Everything starts at the NHL's public API (api-web.nhle.com, no auth required, please be polite). Most of the ingestion layer was prompt-driven from the NHL API documentation: oleander Python was written for each endpoint to conform to the API's actual response shape, then checked against representative responses and iterated until the data landed cleanly. The result is a series of inspectable raw tables - raw_players, raw_schedule, raw_skater_game_logs, raw_goalie_game_logs, the NHL Edge tables, raw_pbp_shots, and raw_shifts - rather than one prematurely modeled dataset. A single entrypoint.py orchestrates the PySpark jobs that populate those tables in the oleander.nhl_2526 namespace; the joins, feature engineering, and models begin downstream of that raw boundary.

The ingestion pattern is deliberately boring. The driver fetches JSON with a 0.3-0.5s delay between requests and exponential backoff on HTTP 429, then Spark takes over for everything downstream:

# nhl/ingest_pbp.py (abridged)
for game_id in game_ids:
    payload = fetch_with_backoff(
        f"https://api-web.nhle.com/v1/gamecenter/{game_id}/play-by-play"
    )
    rows.extend(parse_shot_events(payload))

df = spark.createDataFrame(rows, schema=PBP_SHOT_SCHEMA)
df.writeTo("oleander.nhl_2526.raw_pbp_shots").createOrReplace()

Deployment to oleander's managed Spark is a Makefile and two commands. make zips the nhl module and packs a virtualenv, then:

# the agent ran this, not me
oleander spark jobs upload entrypoint.py \
  --py-files out/pyfiles.zip., \
  --virtualenv out/environment.tar.gz

oleander spark jobs submit entrypoint.py \
  --namespace nhl_2526 --name nhl-full-pipeline

Agents and the oleander MCP built and ran this entire workflow looping through a couple of errors until it happily validated full game and roster collections. Because every job runs through oleander, each table write emits lineage events automatically, so when a downstream WAR number looks suspicious you can walk the lineage graph from skater_percentiles back through rapm_skater_results to the exact ingest_shifts run that fed it, with logs and traces attached.

The full pipeline is five phases, and order matters, because everything compounds:

Phase 1: Ingest          rosters, schedule, game logs, Edge, PBP, shifts
Phase 2: xG + stints     build_xg (expected goals), build_stint_events, game context
Phase 3: RAPM            ridge regression on 5v5 stints
Phase 4: Player models   GAx, (goals above expected) GSAx, (goals saved above expected, for goalies) WAR components (wins above replacement)
Phase 5: Rankings        percentiles within F / D / G

Stints: the atomic unit of on-ice credit

A shift is one player's trip over the boards. A stint is the intersection of everybody's shifts: a maximal slice of game time during which the ten skaters and two goalies on the ice do not change. Any substitution, even one winger swapping on the fly, ends the current stint and starts a new one. The ingest_shifts job takes 986,183 raw shift records, turns every shift start and end into a boundary, and slices each game along all of them.

The 2025-26 season decomposes into 518,282 stints, of which 445,916 are 5v5, totaling 1,064 hours of even-strength hockey. Their durations are wildly skewed:

The median 5v5 stint is 3 seconds long. That sounds broken until you watch a line change: players trickle on and off one at a time, so a single change produces several short-lived intermediate units. The regression handles this gracefully because every stint is weighted by its duration; a 2-second stint contributes almost nothing, and it also gets a 30-second minimum exposure when converting counts to rates, so one lucky bounce in a 2-second stint can't produce an infinite xGF/60.

One shift under the microscope

Here is Macklin Celebrini's fifth shift from San Jose at Detroit on January 16, 2026. He hops over the boards at 9:04 of the first period with Will Smith and Pavol Regenda, Dmitry Orlov and John Klingberg behind them. The Sharks keep the same five skaters out for 50 seconds while Detroit changes all three forwards halfway through:

San Jose Sharksaway · attacking unit

Pavol Regenda

LW

Macklin Celebrini

C

Will Smith

C

Dmitry Orlov

D

John Klingberg

D

Yaroslav Askarov

G

VS

Detroit Red Wingshome · final defending unit

Patrick Kane

RW

Andrew Copp

C

Alex DeBrincat

RW

Moritz Seider

D

Simon Edvinsson

D

John Gibson

G

The puck stays at one end. Celebrini tips the first look wide at 9:06. Smith misses a sharp-angle backhand six seconds later. After Detroit changes forwards, Celebrini puts a snap shot on John Gibson at 9:50. Four seconds later Smith taps in a one-foot backhand, Celebrini gets the primary assist, and the shift ends. The rink view reconstructs all four attempts:

The official NHL clip shows how the sequence ends:

Four attempts, 1.138 expected goals, and a proper beginning-middle-end. The first three chances total only 0.188 xG. Smith's final backhand is worth 0.950 by itself: one foot from the line, no angle, and no plausible world in which the goalie should recover. Celebrini is on the ice for all four attempts, takes two, and creates the last one.

This is also a clean demonstration of why a shift is not a stint. Celebrini's 50 seconds become two 25-second regression observations because Detroit changes from Larkin–Raymond–van Riemsdyk to Kane–Copp–DeBrincat at 9:29. The final stint's corrected ledger is below; columns are from Detroit's home perspective:

duration_seconds = 25          strength_state = 5v5
xgf_home = 0.0000              xga_home = 1.0069
gf_home  = 0                   ga_home  = 1
cf_home  = 0                   ca_home  = 2

The 1.0069 is Celebrini's 0.0569 snap shot plus Smith's 0.9500 goal. Nothing is estimated at the stint level; it is aggregation of the shot model. The word “corrected” matters because the goal occurs exactly at the 9:54 boundary. A half-open interval assigns it to the post-goal lineup unless goals get an end-boundary tie-break. We apply that tie-break here, then fix it in the pipeline below.

The xG math

Every unblocked shot attempt gets a goal probability from a Spark ML logistic regression (xg_model = logistic_v2 in the current build), trained on all attempts with a goalie in net:

P(goal | shot) = σ(β₀ + β·f)     where σ(z) = 1 / (1 + e^(-z))

f = [ log(1 + distance_ft),      # distance decays roughly exponentially
      angle_deg,                 # sharper angle, smaller net
      score_diff (clamped ±3),   # score effects on shot selection
      period,                    # period of game - more likely to score in 2nd period due to longer shift change
      is_rebound,                # secondary shot
      is_rush,                   # in transition, odd numbered sitution
      crossed_royal_road,        # lateral pass before shot
      is_behind_net,
      is_tip,                    # shot that is tipped in front of the goal, changing direction
      is_slap,
      is_wraparound,             # shooting with body from behind the net
      is_pp,                     # is power play
      is_sh ]

Celebrini's snap shot and Smith's one-foot backhand show the tradeoff between angle and distance. Celebrini's attempt is nearly straight on at 1.3 degrees, which helps, but it is still 44 feet from the goal. The model sees log(1 + 44.0) = 3.81, so the shot stays in low-probability territory at 0.057 xG. Smith's finish is the opposite case: no angle penalty, no meaningful distance, and the puck already at the edge of the crease. After calibration, that kind of chance is allowed to hit the model ceiling at 0.950 xG.

Raw logistic output isn't the final xG. It gets calibrated twice:

1. Per strength bucket (5v5, 5v4, 4v5, other): probabilities are scaled so each bucket's total xG equals its total goals. Power plays finish better than the same shot at evens, and the base model underrates that.
2. Globally on shots-on-goal with a goalie in net, so that league-wide Goals Saved Above Expected sums to approximately zero. If you skip this, every goalie in the league looks collectively brilliant or collectively cursed, and goalie WAR inherits the bias.

The receipts, by decile:

Empty-net situations skip the model entirely and use a multiplicative heuristic (xg_model = heuristic_fallback), because a logistic regression trained on goalie-in-net shots has no idea what to do when there's no goalie:

xg_raw = exp(-distance * 0.04) * angle_factor * shot_type_factor
       * rebound_factor * rush_factor * royal_road_factor

Only 1,822 of 112,888 attempts take this path. And where do the other 111,066 come from? Everywhere, but not uniformly:

Toggle that chart to "finishing" and the xG model stops being abstract: shot volume piles up at the blue line where little goes in, and conversion concentrates in a tight cone in front of the crease. Smith's goal sits at the brightest edge of that cone; Celebrini's preceding snap shot came from the high-volume, low-yield interior. Distance and angle carry most of the model because distance and angle carry most of hockey.

The RAPM iteration

Now the main event. RAPM asks: holding everyone else constant, how does xG flow change when a given player is on the ice? Mechanically, every stint becomes two rows in a design matrix, one per attacking perspective. Each row is really just a roster: five attacking skaters, five defending skaters, the defending goalie, team context, and home ice. The final 25 seconds of this shift say that San Jose generated 120.8 xGF/60 with the unit shown above, while those same twelve players pointed the other way produced zero. The targets are per-60 rates, and this stint is short enough that the safety valve in the conversion fires:

rate_factor = 3600 / max(duration, 30)         # 30s floor for tiny stints

duration = 25s → rate_factor = 3600 / 30 = 120   (floor applied)

Row 1 (SJS attacking):  y_xgf = 1.0069 × 120 = 120.8 xGF/60
                        y_gf  = 1      × 120 = 120.0 GF/60
Row 2 (DET attacking):  y_xgf = 0      × 120 = 0.0 xGF/60
                        y_gf  = 0.0

weight w = 25 (stint duration in seconds)

Without the floor, the goal would extrapolate to 144 GF/60, as if this unit scores every 25 seconds forever. Even with it, 120.8 xGF/60 is absurd; a good team sustains about 2.7 over a full game. That is fine. This row carries 25 seconds of weight in a regression built from 1,064 hours of hockey. Individually it is an anecdote; collectively, 891,832 of these rows are a season.

Assemble every 5v5 stint the same way and you get the full system: X is 891,832 rows by 2,043 columns (940 skaters twice, once as offense and once as defense, plus 98 goalies, 32 teams twice, and one home-ice column), W is a diagonal matrix of stint durations, and y is the per-60 rate vector. Ordinary least squares dies here, because teammates who always play together are nearly collinear and low-minute players produce wild coefficients. So we solve the ridge-regularized weighted normal equations instead:

β = (XᵀWX + λI)⁻¹ XᵀWy

λ = 25   for the xGF/60 target
λ = 80   for the GF/60 target
λ × 0.10 for team and home-ice columns

Each choice is doing a specific job:

• W (duration weights): a 25-second observation deserves more trust than a 2-second one, and a 60-second one deserves more than either. Longer stints are also less contaminated by line-change noise.
• λ = 25 vs λ = 80: shrinkage strength tracks noise. Shot generation (xGF) is a high-frequency, fairly repeatable process; goals (GF) are rare events ruled by shooting and save percentage luck, so the GF model gets more than three times the shrinkage. Smith's actual goal is worth 120 GF/60 in row 1, but λ = 80 ensures no player's reputation rides on which bounces happen to go in.
• λ × 0.10 on team/home controls: the nuisance columns get a lighter penalty on purpose. If team effects were shrunk as hard as player effects, systematic team context (coaching, systems, schedule) would leak into every player's coefficient instead of being absorbed by the team column.

After solving, offensive and defensive coefficient blocks are centered on their means, defense is sign-flipped so positive always means good, and each skater gets:

rapm_off_xgf_60 = β_off[i] − mean(β_off)       # xG created above average
rapm_def_xgf_60 = −(β_def[i] − mean(β_def))    # xG suppressed above average
rapm_xgf_60     = rapm_off_xgf_60 + rapm_def_xgf_60

Drive the multiplication yourself

Reading about a dot product is one thing; turning the crank is better. The widget below computes one row of Xβ for our stint using each player's actual published 2025-26 coefficient. Step through the terms one at a time, click any player to bench them (their x flips from 1 to 0 and their term drops out), and flip the perspective to compute the other row. The published coefficients are centered on the league mean, so the team, home-ice, and centering constants are folded into a single baseline term of 2.61, the league-average 5v5 xGF/60:

RAPM contribution waterfall

SJS attackingDET attacking

How the lineup moves the prediction

Start at league context, then add each player's isolated RAPM effect.

ŷ = 0.97 xGF/60

raises predictionlowers predictionClick a player bar to toggle their X entry

Show the underlying matrix math

X (1 × 12) × β (12 × 1) = ŷ (1 × 1)

One attacking row of the design matrix multiplied by its coefficient vector.

0.97 xGF/60

X: design row

Baseline

Regenda

Celebrini

Smith

Orlov

Klingberg

Kane

Copp

DeBrincat

Seider

Edvinsson

Gibson

1

1

1

1

1

1

1

1

1

1

1

1

×

β: coefficients

Baseline

Regenda

Celebrini

Smith

Orlov

Klingberg

Kane

Copp

DeBrincat

Seider

Edvinsson

Gibson

+2.61

−0.03

−0.31

+0.15

−0.05

−0.09

−0.14

−0.48

−0.28

−0.41

−0.08

+0.08

=

ŷ

0.97

xGF/60

The production row has 2,043 columns and 14 non-zero entries. This compact view shows the ten player columns, goalie column, and one baseline term that folds together the two team columns, home ice, and centering constants.

Lineup coefficients · click a player to toggle x between 1 and 0

β₀Baselineleague-average context: team, home ice, and centering constantsx = 12.61+2.61RegendaLW · offense columnx = 1β_off −0.03−0.03CelebriniC · offense columnx = 1β_off −0.31−0.31SmithC · offense columnx = 1β_off +0.15+0.15OrlovD · offense columnx = 1β_off −0.05−0.05KlingbergD · offense columnx = 1β_off −0.09−0.09KaneRW · defense columnx = 1β_def +0.14−0.14CoppC · defense columnx = 1β_def +0.48−0.48DeBrincatRW · defense columnx = 1β_def +0.28−0.28SeiderD · defense columnx = 1β_def +0.41−0.41EdvinssonD · defense columnx = 1β_def +0.08−0.08GibsonG · goalie columnx = 1β_g −0.08+0.08

Multiply next termShow allReset

ŷ = 0.97 xGF/60

model expectation ŷ

1.0

observed y in this stint (w = 25 of 3.8M seconds)

120.8

With everyone on the ice, the model expects this San Jose unit to generate about 0.97 xGF/60. Start at the 2.61 league baseline: the five Sharks' offensive terms subtract 0.33, Detroit's defensive terms subtract another 1.40, and Gibson gives 0.08 back. The observed row is 120.8. That residual looks ridiculous because it is, and then it is multiplied by only 25 seconds of weight out of roughly 3.8 million. One spectacular shift argues with the model; it does not win.

So what did the final 25 seconds in Detroit actually do? They nudged five Sharks offensive coefficients and five Red Wings defensive coefficients toward “ice tilts when these ones are out there,” by an amount proportional to 25 seconds of weight and inversely proportional to everything else those ten players did across the rest of the season. One stint moves a coefficient by a rounding error. A thousand shifts of consistently winning your minutes moves it to the top of this chart:

What the model actually thinks

A few things that are interesting in the results, in the spirit of arguing with a spreadsheet and ML model:

There are two species of RAPM monster. Brady Tkachuk (+1.30) and Blake Lizotte (+1.21) get there by creating: their offensive component dwarfs their defensive one. J.J. Moser gets to functionally the same total (+1.28) with a nearly inverted profile: +1.15 of his value is chance suppression across 1,407 5v5 minutes. Classic plus-minus could never distinguish these players from their linemates; the whole point of the design matrix is that Moser's defensive column gets credit only for suppression that persists after his goalie, his partners, and his team context are accounted for.

The leaderboard is refreshingly weird. Alongside Hintz, Kyrou, and Hagel sit Linus Karlsson, Pontus Holmberg, and Mavrik Bourque. Some of that is real (middle-six players who quietly win their minutes), and some of it is single-season RAPM being noisy for players with sheltered deployment. This is why the WAR model downstream blends RAPM at 55% with individual shot creation rather than swallowing it whole, and why many sports oriented data scientists pools multiple seasons. Season long RAPM can often bias to a very specific set of linemates, coaching sprategies, deployments, and can miss nuances.

xGF and GF RAPM disagree, and that's the design working. Moser's rapm_gf_60 (+2.75) is more than double his xG-based number. Somewhere between shot quality and goals, Tampa got great goaltending and shooting luck in his minutes. The paired targets with different λ values let you see the luck rather than average it away.

From coefficients to player cards

RAPM is deliberately narrow: it estimates how 5v5 chance flow changes with a player on the ice. WAR asks a broader question: how many wins did the player add relative to a replacement-level skater? In this build, replacement level for skaters is the 25th-percentile component rate within position group, forwards and defensemen separately, among skaters with at least 10 games.

The version here is a small, transparent build rather than a claim to have solved hockey. It turns several player signals into the same currency, wins, then adds them:

Each raw signal is compared within a position group, translated into goal value, then divided by 5.15 goals per win. So a component worth roughly +2.6 goals becomes about +0.5 WAR, and a +5.15-goal edge becomes +1.0 WAR. After that conversion, the 25th-percentile skater baseline is shifted to zero, so positive numbers are wins above replacement and negative numbers are wins below that replacement line. That is why a player can have a soft 5v5 RAPM coefficient and still grade out as excellent if the other columns are doing enough work.

For the three Sharks in the featured shift, the accounting looks like this:

Additive WAR ramp

Every player starts at zero WAR. Each component moves the running total; teal steps add value, red steps subtract value.

positive WARnegative WAR

Macklin Celebrini

C · 82 games

+5.22 WAR

final total

Will Smith

C · 69 games

+2.59 WAR

final total

William Eklund

LW · 78 games

+2.42 WAR

final total

Each step is already measured in wins above replacement, so the components can be added left to right. Celebrini's 5.22 WAR is not a claim that his 5v5 RAPM was elite. It was not. The model credits a 99th-percentile season because elite playmaking, power-play value, finishing, transition, defense, and penalty impact overwhelm a soft even-strength offensive coefficient. Eklund's card tells the inverse finishing story: strong two-way and play-driving value survives a negative shooting component.

Choose any skater from the 2025–26 season to render their card:

Choose a playerMCMacklin Celebrini· SJS · C

Loading the full 2025–26 player list…

Macklin Celebrini

SJS · C · 82 games

5.22

WAR

99

WAR percentile

92

Composite

EV offense+0.12 WAR

EV defense+0.57 WAR

Transition+0.40 WAR

Power play+1.22 WAR

Playmaking+1.71 WAR

Shooting+0.93 WAR

Penalties+0.27 WAR

−2.0 WAR0+2.0 WAR

Shared scale for every component and every player

A 99th-percentile WAR season built less by even-strength RAPM than by playmaking, power-play creation, finishing, and drawing more penalties than he takes.

Goalie cards

Goalies need their own accounting. A skater's value can be split across offense, defense, transition, power play, playmaking, shooting, and penalties. A goalie's job is narrower and stranger: face a sequence of shots, inherit the chaos of screens and rebounds and rushes, and turn expected goals against into actual goals against.

The goalie card starts with GSAx, Goals Saved Above Expected. For every shot with a goalie in net, the xG model estimates the chance of a goal. Add those probabilities up and compare them with the goals actually allowed:

GSAx = expected goals against − actual goals against

Then the model shrinks GSAx toward zero based on shots faced, because a goalie with 200 shots faced should not get the same certainty as one with 1,800. That shrunk GSAx is divided by 5.15 goals per win, shifted above the 35th-percentile replacement goalie, and blended with goalie RAPM. This is intentionally different from the skater baseline: skaters use the 25th percentile within position group, while goalies use the 35th percentile among qualified goalies. The actual sequence is:

GSAx WAR        = p35_shift(shrunk GSAx / 5.15)
goalie RAPM WAR = clamp((rapm_xga_60 − p35_rapm_rate) × TOI_hours / 5.15)

goalie WAR raw  = 0.80 × GSAx WAR + 0.20 × goalie RAPM WAR
goalie WAR      = goalie WAR raw − p35(goalie WAR raw)

The RAPM piece asks whether expected goals against changes when that goalie is in net after accounting for the skaters, teams, and home ice. The GSAx piece still does most of the work, but RAPM keeps the card from treating every environment as equally difficult.

Choose a goalieISIlya Sorokin· NYI · G

Loading the full 2025–26 goalie list…

Ilya Sorokin

NYI · G · 55 games

4.33

Goalie WAR

100

WAR percentile

90

Composite

GSAx WAR+4.98 WAR

RAPM WAR+2.00 WAR

How the WAR is built

0.80 × +4.98 + 0.20 × +2.00 = 4.39

4.39 − 0.06 p35 shift = 4.33 WAR

Adj SV%

.907

Raw SV%

.906

GSAx

+28.0

xGA

172.0

Ilya Sorokin's goalie card starts with shot quality faced, shrinks GSAx for sample size, blends in goalie RAPM, and then asks how far above a replacement-level goalie that result sits.

This work is inspired by the excellent content at Evolving Hockey and HockeyStats, simplified to a single season. RAPM's regularized-regression lineage traces back to basketball's adjusted plus-minus work, and the shrinkage story is the same in any sport: when your signal is rare events distributed across shared credit, the honest move is to pull every estimate toward zero and let sustained evidence drag it back out. Celebrini's 50-second shift produced four attempts, 1.138 expected goals, and the winning geometry of Smith's one-foot finish. The model enjoyed the goal as much as anyone. It just refuses to plan around it.