# Data manipulation
import pandas as pd
import numpy as np
# We'll need to solve a system
from scipy.linalg import solve
# Graphics
import plotly.express as px
import plotly.graph_objects as go
import graphviz as gvMassey ratings
In this notebook, we’re going to outline a technique based in linear algebra to rate teams playing games in a common league based on the score differences that we see in the games. The technique dates back to Ken Massey’s 1997 senior honors thesis at Bluefield College. My understanding of the technique is based mostly on Chapter 2 of Langville and Meyer’s text on sport’s ranking.
Imports
First, here are a few of the libraries we’ll use:
The general problem
In general, we’ve got \(m\) games played between \(n\) teams. Let’s suppose that in the \(k^{\text{th}}\) game, team \(i\) beat team \(j\) by \(y_k\) points. We’d like the teams ratings \(r_i\) and \(r_j\) to reflect this fact, i.e. we’d like \[r_i - r_j = y_k.\] To put it another way, we’d like our ratings of teams \(i\) and \(j\) to indicate the expected score difference between the teams, if they play again.
If \(\vec{r}\) is the vector of ratings and \(\vec{y}\) is the vector of point differences, then we can express this as the matrix equation \[X\vec{r} = \vec{y},\] where \(X_{ki} = 1\), \(X_{kj} = -1\), and all other entries are zero.
Now, there are generally many more games than teams so that the system we’ve written down is likely overdetermined. Nonetheless, we can attempt to find an approximate solution in the least squares sense and an efficient way to do so is to solve the so-called normal equations obtained by multiplying both sides by the transpose of \(X\): \[X^T X\vec{r} = X^T \vec{y},\]
The matrix \(M = X^T X\) formed in this process is called the Massey matrix and has some interesting properties. Note that each entry in column \(i\) of the matrix \(X\) is either \(1\) (if team \(i\) won that game), \(-1\) (if team \(i\) lost that game), or zero (if team \(i\) didn’t play in that game). As a result, \(M_{ii}\) (which is obtained by computing the dot product of the \(i^{\text{th}}\) column of \(X\) with itself) denotes the number of games played by team \(i\). Similarly, \(M_{ij}\) is the opposite of the number of games played by teams \(i\) and \(j\). Thus, the sum of each row is zero. Put another way, the columns sum to the zero vector and are therefore linearly dependent. This implies that \(M\) is singular so that the normal equations don’t have a unique solution.
A way to get around this problem is to set all the elements in some row of \(M\) equal to 1. This is conventionally done in the last row. We similarly replace the last element of \[\vec{p} = X^T y\] to zero. This implies the constraint that the sum of the ratings must be zero. Furthermore, the resulting system \[M \vec{r} = \vec{p}\] can be shown to have a unique solution.
A simple working example
As a simple example, let’s suppose that we have three teams labeled \(A\), \(B\), and \(C\) and that each team plays the other teams twice each. We might visualize the collection of games via a labeled, directed graph like so:
Code
# Code to generate digraph
dot = gv.Digraph()
dot.graph_attr.update({'rankdir': 'LR'})
dot.edge("A", "B", label="2")
dot.edge("A", "B", label="3")
dot.edge("B", "C", label="1")
dot.edge("C", "B", label="2")
dot.edge("A", "C", label="3")
dot.edge("C", "A", label="1")
dot
Note that a directed edge from team \(i\) to team \(j\) indicates that team \(i\) defeated team \(j\) in that game; the edge labels indicates the score differences. Given the directed graph shown above, how would you rank the teams?
To rank the teams the Massey way, we start with the general equation \[X \vec{r} = \vec{y}\] which takes the specific form
\[ \begin{pmatrix}1&-1&0\\1&-1&0\\0&1&-1\\0&-1&1\\1&0&-1\\-1&0&1\end{pmatrix} \begin{pmatrix}r_0\\r_1\\r_2\end{pmatrix} = \begin{pmatrix}3\\2\\1\\2\\3\\1\end{pmatrix}. \]
If we multiply both sides by \(X^T\) and modify the last rows according Massey’s strategy, the equation becomes
\[ \begin{pmatrix}4&-2&-2\\-2&4&-2\\1&1&1\end{pmatrix} \begin{pmatrix}r_0\\r_1\\r_2\end{pmatrix} = \begin{pmatrix}7\\-6\\0\end{pmatrix}. \]
The solution of this system is
\[ \begin{pmatrix}r_0\\r_1\\r_2\end{pmatrix} = \begin{pmatrix}7/6\\-1\\-1/6\end{pmatrix}. \]
In particular, note that \(r_1 < r_2 < r_0\), which is probably in accordance with the order you expected.
The Massey equation
It’s worth noting that the Massey equation \[M\vec{r} = \vec{p}\] can itself be written down and interpreted without reference to \(X\) or \(y\). As we’ve already noted, \(M_{ii}\) denotes the number of games played by team \(i\) and \(M_{ij}\) is the opposite of the number of games played by teams \(i\) and \(j\). Similarly, the entries of \(\vec{p}\) is the accumulated set of score differences for each team.
For example, if we expand the first term in the product \(M\vec{r}\) and set that equal to the first term of \(p\), we get \[4r_0 - 2r_1 - 2r_2 = 7.\] This indicates that team \(0\) played four games, two against each other team, and that it outscored its opponents by a combined score of seven. If we solve for \(r_0\), we get \[r_0 = \frac{1}{2}r_1 + \frac{1}{2}r_2 + \frac{7}{2}.\] In general, the rating of team \(i\) is a combination of the ratings of the other teams together with the number of points that team \(i\) outscored those other teams by.
Implementation
Here’s an implementation of the technique for this example in Python:
X = np.array([
[1, -1, 0],
[1, -1, 0],
[0, 1, -1],
[0, -1, 1],
[1, 0, -1],
[-1, 0, 1]
])
M = X.transpose().dot(X)
M[-1] = np.ones(len(M))
y = np.array([[3], [2], [1], [2], [3], [1]])
y = np.array([3, 2, 1, 2, 3, 1])
y = X.transpose().dot(y)
y[-1]=0
solve(M,y)array([ 1.16666667, -1. , -0.16666667])This is exactly the [7/6,-1,-1/6] ratings that we mentioned above.
The whole NCAA
Let’s take a look at a somewhat bigger example - namely, the 2023/24 Men’s NCAA. We’ll use our CSV file that lists the results of all \(m=5607\) games:
games = pd.read_csv('https://marksmath.org/data/CombinedCondensedNCAAData2024.csv')
m = len(games)
games| WTeamName | WIDX | WScore | LTeamName | LIDX | LScore | |
|---|---|---|---|---|---|---|
| 0 | Abilene Chr | 0 | 64 | Oklahoma St | 221 | 59 |
| 1 | Akron | 2 | 81 | S Dakota St | 249 | 75 |
| 2 | Alabama | 3 | 105 | Morehead St | 183 | 73 |
| 3 | Arizona | 9 | 122 | Morgan St | 184 | 59 |
| 4 | Ark Little Rock | 11 | 71 | Texas St | 302 | 66 |
| ... | ... | ... | ... | ... | ... | ... |
| 5602 | Auburn | 16 | 86 | Florida | 89 | 67 |
| 5603 | Duquesne | 76 | 57 | VCU | 332 | 51 |
| 5604 | Illinois | 120 | 93 | Wisconsin | 355 | 87 |
| 5605 | UAB | 309 | 85 | Temple | 296 | 69 |
| 5606 | Yale | 360 | 62 | Brown | 30 | 61 |
5607 rows × 6 columns
We can extract and sort all \(n=362\) teams:
teams = games.LTeamName.drop_duplicates().sort_values()
n = len(teams)
teams209 Abilene Chr
70 Air Force
700 Akron
872 Alabama
166 Alabama A&M
...
182 Wright St
455 Wyoming
339 Xavier
189 Yale
96 Youngstown St
Name: LTeamName, Length: 362, dtype: objectWe can now programmatically set up and solve the Massey system of equations and simultaneously set up the directed graph visualization
X = np.zeros([m,n])
y = np.zeros([m])
X.transpose().dot(y)
for i,r in games.iterrows():
X[i, r.WIDX] = 1
X[i, r.LIDX] = -1
y[i] = r.WScore - r.LScore
M = X.transpose().dot(X)
M[-1,:] = 1
p = X.transpose().dot(y)
p[-1] = 0
r = solve(M,p)
rarray([-3.44800331e+00, -3.96214497e+00, 3.38656310e+00, 2.14693447e+01,
-1.42910781e+01, -1.13623436e+01, -1.12216367e+01, -8.41941895e+00,
6.80294991e+00, 2.57083745e+01, 5.24278134e+00, -2.58561054e+00,
-1.52541184e+01, 6.73854210e+00, 8.98091753e-01, -1.40739773e+01,
2.36914293e+01, -4.87069887e+00, 2.07265696e+01, -6.10134029e+00,
2.05852606e+01, -1.10502701e+01, 3.29439289e+00, -1.17983432e+01,
-8.73830067e+00, 1.41974194e+01, 9.84440618e+00, -9.29918078e+00,
-4.42542558e+00, 1.00494549e+01, -3.00410689e+00, -2.47534474e+00,
-7.73284693e+00, -1.45032474e+01, 1.25921759e+01, -7.46411131e+00,
-4.68926078e+00, -4.37680648e+00, -4.07560987e+00, -1.00058338e+01,
-2.00938354e+00, -1.28893322e+01, 5.60270582e+00, -1.00677801e+01,
-7.35547550e+00, -1.49595593e+01, -7.65488270e+00, -1.15039052e+01,
4.74526096e+00, 3.38492217e-01, -9.34786582e+00, 1.59359853e+01,
-6.05955553e+00, 1.52859383e+01, -2.53364284e+00, -9.81570182e+00,
4.57246121e+00, 4.40198397e-01, 1.63838330e+01, 1.38369732e+01,
-4.86438688e+00, 2.50806036e+01, -2.10651834e+01, 3.99942924e+00,
2.03892159e+01, -1.35171367e+01, 4.21504051e+00, 1.39266964e+01,
-5.97462178e+00, -1.10401382e+00, -1.01461548e+01, -6.64945493e+00,
-1.63369083e+01, 1.15668621e+01, 1.77265747e+00, 2.06052334e+01,
7.90444234e+00, -1.09139969e+01, -4.73628257e+00, -1.21555519e+01,
-2.14799073e-01, -2.56135514e+00, 1.77915989e-02, -1.08558324e+01,
-3.01258172e+00, -1.38894304e+01, 1.42594388e+01, -5.59816325e+00,
-1.18527859e+00, 1.58457369e+01, -1.53997119e+01, -7.39680565e+00,
9.79552217e+00, 9.50693538e-02, -2.46108931e+00, 7.39252413e-01,
-4.36071216e-01, -7.77752420e+00, -4.04444378e+00, 6.99314722e+00,
9.08806681e-01, 8.87254150e+00, -3.88861880e+00, 4.91946966e+00,
1.76755906e+01, -9.79992302e+00, 1.03905517e+01, -1.52639947e+01,
-3.98347977e+00, 6.60909876e-01, 2.60347272e+00, 2.82502134e+00,
-1.65466412e+01, 2.67197310e+01, -1.85153339e+01, -8.11206766e+00,
-1.07251920e+00, -2.07227902e+01, -1.04994975e+01, -6.66309719e+00,
2.02341314e+01, -8.51480174e-01, -1.56428661e+01, 8.39766703e+00,
1.31361947e+01, -4.51824745e+00, 1.36061399e+01, 2.33809631e+01,
-1.11186211e+01, -9.14587776e+00, -2.93095638e+00, 9.04260634e+00,
1.79889605e+01, 1.21862953e+01, -7.32859215e+00, -1.41102276e-01,
1.79789935e+01, -1.73621287e+01, 8.99973068e+00, -1.04744110e+00,
-1.05685169e+01, -4.66772164e+00, -1.05090015e+01, -6.89962874e+00,
2.14909994e+00, -1.80539637e+01, -1.23977872e+00, -6.11130563e-01,
-1.11334847e+00, 5.11872725e-02, 5.42179262e+00, 7.46434554e-01,
4.47070600e-01, -1.56152412e+01, 7.44229605e+00, -7.28711794e-01,
-1.69088502e+01, -2.57887432e+01, -6.86084987e+00, -5.79726018e+00,
-1.46791709e+01, -4.96568434e+00, 1.99378706e+01, -4.74613756e+00,
1.11955343e+01, 8.17970620e+00, 7.52481719e+00, 9.59672987e+00,
-4.13330305e+00, -4.73858803e+00, 9.90786815e+00, -5.83447839e+00,
5.94779377e+00, 1.82584036e+01, 1.02687717e+01, 8.89349601e+00,
1.52965957e+01, 3.86364019e+00, -4.60926816e+00, 1.49465006e+00,
-2.62210839e+00, -1.02333340e-01, -4.20511903e+00, 1.06389081e+00,
-1.55381765e+01, -7.14712630e+00, 9.29439229e-01, -3.62531210e+00,
-6.87285570e+00, -8.59293184e+00, -3.02233175e+00, -1.57797100e+01,
-7.36400755e+00, 1.19353410e+01, -6.61100850e+00, -1.41073338e+01,
-1.10903185e+01, 1.54899393e+01, 1.37367099e+01, -5.61452483e+00,
1.63804168e+01, -7.95141190e+00, -1.59675214e+01, -7.76576394e+00,
-6.28361672e+00, -5.92006069e+00, -6.40271463e+00, 2.13593658e+01,
-5.54168274e+00, -6.26262161e+00, 9.17376335e+00, -5.54516942e+00,
-1.14868286e+01, 4.33394875e+00, 1.31122677e+01, -1.35994927e+01,
5.31026384e+00, 1.97526704e+00, 1.10877939e+00, 1.44682600e+01,
1.55261109e+01, 7.16374100e+00, -7.65790547e+00, -6.71980671e+00,
1.26187086e+01, 3.90904582e+00, 1.84342498e-01, -1.68277987e+01,
-2.69511655e+00, 9.80645693e+00, -1.56289871e+00, 1.51662436e+01,
-7.64390328e+00, -6.32216802e+00, -1.50248175e+01, -7.65572629e+00,
7.30964473e+00, 1.30739963e+01, 2.51384479e+01, -8.06538639e+00,
-2.38524163e+00, -5.11688662e+00, -1.64858929e+00, -2.56528347e+00,
8.62527382e+00, -5.06360496e+00, -8.95604905e+00, 8.03493182e+00,
-1.13955754e+01, 1.08049633e+00, 4.62294972e+00, -8.74326004e+00,
-1.04907897e+01, -1.76304709e+01, -5.69326991e-01, -8.43937623e+00,
1.10534199e+01, -6.09235464e+00, -9.38300335e+00, 2.27732518e-01,
6.62553229e+00, -4.71413666e+00, 1.56762576e+01, 1.00602936e+01,
-2.99064385e+00, 5.98084074e+00, 3.71971725e+00, 1.17277194e+01,
-1.85461423e+01, -5.26942831e+00, 1.19178405e+01, -1.18518312e+01,
7.75188067e+00, -1.30711958e+01, -5.26000529e+00, -8.29321133e+00,
-5.34340628e+00, 7.42610434e+00, -1.81313244e+01, 1.66547620e+01,
7.14051326e+00, -1.01068579e+00, 1.73228523e+01, -3.70202708e+00,
-2.54041878e-01, 8.01207895e+00, -4.13288791e+00, -1.96267839e+01,
-2.05563143e+00, 8.68491758e+00, -3.29501388e+00, 1.60295228e+01,
-4.96844622e+00, -1.36508625e+01, -8.70917877e+00, 1.15639452e+00,
-6.23052534e-01, 2.21523176e+01, -9.03731368e+00, -1.44162764e+01,
1.64152601e+01, 1.34014902e+01, -3.02461327e+00, 1.63455811e+01,
2.37493864e+00, 7.67192499e-01, 1.86819981e+00, 4.33066346e+00,
-5.53670616e-01, 5.57112903e+00, -5.62647492e-01, 7.67408872e+00,
-3.56210789e+00, 3.17310957e+00, -2.43359577e+00, 1.31621322e+01,
8.09464096e+00, -9.35895107e+00, -8.11401334e+00, -1.29470035e+00,
2.74028560e-01, 1.81761988e+00, 9.28639394e+00, 9.96982957e+00,
2.76608270e+00, -4.59579431e+00, -9.17425404e-01, -1.29788373e+01,
1.53572858e+01, 1.25310972e+01, -7.71666611e+00, -9.80123387e-01,
9.14598079e+00, -2.00663196e+01, -7.78886717e+00, 8.17922055e-01,
3.29791704e+00, 1.51526116e+01, 1.13470274e+01, 1.21537913e+01,
2.51265545e+00, -8.09476923e+00, -8.68951090e+00, -4.98787839e+00,
-4.40408292e+00, 1.80359123e+00, -1.04313911e+01, 1.55234516e+01,
1.24477313e+01, 1.46496609e+01, 3.46168443e-01, 4.61321056e+00,
2.69079300e+00, -1.07568176e+01, -1.42021929e+00, 1.71678392e+01,
-3.34000965e+00, 3.78057190e-01, 1.33258493e+00, 1.26120476e+01,
5.21484520e+00, 1.57686732e+00])That won’t mean much, unless we associate the ratings with the teams:
ratings = pd.DataFrame()
ratings['team'] = teams
ratings['r'] = r
ratings = ratings.sort_values('r', ascending=False)
ratings.index = range(1,len(teams)+1)
ratings[:25]| team | r | |
|---|---|---|
| 1 | Houston | 26.719731 |
| 2 | Arizona | 25.708375 |
| 3 | Purdue | 25.138448 |
| 4 | Connecticut | 25.080604 |
| 5 | Auburn | 23.691429 |
| 6 | Iowa St | 23.380963 |
| 7 | Tennessee | 22.152318 |
| 8 | Alabama | 21.469345 |
| 9 | North Carolina | 21.359366 |
| 10 | BYU | 20.726570 |
| 11 | Duke | 20.605233 |
| 12 | Baylor | 20.585261 |
| 13 | Creighton | 20.389216 |
| 14 | Illinois | 20.234131 |
| 15 | Marquette | 19.937871 |
| 16 | Michigan St | 18.258404 |
| 17 | Kansas | 17.988960 |
| 18 | Kentucky | 17.978994 |
| 19 | Gonzaga | 17.675591 |
| 20 | St Mary's CA | 17.322852 |
| 21 | Wisconsin | 17.167839 |
| 22 | St John's | 16.654762 |
| 23 | Texas | 16.415260 |
| 24 | Colorado | 16.383833 |
| 25 | New Mexico | 16.380417 |
Looks good!