html`<div style="display:flex;flex-wrap:wrap;gap:12px">${OUT.map(v => {
const d = ppc.distributions[v];
const bins = d.observed.map((c, k) => ({x1: d.edges[k], x2: d.edges[k + 1], c}));
const mid = k => d.log_x ? Math.sqrt(d.edges[k] * d.edges[k + 1]) : (d.edges[k] + d.edges[k + 1]) / 2;
const reps = d.replicates.flatMap((r, j) => r.map((c, k) => ({rep: j, x: mid(k), c})));
return Plot.plot({
width: panelW, height: 200, marginLeft: 44, title: d.label,
x: {type: d.log_x ? "log" : "linear", label: null, tickFormat: d.log_x ? (t => t.toFixed(2)) : undefined, ticks: 5},
y: {label: "↑ seasons", grid: true},
marks: [
Plot.rectY(bins, {x1: "x1", x2: "x2", y: "c", fill: "#c9c8c3", inset: 0.5}),
Plot.line(reps, {x: "x", y: "c", z: "rep", stroke: REP, strokeOpacity: 0.25, strokeWidth: 1}),
Plot.ruleY([0])
]
});
})}</div>`The Soccer Factor Model, explained
How Andorra and Göbel separate a forward from his team, and what we changed to fit their model to public data
An account of the Soccer Factor Model (Andorra & Göbel, 2024), our four-skill season-level version of it, how we checked it and what it found.
This page explains the model behind the explorer that Wyatt and I built. The data dictionary describes the data. The code is in src/sfm/ and scripts/. Given the Kaggle Understat dump, the scripts rebuild the fit, the tables and the figures.
Skill and the team around it
A striker’s goal count measures the striker and his team at once. A team that wins the ball high and plays most games against weaker sides feeds its forward more chances, and better ones, than a team fighting relegation. Andorra and Göbel (2024) call the raw count observed performance, OP. They model it as a function of the player’s own skill, \(\alpha\), and his team’s performance, TP:
\[\text{OP} = g(\alpha, \text{TP}; \theta).\]
A scout wants to know whether the player carries the team or the team carries the player. To answer that, you have to pull \(\alpha\) out of OP.
The idea comes from finance. A fund manager who beats the market may have picked good stocks. Or he may have held the kinds of stock that rose that year, small firms say, or cheap ones. Finance separates the two by regressing the manager’s returns on these market-wide factors. The intercept, the return the factors cannot explain, is the manager’s alpha. The Soccer Factor Model treats team strength as the factors. Whatever is left of the player’s performance is his alpha.
The model in the paper
The paper follows 144 strikers in England, Germany, Italy and Spain from 2000/01 to 2023/24. That is about 33,000 player-matches, scraped from the German magazine kicker.
The outcome is goals per player per match. Only 1.5% of player-matches produce three goals or more, so the authors cap the count at three. That leaves four ordered categories, 0, 1, 2 and 3+. They model these with an ordered logistic likelihood, which turns a single latent score into a probability for each category by way of a set of cutpoints. The first cutpoint is fixed so that the model is identified.
Two factors stand for team strength. One marks home games. The other is the gap in league points between the player’s team and its opponent so far that season. Both are known before kickoff. The authors tried more factors and found that the points gap had already absorbed them.
Skill in the paper has three parts:
\[\alpha_{i,m,s} = \underbrace{\mu_b + \delta_{b,i}}_{\text{talent}} + \underbrace{h^W(D_{i,m,s})}_{\text{within-season}} + \underbrace{h^C(N^S_{i,s})}_{\text{cross-season}}.\]
The first part is talent: a common level \(\mu_b\) plus a player effect \(\delta_{b,i}\), with the player effects summing to zero. The second, \(h^W\), tracks form within a season. Its input is the number of days since the player’s first match of the season. It adds a short-lengthscale curve to a medium one, so it can pick up a hot streak as well as slow fatigue. The third, \(h^C\), follows the career. Its input is the number of seasons the player has been in the sample.
Each curve is a Gaussian process. An exact Gaussian process slows to a crawl with tens of thousands of matches. The paper therefore uses a Hilbert-space approximation, or HSGP, which builds each curve from a fixed set of basis functions under a Matérn-5/2 kernel.
To compare players, the authors borrow baseball’s wins above replacement. Performance above replacement, PAR, sets a player’s predicted goals against those of a hand-picked pool of ordinary strikers. Skill above replacement, SAR, makes the same comparison with every team factor set to zero, as though all teams were equally strong. The difference, SAR − PAR, measures how far the team inflated or deflated the player’s numbers.
Messi comes out as the most skilled striker of the period. Haaland’s estimate is high and very uncertain, because his career is short. The three most skilled players, Messi, Haaland and Ronaldo, also gained most from their teams. Skill peaks about four seasons into a career. In leave-one-out cross-validation the full model predicts better than two simpler versions of itself, though the authors call its overall predictive performance “not outstanding”.
What we changed
We wanted to fit the paper’s ideas to data anyone can download. We used Understat for the top five leagues from 2014/15 to 2023/24, taken from a Kaggle archive, and match results from football-data.co.uk.
The archive holds every match and every shot, but no lineups. When a forward plays without taking a shot or making a key pass, his appearance leaves no trace. That happens in 29% of appearances. The paper’s per-match model needs those empty games, so we could not fit it as written. Season totals, minutes included, are complete. We fitted the model to player-seasons instead, and used the change to measure four skills where the paper measures one.
Four skills
Goals per 90 minutes split exactly into three ratios:
\[\frac{\text{goals}}{90} = \frac{\text{shots}}{90}\times\frac{\text{xG}}{\text{shot}}\times\frac{\text{goals}}{\text{xG}}.\]
Expected goals, or xG, is the probability that a shot like the one taken goes in. We treat each ratio as a separate ability, though each also carries the team’s influence and luck. The first counts how often a forward gets a shot away. Wyatt calls a forward who is good at this a Shooter. The second measures where the shots come from: the Fox. The third asks whether he scores more than an average player would from the same chances: the Sniper. We added a fourth skill, the Playmaker, which counts the shots a forward sets up for team-mates. It adds nothing to his own goals, but it is plainly part of a forward’s game.
The team factor, and why it needs help
For each player-season we measure team strength as the paper does, before kickoff. For every league match we take the gap in points per game between the team and its opponent, using only earlier matches that season. We then average those gaps over the team’s season. Call the result TS. A side that sits about one point per game clear of its opponents on average, roughly the strongest team in a league, has TS near 1.
A season-level measure brings a problem the paper does not face in the same form. Good forwards do not land at clubs at random. The best of them play for the best teams. If the model assumes a forward’s skill is unrelated to the strength of his teams, it cannot tell “this team makes forwards score” from “this team hires forwards who score”. The team coefficient then absorbs part of the players’ skill, and the best players’ skill is understated.
We deal with this using a device from panel econometrics, the Mundlak term. For each player we also include his average team strength across all his seasons, \(\overline{\text{TS}}_p\). The average soaks up the systematic link between a player and the kind of team he plays for. What is left for the team coefficient comes from changes within a career: the same forward on a stronger or weaker team in different seasons. That is the comparison we want, because the player’s ability is held fixed in it.
For player-season \(r\) and skill \(k\), the linear predictor is
\[\eta_{k,r} = \mu_k + \lambda_{k,\ell(r)} + \beta_k\,\text{TS}_r + \gamma_k\,\overline{\text{TS}}_{p(r)} + h_k(N^S_r) + \alpha_{p(r),k}.\]
Here \(\lambda_{k,\ell}\) is a league effect, and the five league effects sum to zero. \(\beta_k\) is the within-player team effect. \(\gamma_k\) is the extra association between a player’s skill and the strength of the teams he plays for. \(h_k\) is a Matérn-5/2 HSGP over seasons observed. The four \(\alpha\)s are correlated, with an LKJ prior on their correlations:
\[\alpha_p \sim \text{MVN}(0, \Sigma), \qquad \Sigma = \text{diag}(\sigma)\,\Omega\,\text{diag}(\sigma), \quad \Omega \sim \text{LKJ}(2).\]
A player’s skill is everything about him that does not change with his team in a given season: \(\theta_{p,k} = \alpha_{p,k} + \gamma_k\,\overline{\text{TS}}_p\). That is the quantity the explorer shows and the one SAR is built on.
Likelihoods
| Skill | Data | Likelihood |
|---|---|---|
| Shooter | non-penalty shots | \(\text{NegBin}\big(\tfrac{\text{min}}{90}e^{\eta_\text{shoot}},\ \phi_s\big)\) |
| Fox | non-penalty xG per shot | \(\log(\text{xG}/\text{shot}) \sim N\big(\eta_\text{fox},\ \sigma_f/\sqrt{\text{shots}}\big)\) |
| Sniper | each non-penalty shot | \(\text{Bernoulli}\big(\text{logit}^{-1}(\text{logit}\,\text{xG}_i + \eta_\text{sniper})\big)\) |
| Playmaker | key passes | \(\text{NegBin}\big(\tfrac{\text{min}}{90}e^{\eta_\text{create}},\ \phi_c\big)\) |
The standard deviation for the Fox shrinks with the square root of the number of shots. A season of 120 shots counts for more than a season of 12.
Finishing from three numbers a season
The Sniper likelihood works shot by shot. Shot \(i\) has xG \(p_i\), and it goes in with probability \(\text{logit}^{-1}(\text{logit}\,p_i + \eta)\), where \(\eta\) is the forward’s finishing in that season on the log-odds scale. At \(\eta = 0\) the forward converts exactly at xG. At \(\eta = 0.1\) a 10% chance becomes a 10.9% chance.
Summed over a season’s shots, with \(y_i = 1\) for a goal, the log-likelihood is
\[\ell(\eta) = G\,\eta - \sum_i \log\!\big(1 + e^{\,\text{logit}\,p_i + \eta}\big) + \text{const}, \qquad G = \sum_i y_i .\]
The outcomes enter only through \(G\), the season’s goals. Which shots went in makes no difference once the xG of every shot is known. The xGs still matter, through the sum of logs, and that sum is what made the exact fit slow. It has a term for each of 182,391 shots, and the sampler evaluates it thousands of times. A fit took close to an hour.
The first derivative of the sum at \(\eta = 0\) is \(X = \sum_i p_i\), the season’s xG. The second is \(V = \sum_i p_i(1-p_i)\). Expanding to second order gives
\[\ell(\eta) \approx (G - X)\,\eta - \tfrac{1}{2}\,V\,\eta^2 + \text{const}.\]
A season now needs three numbers: goals \(G\), xG \(X\) and \(V\). They are the columns fin_goals, fin_xg and fin_var in the model table, and they replace the shot list. A fit takes under five minutes.
Completing the square shows what the approximation says. It is the likelihood of a single normal observation,
\[z = \frac{G - X}{V} \sim N\!\big(\eta,\ 1/\sqrt{V}\big).\]
So \(z\) is the finishing that one season would suggest on its own, and \(1/\sqrt{V}\) is its standard error. Messi in 2016/17 scored 31 non-penalty goals from 21.7 xG, with \(V = 14.9\): \(z = 0.63\), with a standard error of 0.26. Haaland in 2022/23 scored 29 from 27.4, with \(V = 16.3\): \(z = 0.10\), standard error 0.25. The model does not take either number at face value. It pools each forward’s seasons into one skill, and pulls that skill toward the average of all forwards, the more so when \(V\) is small.
\(V\) measures how much a season can say about finishing. Each shot adds \(p(1-p)\) to it. A 50% chance adds 0.25, a 90% tap-in 0.09, a 5% shot from distance 0.05. Over all shots in the data, \(V\) is 0.66 of the xG. It is lower for a forward whose xG comes from close range: about 0.6 in Haaland’s seasons against 0.67 in Messi’s. A goal expected from a tap-in says less about finishing than a goal expected from a harder chance, because a tap-in is scored by nearly everyone.
\(V\) also explains an error in an earlier version of the model. That version gave goals a Poisson likelihood with xG as the exposure. A Poisson takes the variance of a season’s goals to be \(X\). The shots themselves imply \(V\), about two-thirds as much. The Poisson model expected more luck than the shots could produce, and it put the real differences between finishers down to luck. It estimated the spread of finishing skill at 0.05. The model here puts it at 0.14.
For a typical season of about 40 shots, the approximation differs from the exact likelihood by about 0.002 nats at a finishing effect of ±0.2, and by about 0.02 at ±0.4. The error grows with the number of shots. In the busiest seasons, around 130 shots, it reaches about 0.008 and 0.06. Adding the third-order term would be more accurate near zero, but the resulting likelihood rises without limit for large \(\eta\), so it cannot be used. The predictive checks below use the exact shot-by-shot likelihood. In the pipeline the quadratic enters the model as a pm.Potential. The notebook writes it as the normal observation of \(z\). That has the same dependence on \(\eta\) and lets PyMC simulate data from the model.
Skill above replacement
SAR asks one question of each forward. How many more non-penalty goals per 90 minutes would he score than an ordinary forward, if both played in an average league, for an average team, at an average point in their careers? Answering it takes three steps: a scoring rate on that level field, a replacement level, and the difference between them.
The level-field rate. Set the league effect, team strength and career curve to zero. What is left of each skill’s linear predictor is \(\mu_k + \theta_{p,k}\), and the rate is
\[\text{rate}_p = \underbrace{e^{\mu_\text{shoot} + \theta_{p,\text{shoot}}}}_{\text{shots per 90}} \times \underbrace{e^{\mu_\text{fox} + \theta_{p,\text{fox}}}}_{\text{xG per shot}} \times \underbrace{c\big(\mu_\text{sniper} + \theta_{p,\text{sniper}}\big)}_{\text{goals per xG}} .\]
The first two factors follow from the likelihoods. The third needs a choice. Finishing shifts each shot’s log-odds, and the same shift adds more goals per xG to long shots than to tap-ins. So goals per xG depends on the kind of chances taken. We fix one mix of chances for every forward: all the non-penalty shots in the data. Then
\[c(\eta) = \frac{\text{mean}_i\ \text{logit}^{-1}(\text{logit}\,p_i + \eta)}{\text{mean}_i\ p_i} \approx 1 + \frac{V}{X}\,\eta = 1 + 0.66\,\eta ,\]
and the linear form is the same second-order approximation as before. A finishing skill of +0.1 is worth about 7% more goals per xG. Fixing the mix values finishing the same way for every forward. The chances a forward actually takes are counted once, in the Fox.
The table takes the rate apart for three forwards and for the replacement pool, using posterior means.
| Shots per 90 | xG per shot | Goals per xG | Goals per 90 | |
|---|---|---|---|---|
| Erling Haaland | 3.14 | 0.191 | 0.98 | 0.59 |
| Lionel Messi | 4.73 | 0.107 | 1.12 | 0.57 |
| Harry Kane | 3.70 | 0.132 | 1.08 | 0.53 |
| Replacement pool | 2.33 | 0.114 | 0.99 | 0.26 |
Haaland and Messi reach almost the same rate by different routes. Haaland takes fewer shots from much better places and finishes them at about xG. Messi shoots half as often again from ordinary places and converts 12% above xG. (The last column is the mean of the product over posterior draws, so it can differ slightly from the product of the means.)
The replacement level. The pool holds the 119 forwards whose level-field rates fall between the 40th and 60th percentiles of the 593. Their average rate is 0.255 goals per 90. The paper chose its replacement strikers by hand; we chose ours by rule. Every forward in the data has at least 6,000 minutes in a top-five league, so “replacement” here means a typical regular forward in those leagues, not a fringe player. Against a weaker pool every SAR would be larger. A different pool moves every forward’s SAR by the same amount in each draw, so it does not change the order. The explorer lets you choose another pool.
The difference. \(\text{SAR}_p = \text{rate}_p - \overline{\text{rate}}_\text{pool}\), computed draw by draw, so the interval carries uncertainty about the forward and about the pool. For Haaland it is 0.59 − 0.26 = 0.33 goals per 90. Over 34 full matches, about a season for a regular starter, that is 11 non-penalty goals more than a replacement forward would score on the same field.
SAR leaves out penalties, creating chances for others, and the forward’s real league, teams and age. It measures skill at scoring from open play, not a forward’s whole value.
PAR and the team’s share. PAR replaces the level field with the forward’s own teams. For each of his seasons we add \(\beta_k\,\text{TS}\) to the three scoring skills, compute the rate, and average the rates by minutes. League and career stage stay at their averages, so the team factor is the only difference from SAR, as in the paper. The pool’s rate is recomputed the same way, with its own teams. The team lift is the difference between the two rates: for Haaland, 0.75 − 0.59 = 0.16 goals per 90. The pool’s average lift is 0.007, because its forwards played for teams close to average. So
\[\text{SAR} - \text{PAR} = -\big(\text{lift}_p - \overline{\text{lift}}_\text{pool}\big),\]
and a negative value means a forward’s teams helped him more than the pool’s teams helped them. Haaland’s is −0.16.
The table sets out where we depart from the paper.
| Paper | Here | |
|---|---|---|
| Unit | player-match | player-season (no public lineups) |
| Outcome | goals, ordered logistic | shots, xG/shot, goals given each shot’s xG, key passes |
| Skill | one \(\alpha\) per player | four correlated skills (LKJ prior) |
| Player effects | sum to zero | correlated, MVN with mean zero |
| Team factor | home + pre-match points gap | pre-match points-per-game gap, averaged over the season, plus a Mundlak term |
| Career curves | per player, within- and cross-season HSGPs | one population HSGP per skill, over seasons |
| Penalties | included | excluded |
| Players | 144 hand-picked strikers | 593 forwards: any player Understat lists at F for most of his minutes, ≥ 6,000 minutes, plus a few of Wyatt’s picks |
| Replacement pool | 102 hand-picked strikers | forwards between the 40th and 60th percentile of level-field skill |
Checking the model
We checked the priors before fitting. Seasons simulated from the priors alone had a median of 39 shots, against 42 in the data, with a 90% range of 6 to 176 against 10 to 100. Key passes and xG per shot looked the same: the right scale, wider than the data. Weak priors should behave like that.
Next we simulated a whole dataset from the model and fitted the model to it. The simulation had known skills, correlations, team effects, Mundlak terms and a humped career curve. We set the true spreads near the real estimates, including 0.13 for finishing. The fit produced no divergences. The 83% intervals covered 81% of the true Shooter values, 85% for the Fox, 80% for the Sniper and 83% for the Playmaker. The spreads came back within 0.01 of the truth. The team and Mundlak coefficients came back within 0.02.
Finishing was again the hardest skill to recover. Its posterior means correlated 0.55 with the true values, against about 0.95 for each of the other skills. An earlier run with a smaller finishing spread, 0.05, did worse: a correlation of 0.17. A forward’s finishing shows only in a few dozen goals a season, and the model sees that.
The real fit ran 4 chains of 1,500 draws in nutpie, with target_accept = 0.95. It produced one divergence in 6,000 draws. The largest \(\hat R\) was 1.007, and the smallest effective sample size 368, both on an element of the Cholesky factor behind the skill correlations. The skills themselves have \(\hat R\) at most 1.004 and effective sample sizes above 3,900.
Posterior predictive checks
Convergence shows that the sampler explored the posterior. It says nothing about whether the model fits. For each check we draw replicated data from the posterior predictive distribution and set them beside what the players did.
Checks on the data a model was fitted to flatter it, so for calibration we use leave-one-out versions. Pareto-smoothed importance sampling reweights the posterior draws as if each player-season had been left out. For each outcome, at most six of the 3,950 seasons have weights too unstable to trust (Pareto \(k\) above 0.7).
The first check compares whole distributions. The grey bars show the observed player-seasons. Each thin blue line is one of 40 datasets the model generated. For goals, the replicates keep each season’s real shots and redraw which of them went in.
The replicated datasets match the observed distributions for all four outcomes, including the long right tails of shots and goals.
The second check asks whether each season looks like a typical draw from its own leave-one-out prediction. For every player-season we compute the share of predicted values below the observed one, the probability integral transform, or PIT. A calibrated model gives PIT values spread evenly between 0 and 1. If the predictions are too narrow, values pile up at the ends. If they are too wide, the ends thin out. The dashed line marks an even spread, and the band shows how far a calibrated model would wander from it by chance.
The leave-one-out histograms sit mostly inside the band. Shots, key passes and xG per shot are well calibrated: their 83% leave-one-out intervals hold 83% to 86% of seasons. The in-sample figures are four to six points higher, which is the flattery the leave-one-out version removes. Goals are slightly over-covered, at about 90%. There are a few more low PIT values than an even spread allows, so some seasons scored fewer goals than the model expected from their chances.
The third check looks for patterns the model leaves behind across team strength, leagues and stages of a career. Here goals are simulated through the whole chain: a replicated number of shots, then a replicated average xG per shot, then replicated goals. This tests how team strength raises the number and quality of chances as well as their conversion. To keep it fast, the chain draws goals as a binomial with the season’s average xG, ignoring the spread of chance quality within a season. Each panel compares observed non-penalty goals per 90 (dark dots) with the model’s 83% predictive interval (blue bars).
The last check looks at players one by one. For each forward we add up his replicated goals across all his seasons, again through the whole chain, and compare the result with his actual career total. Points on the diagonal match exactly. The grey bars are 83% intervals, and orange marks a player whose total falls outside his.
The checks found no large misfit. Shots, key passes and chance quality are well calibrated out of sample, and goals given shots are slightly over-covered. The full chain reproduces scoring across team strength, leagues and career stage, with intervals that run wide.
Findings
How much the team matters
These are the team coefficients for a team one point per game clear of its opponents before kickoff, with 83% intervals. \(\beta\) is the within-player effect. \(\gamma\) is the Mundlak term.
| Skill | \(\beta\) (within player) | \(\gamma\) (Mundlak) |
|---|---|---|
| Shooter | 0.087 [0.073, 0.101] | 0.167 [0.119, 0.214] |
| Fox | 0.166 [0.150, 0.182] | −0.007 [−0.057, 0.043] |
| Sniper | 0.096 [0.062, 0.129] | 0.059 [0.009, 0.108] |
| Playmaker | 0.093 [0.076, 0.110] | 0.313 [0.260, 0.365] |
The Mundlak term matters for shooting and creating. Forwards who shoot and create a lot cluster at strong clubs, and without the term the model would have credited that to the clubs. An earlier version of this model, without it and with the season’s final points as the team factor, put the shooting coefficient at 0.150 and the creating coefficient at 0.173. Within a player’s own career both are a little over half that.
Chance quality is different. Its Mundlak term is close to zero, so better teams do not hire forwards who find better chances. The same forward gets better chances on a better team.
Put together, a forward moving to a team one point per game stronger scores about 1.37 times as many non-penalty goals. Chance quality accounts for about half of that, shot volume for under a third, finishing for the rest. These are still associations. Our pre-match measure, like the paper’s, includes the player’s own goals in earlier matches that season. The Mundlak term also assumes that a player’s skill relates to team strength only through his career average.
Skill spreads
The spread of each skill across forwards (the standard deviation of \(\theta\) on the log or log-odds scale):
| Shooter | Fox | Sniper | Playmaker |
|---|---|---|---|
| 0.29 | 0.30 | 0.14 | 0.33 |
A forward one standard deviation above average at shooting takes about 33% more shots. One standard deviation above average at finishing scores about 9% more goals than his chances would give an average forward. Finishing varies less than the other skills, but it is not negligible. The Poisson version of the model, described above, put its spread at 0.05.
The model is still unsure about most individual finishers. The explorer gives a forward three Sniper stars unless his 83% interval excludes zero, and 536 of the 593 keep three. Messi and Son Heung-min get five. Harry Kane gets four. Haaland scored 103 non-penalty goals from 94.9 xG and keeps three: his interval runs from −0.14 to 0.11.
Skills that go together
Posterior correlations between the player skills \(\theta\), with 83% intervals:
| Fox | Sniper | Playmaker | |
|---|---|---|---|
| Shooter | +0.23 [+0.20, +0.26] | −0.07 [−0.17, +0.03] | −0.12 [−0.15, −0.10] |
| Fox | −0.40 [−0.49, −0.30] | −0.57 [−0.60, −0.55] | |
| Sniper | +0.45 [+0.36, +0.54] |
Creators take worse shots. The correlation between Fox and Playmaker, −0.57, is the most precisely estimated in the table. A penalty-box striker lives on tap-ins and sets up little. A forward who drops deep to create shoots from further out.
Forwards with very good chances convert slightly below their xG, and creators convert slightly above it. Both patterns could be skill. They could also come from Understat’s xG model, if it rates very good chances a little too high and long shots a little too low. A season-level model cannot tell which. The best-rated finishers include James Rodríguez, Suso and Ruslan Malinovskiy, all creators who often shoot from distance.
The most skilled forwards
The top of the SAR ranking, in non-penalty goals per 90 above the replacement pool on a level field. SAR − PAR is minus the team’s share, as defined above.
| SAR | 83% interval | SAR − PAR | |
|---|---|---|---|
| Erling Haaland | 0.33 | 0.24 – 0.44 | −0.16 |
| Lionel Messi | 0.31 | 0.25 – 0.38 | −0.20 |
| Kylian Mbappé | 0.31 | 0.23 – 0.39 | −0.20 |
| Darwin Núñez | 0.31 | 0.18 – 0.45 | −0.11 |
| Cristiano Ronaldo | 0.31 | 0.24 – 0.37 | −0.18 |
| Robert Lewandowski | 0.30 | 0.24 – 0.37 | −0.21 |
| Sergio Agüero | 0.28 | 0.21 – 0.37 | −0.17 |
| Luis Suárez | 0.27 | 0.20 – 0.35 | −0.18 |
The intervals overlap almost completely. Read the table as a group of elite forwards. The order within it is mostly noise. Núñez’s interval is the widest, because he played fewer than 4,000 minutes. The paper made the same point about Haaland: a short career earns a wide interval.
The Mundlak term moved this table. The model no longer credits these forwards’ teams with part of their skill. Haaland, Messi and Ronaldo, the paper’s top three, are all now in the top five.
Lewandowski, Mbappé, Messi, Suárez and Ronaldo have the largest team lift, 0.18 to 0.21 goals per 90. They also played for the strongest teams. Team strength multiplies the scoring rate, so the same team adds more goals to a better forward; the model builds that in. At the other end, the model attributes the shortfall of Simy, Charlie Austin, Aleksandar Mitrović and Jermain Defoe to weak teams.
Messi’s creating, 2.0 times a typical forward’s, adds nothing to his goals, and we leave out penalties. Both count against him here. Our window also includes his two seasons at Paris Saint-Germain, at 34 and 35. The model has no individual career curves to set those seasons apart.
Careers
The average multiplier for each skill, by season observed. Season 1 is a player’s first season in the 2014/15–2023/24 window.
| Season 1 | Season 3 | Season 5 | Season 10 | |
|---|---|---|---|---|
| Shooter | 1.04 | 1.02 | 1.00 | 0.93 |
| Fox | 0.95 | 0.98 | 1.02 | 1.05 |
| Sniper | 1.04 | 1.05 | 0.98 | 0.87 |
| Playmaker | 1.00 | 1.00 | 1.00 | 1.01 |
Across their seasons in the window, forwards shoot less often and from better places, and their finishing fades. Scoring peaks around the third season observed, about 3% above average, and ends about 11% below it.
The curve has a flaw. Seasons are counted from the start of our window, so Messi’s first season is 2014/15, his eleventh as a professional. The curves show how forwards changed over one decade. They are not age curves.
What remains
Understat’s per-match lineups would let us fit the paper’s model as written: goals per match, an ordered logistic likelihood, the home and pre-match points factors, and form within the season. We already hold the match ids we would need.
Age in place of seasons observed would fix the career curves. Curves for individual players would stop a veteran’s decline years from dragging down his estimated skill. That was the paper’s main reason for fitting them. We have not yet repeated the paper’s cross-validation against simpler models, though the saved log-likelihood now makes that straightforward. And the Fox and Sniper skills rest on Understat’s xG model, so they inherit its errors, possibly including the correlations between them.
Reproducing it
The quickest route is reproduce.ipynb. It fits the model in plain PyMC from a single 0.4 MB table, data/model_table.csv, and runs the prior check, the fit, SAR and a parameter-recovery test with pm.do.
The full pipeline, from football_models_with_wyatt/ with the conda environment in environment.yml:
python scripts/01_fetch.py odds # football-data.co.uk CSVs
python scripts/02_team_matches.py # team factors
python scripts/03_build_panel.py # player tables (needs the Kaggle Understat dump in data/understats/)
python scripts/04_fit_skills.py prior # prior predictive summary -> results/
python scripts/04_fit_skills.py recovery # simulation recovery -> results/
python scripts/04_fit_skills.py real # the fit -> idata/
python scripts/05_sar.py
python scripts/06_export_web.py # web/*.json, including the predictive checks
python scripts/07_model_table.py # data/model_table.csv for reproduce.ipynb
python -m pytest tests -qAndorra, A. and Göbel, M. (2024). Unveiling True Talent: The Soccer Factor Model for Skill Evaluation. arXiv:2412.05911.