(item_response_nba)=
:::{post} Apr 17, 2022 :tags: hierarchical model, case study, generalized linear model :category: intermediate, tutorial :author: Austin Rochford, Lorenzo Toniazzi :::
This tutorial shows an application of Bayesian Item Response Theory {cite:p}fox2010bayesian to NBA basketball foul calls data using PyMC. Based on Austin Rochford's blogpost NBA Foul Calls and Bayesian Item Response Theory.
Our scenario is that we observe a binary outcome (a foul being called or not) from an interaction (a basketball play) of two agents with two different roles (the player committing the alleged foul and the player disadvantaged in the play). Moreover, each committing or disadvantaged agent is an individual which might be observed several times (say LeBron James observed committing a foul in more than one play). Then it might be that not only the agent's role, but also the abilities of the single individual player contribute to the observed outcome. And so we'd like to estimate the contribution to the observed outcome of each individual's (latent) ability as a committing or disadvantaged agent. This would allow us, for example, to rank players from more to less effective, quantify uncertainty in this ranking and discover extra hierarchical structures involved in foul calls. All pretty useful stuff!
So how can we study this common and complex multi-agent interaction scenario, with hierarchical structures between more than a thousand individuals?
Despite the scenario's overwhelming complexity, Bayesian Item Response Theory combined with modern powerful statistical software allows for quite elegant and effective modeling options. One of these options employs a {term}Generalized Linear Model called Rasch model, which we now discuss in more detail.
We sourced our data from the official NBA Last Two Minutes Reports with game data between 2015 to 2021. In this dataset, each row k is one play involving two players (the committing and the disadvantaged) where a foul has been either called or not. So we model the probability p_k that a referee calls a foul in play k as a function of the players involved. Hence we define two latent variables for each player, namely:
theta: which estimates the player's ability to have a foul called when disadvantaged, andb: which estimates the player's ability to have a foul not called when committing.Note that the higher these player's parameters, the better the outcome for the player's team. These two parameters are then estimated using a standard Rasch model, by assuming the log-odds-ratio of p_k equals theta-b for the corresponding players involved in play k. Also, we place hierarchical hyperpriors on all theta's and all b's to account for shared abilities between players and largely different numbers of observations for different players.
Our analysis gives an estimate of the latent skills theta and b for each player in terms of posterior distributions. We analyze this outcome in three ways.
We first display the role of shared hyperpriors, by showing how posteriors of players with little observations are drawn to the league average.
Secondly, we rank the posteriors by their mean to view best and worst committing and disadvantaged players, and observe that several players still rank in the top 10 of the same model estimated in Austin Rochford blogpost on different data.
Thirdly, we show how we spot that grouping payers by their position is likely to be an informative extra hierarchical layer to introduce in our model, and leave this as an exercise for the interested reader. Let us conclude by mentioning that this opportunity of easily adding informed hierarchical structure to a model is one of the features that makes Bayesian modelling very flexible and powerful for quantifying uncertainty in scenarios where introducing (or discovering) problem-specific knowledge is crucial.
The analysis in this notebook is performed in four main steps:
We first import data from the original data set, which can be found at this URL. Each row corresponds to a play between the NBA seasons 2015-16 and 2020-21. We imported only five columns, namely
committing: the name of the committing player in the play.disadvantaged: the name of the disadvantaged player in the play.decision: the reviewed decision of the play, which can take four values, namely:
CNC: correct noncall, INC: incorrect noncall, IC: incorrect call, CC: correct call.committing_position: the position of the committing player which can take values
G: guard, F: forward, C: center, G-F, F-G, F-C, C-F.disadvantaged_position: the position of the disadvantaged player, with possible values as above.We note that we already removed from the original dataset the plays where less than two players are involved (for example travel calls or clock violations). Also, the original dataset does not contain information on the players' position, which we added ourselves.
We now process our data in three steps:
df by removing the position information from df_orig, and we create a dataframe df_position collecting all players with the respective position. (This last dataframe will not be used until the very end of the notebook.)df, called foul_called, that assigns 1 to a play if a foul was called, and 0 otherwise.Finally, we display the head of our main dataframe df along with some basic statistics.
We denote by:
We assume that each disadvantaged player is described by the latent variable:
and each committing player is described by the latent variable:
Then we model each observation $y_k$ as the result of an independent Bernoulli trial with probability $p_k$, where
for $k=1,2,...,K$, by defining (via a non-centered parametrisation)
\begin{align*} \theta_{i}&= \sigma_\theta\Delta_{\theta,i}+\mu_\theta\sim \text{Normal}(\mu_\theta,\sigma_\theta^2), &i=1,2,...,N_d,\ b_{j}&= \sigma_b\Delta_{b,j}\sim \text{Normal}(0,\sigma_b^2), &j=1,2,...,N_c, \end{align*}
with priors/hyperpriors
\begin{align*} \Delta_{\theta,i}&\sim \text{Normal}(0,1), &i=1,2,...,N_d,\ \Delta_{b,j}&\sim \text{Normal}(0,1), &j=1,2,...,N_c,\ \mu_\theta&\sim \text{Normal}(0,100),\ \sigma_\theta &\sim \text{HalfCauchy}(2.5),\ \sigma_b &\sim \text{HalfCauchy}(2.5). \end{align*}
Note that $p_k$ is always dependent on $\mu_\theta,,\sigma_\theta$ and $\sigma_b$ ("pooled priors") and also depends on the actual players involved in the play due to $\Delta_{\theta,i}$ and $\Delta_{b,j}$ ("unpooled priors"). This means our model features partial pooling. Morover, note that we do not pool $\theta$'s with $b$'s, hence assuming these skills are independent even for the same player. Also, note that we normalised the mean of $b_{j}$ to zero.
Finally, notice how we worked backwards from our data to construct this model. This is a very natural way to construct a model, allowing us to quickly see how each variable connects to others and their intuition. Meanwhile, when instantiating the model below, the construction goes in the opposite direction, i.e. starting from priors and moving up to the observations.
We now implement the model above in PyMC. Note that, to easily keep track of the players (as we have hundreds of them being both committing and disadvantaged), we make use of the coords argument for {class}pymc.Model. (For tutorials on this functionality see the notebook {ref}data_container or this blogpost.) We choose our priors to be the same as in Austin Rochford's post, to make the comparison consistent.
We now plot our model to show the hierarchical structure (and the non-centered parametrisation) on the variables theta and b.
We now sample from our Rasch model.
We plot below the energy difference of the obtained trace. Also, we assume our sampler has converged as it passed all automatic PyMC convergence checks.
Our first check is to plot
These plots show, as expected, that the hierarchical structure of our model tends to estimate posteriors towards the global mean for players with a low amount of observations.
As we successfully estimated the skills of disadvantaged (theta) and committing (b) players, we can finally check which players perform better and worse in our model.
So we now plot our posteriors using forest plots. We plot the 10 top and bottom players ranked with respect to the latent skill theta and b, respectively.
By visiting Austin Rochford post and checking the analogous table for the Rasch model there (which uses data from the 2016-17 season), the reader can see that several top players in both skills are still in the top 10 with our larger data set (covering seasons 2015-16 to 2020-21).
A natural question to ask is whether players skilled as disadvantaged players (i.e. players with high theta) are also likely to be skilled as committing players (i.e. with high b), and the other way around. So, the next two plots show the theta (resp. b) score for the top players with respect to b ( resp.theta).
These plots suggest that scoring high in theta does not correlate with high or low scores in b. Moreover, with a little knowledge of NBA basketball, one can visually note that a higher score in b is expected from players playing center or forward rather than guards or point guards.
Given the last observation, we decide to plot a histogram for the occurrence of different positions for top disadvantaged (theta) and committing (b) players. Interestingly, we see below that the largest share of best disadvantaged players are guards, meanwhile, the largest share of best committing players are centers (and at the same time a very small share of guards).
The histograms above suggest that it might be appropriate to add a hierarchical layer to our model. Namely, group disadvantaged and committing players by the respective positions to account for the role of position in evaluating the latent skills theta and b. This can be done in our Rasch model by imposing mean and variance hyperpriors for each player grouped by the positions, which is left as an exercise for the reader. To this end, notice that the dataframe df_orig is set up precisely to add this hierarchical structure. Have fun!
A warm thank you goes to Eric Ma for many useful comments that improved this notebook.
:::{bibliography} :filter: docname in docnames :::
:::{include} ../page_footer.md :::