chrisberry888.github.io
<!DOCTYPE html>
Expected Value, Win Probability, and why "Common Knowledge" Hurts Sports Teams
Tutorial by Chris Berry
The four major American sports leagues (MLB, NFL, NBA, and NHL) have very long and extensive histories, and established throughout each league's history were traditions and "common knowledge" that many managers and front offices adhere to, to this day. However, many of the rules and practices that teams follow simply do not match up with what they mathematically should be doing.
The goal of every sports team should be to win games, and so every decision made by a coach should be made with the intention of increasing the probability that that team will win the game. For example, if a football coach thinks that kicking an onside kick will make them more likely to win the game, then they should choose to kick an onside kick. Despite this, there are many decisions that coaches make that are made because of deep-rooted traditions in the sport, without considering the data whatsoever. One of the most famous examples of this is punting in football: conventional wisdom states that football teams should almost never go for it on 4th down, but many different sources show that punting is almost never the mathematically correct decision (https://www.footballstudyhall.com/2013/11/15/5105958/fourth-down-pulaski-academy-kevin-kelley#:~:text=Kevin%20Kelley%2C%20Pulaski's%20coach%20never,attempts%20onside%20kicks%20after%20scoring).
In this tutorial I will discuss the ways in which win probability is used and misused in when deciding to pull the goalie in the NHL. In doing so, I will walk through all the steps of analyzing data, from retreiving, cleaning, and processing the data to creating predictive models of different real-life scenareos.
Pulling the Goalie
Before we talk about win probability, we need to talk about a concept called "expected value". In the context of sports, expected value (EV) is a measure of how many points that you can expect will result from a decision. To explain EV, here's an example: say we are playing a game, where whenever I flip a coin, if it lands on heads, I give you 10 dollars, and if it lands on tails, you give me 5 dollars. The equation for EV in this situation is EV = .5(10 dollars) + .5(-5 dollars) = 2.5 dollars; therefore, you can "expect" to win 2 dollars and 50 cents every time I flip the coin.
You might expect that making the highest EV decisions would always result in a higher chance for winning a game, but this is not always the case. Imagine a hockey game between the Boston Bruins and the New York Rangers: there are 2 minutes left in the game, and the Bruins are losing 2-3. With two minutes left, the Bruins decide to substitute their goalie for another skater. Andrew Thomas of Harvard University analyzed data from all NHL games between 2003 and 2007 and found that, in games where a goalie was pulled, the winning team scored on the empty net 34% of the time, while the team with the extra skater scored 30% of the time (https://fliphtml5.com/snvq/xptx/basic). So if the Rangers and Bruins are evenly matched, then the Bruins' expected value of pulling the goalie is EV = .3 (1) + .34 (-1) = -.04.
If the Bruins get negative EV from pulling the goalie, then why do they do it? The answer is because they don't care whether they lose by 1 or lose by 2: they just want to tie the game. Despite being a negative EV decision, deciding to pull the goalie increases the Bruins' chances of winning the game. It would be terribly stupid to pull the goalie while the score is still 0-0 in the first period, but it becomes the smart thing to do if you're losing in the third. This begs the question: when should a team decide to pull their goalie?
In this part, we will be using NHL Game Data courtesy of Martin Ellis (https://www.kaggle.com/martinellis/nhl-game-data?select=game_plays.csv). This database contains game data from almost all NHL games played between the 2000-01 season and the 2019-20 season. We will be extracting data from 4 of the tables in this database: "total_shifts" tells us how much time each player is on the ice; "total_plays" tells us about every single play that happens in each game, be it a shot, penalty, goal, faceoff, or anything else; "games" gives us general stats about each game; and "game_goalie_stats" gives us the stats for each goaltender for every game. The ultimate goal is to learn more about when the optimal time to pull the goalie is in a losing game, based on available data.
The following are all of the python packages that we will be using. They will be used later on for manipulating data and creating models, among other things
import requests
from bs4 import BeautifulSoup
import pandas as pd
import datetime
import numpy as np
import re
from matplotlib import pyplot as plt
import sqlite3
import statistics
from sklearn.linear_model import LinearRegression
from sklearn.preprocessing import PolynomialFeatures
from sklearn.neural_network import MLPRegressor
from sklearn.linear_model import LogisticRegression
from sklearn.svm import SVR
import warnings
import random
warnings.filterwarnings("ignore")
pref = "D:\\Documents\\Computer Programs\\CMSC320\\Final_Project\\"
We start by importing the data into Python. The "drop_duplicates()" command ensures that no game, goal, or player is being double-counted.
total_shifts = pd.read_csv(pref + "NHL\\game_shifts.csv")
total_shifts = total_shifts.drop_duplicates()
total_plays = pd.read_csv(pref + "NHL\\game_plays.csv")
total_plays = total_plays.drop_duplicates()
games = pd.read_csv(pref + "NHL\\game.csv")
games = games.drop_duplicates(subset=["game_id"])
game_goalie_stats = pd.read_csv(pref + "NHL\\game_goalie_stats.csv")
game_goalie_stats = game_goalie_stats.drop_duplicates()
Next, we will find every game in which a team was losing at the end of the game and decided to pull the goalie. Unfortunately, this information is not readily available in the database, so we will have to extrapolate from what we have. So in this step, we will first isolate every goalie's shift. Goalies are the only players in hockey who stay on the ice for more than just a minute or so at a time. All goalies start their shift at the beginning of the third period (2400 seconds into the game), and given that the longest shift by a skater in NHL history is around 310 seconds, we can safely say that if a shift starts at time marker 2400 and ends beyond 2710, then it's a goalie's shift.
We also know that if a goalie finished the game in net (i.e. their shift ends at time marker 3600), then they were not pulled. Therefore, any goalie whose shift ends before the 3600 mark must have been pulled. We store these shifts in the "pulled_goalies" variable.
goalies = total_shifts[ (total_shifts["shift_start"] == 2400) & (total_shifts["shift_end"] > 2710) ]
pulled_goalies = goalies[goalies["shift_end"] < 3600]
There are almost always only 2 goalies playing per game, which means that by dividing the total number of third-period goalie shifts by two, we will get the total number of games in this data set. Dividing the number of games with a pulled goalie by the total number of games yields the proportion of games that have a pulled goalie:
number_of_games = goalies.shape[0] / 2
print(pulled_goalies.shape[0] / number_of_games)
So about 72.93% of games have a situation where the goalie is pulled at the end of the game. While this is a highly accurate figure, it is not completely accurate. The main culprit of this inaccuracy has to do with delayed penalties: when a team gets flagged for a penalty, the other team's goalie will go to the bench as a way of starting their power play early. This means that the 72.9% number is an overestimate, since it counts all delayed penalties like that in its calculation. Despite this inaccuracy, we will continue to use this table since it is accurate enough.
pulled_goalies.head()
This is a small section of the "pulled_goalies" table. From now on, we will use this table to represent the games that have a pulled goalie in them. This table doesn't tell us very much about each game, so we'll use the dataset to add more and more useful information. First, we'll add the season that each game was played in, as well as the time left in the game when the goalie was pulled:
seasons = []
for index, row in pulled_goalies.iterrows():
seasons.append(int(games.loc[games["game_id"] == row["game_id"], "season"].drop_duplicates()))
pulled_goalies["season"] = seasons
pulled_goalies["time_left"] = 3600 - pulled_goalies["shift_end"]
pulled_goalies.head()
Next, we'll figure out what teams played in this game, what the score was when the goalie was pulled, and what happened afterwards. To do this, we'll start by creating another table containing every goal, as well as what game the goal belongs to:
goals = total_plays[total_plays["event"] == "Goal"]
goals["gameTime"] = goals["periodTime"] + (goals["period"] - 1) * 1200
This next chunk of code looks intimidating, but it's simpler than it looks. We want a few bits of information from the "goals" table: for each game we want to figure out who's winning, who's losing, what's the score, how much is the losing team losing by, how many goals did the losing team get after pulling the goalie, and did the winning team score at all.
losing_scores = []
losing_teams = []
winning_scores = []
winning_teams = []
goals_for_losing_team = []
did_winning_team_score = []
for index, row in pulled_goalies.iterrows():
game = row["game_id"]
time = row["shift_end"]
home_team = games[games["game_id"] == game]
home_team = int(home_team["home_team_id"])
away_team = games[games["game_id"] == game]
away_team = int(away_team["away_team_id"])
p3_goals = goals[goals["game_id"] == game]
post_hook_goals = p3_goals[(p3_goals["gameTime"] > time)&(p3_goals["period"] == 3)]
p3_goals = p3_goals[p3_goals["gameTime"] < time]
if p3_goals.shape[0] > 0:
goal_before = p3_goals[p3_goals["gameTime"] == p3_goals["gameTime"].max()]
a = int(goal_before["goals_away"])
h = int(goal_before["goals_home"])
losing_score = ([a, away_team] if a < h else [h, home_team])
winning_score = ([h, home_team] if a < h else [a, away_team])
diff = winning_score[0] - losing_score[0]
phg = 0 #goals scored by losing team before 1) the period ends, 2) the other team scores, or 3) they even the score
dwts = 0
if(len(post_hook_goals) > 0):
g1 = post_hook_goals.iloc[0]
if g1.loc["team_id_for"] == losing_score[1]:
phg += 1
if diff >= 2 and len(post_hook_goals) >= 2:
g2 = post_hook_goals.iloc[1]
if g2.loc["team_id_for"] == losing_score[1]:
phg += 1
if diff >= 3 and len(post_hook_goals) >= 3:
g3 = post_hook_goals.iloc[2]
if g3.loc["team_id_for"] == losing_score[1]:
phg += 1
else:
dwts = 1
else:
dwts = 1
else:
dwts = 1
losing_scores.append(losing_score[0])
losing_teams.append(losing_score[1])
winning_scores.append(winning_score[0])
winning_teams.append(winning_score[1])
goals_for_losing_team.append(phg)
did_winning_team_score.append(dwts)
else:
losing_scores.append(0)
losing_teams.append(home_team)
winning_scores.append(0)
winning_teams.append(away_team)
goals_for_losing_team.append(0)
did_winning_team_score.append(0)
pulled_goalies["losingScore"] = losing_scores
pulled_goalies["losingTeam"] = losing_teams
pulled_goalies["winningScore"] = winning_scores
pulled_goalies["winningTeam"] = winning_teams
pulled_goalies["goalsForLosingTeam"] = goals_for_losing_team
pulled_goalies["didWinningTeamScore"] = did_winning_team_score
pulled_goalies = pulled_goalies[pulled_goalies["losingScore"] != pulled_goalies["winningScore"]]
pulled_goalies["goal_diff"] = pulled_goalies["winningScore"] - pulled_goalies["losingScore"]
pulled_goalies.head()
At this point, we have enough data to do some simple analysis. But what other factors might influence the losing team's ability to score? Well, if the losing team is usually pretty good on offense, then clearly they will be more likely to come back. Likewise, if the winning team is usually pretty good on defense, then a comeback will be less likely. Because of this, we should calculate how good each losing team's offense is, as well as how good each winning team's defense is. How can we quantify this? One way of measuring a team's offense is by using their Goals Per Game (GPG), which is the average number of goals that team scores in an average game. A way of quantifying a team's defense is by using their goalie's Goals Against Average (GAA), which is the average number of goals that they allow in an average game.
seasons = []
for index, row in game_goalie_stats.iterrows():
seasons.append(int(games.loc[games["game_id"] == row["game_id"], "season"].drop_duplicates()))
game_goalie_stats["season"] = seasons
goalies = game_goalie_stats[["team_id", "season"]].drop_duplicates()
game_goalie_stats["goalsAgainst"] = game_goalie_stats["shots"] - game_goalie_stats["saves"]
GAA = []
for index, row in goalies.iterrows():
total_GA = game_goalie_stats.loc[(game_goalie_stats["team_id"] == row["team_id"])&(game_goalie_stats["season"] == row["season"]), "goalsAgainst"].sum()
TOI = game_goalie_stats.loc[(game_goalie_stats["team_id"] == row["team_id"])&(game_goalie_stats["season"] == row["season"]), "timeOnIce"].sum()
GAA.append(total_GA/(TOI/3600))
teams = games[["away_team_id", "season"]].drop_duplicates()
GPG = []
for index, row in teams.iterrows():
total_GA = games.loc[(games["away_team_id"] == row["away_team_id"])&(games["season"] == row["season"]), "away_goals"].sum()
total_GA += games.loc[(games["home_team_id"] == row["away_team_id"])&(games["season"] == row["season"]), "home_goals"].sum()
total_games = int(games.loc[(games["away_team_id"] == row["away_team_id"])&(games["season"] == row["season"])].shape[0]*2)
GPG.append(total_GA/total_games)
goalies["GAA"] = GAA
teams["GPG"] = GPG
We then add this information to the main table:
losing_GPG = []
winning_GAA = []
for index, row in pulled_goalies.iterrows():
losing_GPG.append(float(teams.loc[(teams["away_team_id"] == row["losingTeam"])&(teams["season"] == row["season"]), "GPG"]))
winning_GAA.append(float(goalies.loc[(goalies["team_id"] == row["winningTeam"])&(goalies["season"] == row["season"]), "GAA"]))
pulled_goalies["losingGPG"] = losing_GPG
pulled_goalies["winningGAA"] = winning_GAA
pulled_goalies.head()
We finally have all the information we need. In our analysis, we'll define a few measures of success. The first is simply whether the losing team came back and won the game: 1 represents a win, 0 represents a loss.
did_they_win = []
for index, row in pulled_goalies.iterrows():
succ = 0
game = games.loc[games["game_id"] == row["game_id"]]
if (row["losingTeam"] == game["away_team_id"]).bool() and ("away win" in str(game["outcome"])):
succ = 1
elif (row["losingTeam"] == game["home_team_id"]).bool() and "home win" in game["outcome"]:
succ = 1
did_they_win.append(succ)
pulled_goalies["didTheyWin"] = did_they_win
The next is whether the losing team overcame their defecit, regardless of whether they won or lost: 1 represents if they did, 0 represents if they didn't.
pulled_goalies["didTheyComeBack"] = np.where(pulled_goalies["goal_diff"] - pulled_goalies["goalsForLosingTeam"] == 0, 1, 0)
The third measure of success will be net goals, i.e. the goals they scored with the net empty minus the goals their opponent scored.
pulled_goalies["net_goals"] = pulled_goalies["goalsForLosingTeam"] - pulled_goalies["didWinningTeamScore"]
Before we construct any models of this data, we will first make some general observations. First, one might wonder: how often are teams able to eliminate their defecit when they pull the goalie?
temp = pulled_goalies
print(temp.shape[0]) #sample size
print(temp.loc[temp["didTheyComeBack"] == 1].shape[0]) # number of games where they came back
print(temp.loc[temp["didTheyComeBack"] == 0].shape[0]) # number of games where they didn't come back
print(temp.loc[temp["didTheyComeBack"] == 1].shape[0]/temp.shape[0]) # proportion of games where they came back
According to our data, pulling the goalie works about 9.3% of the time. Now, assuming parity in the league, one would expect teams to win even-scored games about equally as often as they lose them, so one would expect the losing team to win about 4.65% of the time:
print(temp.loc[temp["didTheyWin"] == 1].shape[0]/temp.shape[0])
However, this is not the case; the losing teams win a measley 2.3% of the time when they pull the goalie. Why is this, what is causing this discrepency? There are a few explanations, and the first is with the data itself. The database is not 100% accurate, and as stated at the beginning of this tutorial, there are incorrect data points representing delayed penalties and mid-period goalie changes. These inaccurace data points will make our findings less accurate, but it certainly doesn't explain the above discrepency. Another explanation is that there isn't parity in these games: if a team gets out to a big lead at the beginning of the game, then there's a good chance that they're simply a better team.
Next we'll look at how success changes as the defecit changes:
temp1 = pulled_goalies[(pulled_goalies["goal_diff"] == 1)]
temp2 = pulled_goalies[(pulled_goalies["goal_diff"] == 2)]
temp3 = pulled_goalies[(pulled_goalies["goal_diff"] == 3)]
print(temp1.loc[temp["didTheyComeBack"] == 1].shape[0]/temp1.shape[0])
print(temp2.loc[temp["didTheyComeBack"] == 1].shape[0]/temp2.shape[0])
print(temp3.loc[temp["didTheyComeBack"] == 1].shape[0]/temp3.shape[0])
This indicates that pulling the goalie will help eliminate a 1, 2, and 3 goal defecit 15.79%, 2.04%, and .93% of the time, respectively. This tells us the obvious: the larger your opponent's lead, the harder it will be to claw your way back. In fact, only once in the past 20 years has a team come back to win the game after pulling the goalie with a 3-0 defecit.
Next we'll look at the times that successful and unsuccessful teams decide to take the goalie out:
temp1 = pulled_goalies
temp2 = pulled_goalies[pulled_goalies["didTheyComeBack"] == 1]
temp3 = pulled_goalies[pulled_goalies["didTheyComeBack"] == 0]
q1 = temp2["time_left"].quantile(.25)
q3 = temp2["time_left"].quantile(.75)
mask = temp2["time_left"].between(q1, q3, inclusive=True)
iqr_success = temp2.loc[mask, "time_left"]
q1 = temp3["time_left"].quantile(.25)
q3 = temp3["time_left"].quantile(.75)
mask = temp3["time_left"].between(q1, q3, inclusive=True)
iqr_fail = temp3.loc[mask, "time_left"]
print(temp1.loc[temp1["didTheyComeBack"] == 1, "time_left"].median())
print(temp1.loc[temp1["didTheyComeBack"] == 0, "time_left"].median())
print()
print(iqr_success.mean())
print(iqr_fail.mean())
The above block calculated the median time and interquartile mean time that losing teams were successful/unsuccessful when pulling the goalie: the median successful time is 121 seconds left, while the median unsuccessful time is 116 seconds. Meanwhile, the IQM successful time is a full 34 seconds sooner on average than the IQM unsuccessful time. This gives the indication that teams should pull their goalies sooner than they do now.
NOTE: I decided to use the interquartile mean instead of the regular mean in order to account for the very large and very small outliers in the dataset; the IQM is insensitive to outliers, unlike the regular mean.
Lastly, we'll attempt to model the relationship between the different independent variables with success rate. We'll start by creating a neural network that takes 3 inputs (time left in the game, the losing team's GPG, and the winning team's GAA) and returns one output, which we can interpret as the probability that a team will come back if a goalie is pulled in a given situation.
data = pulled_goalies[pulled_goalies["goal_diff"] == 1]
X = data[["time_left", "losingGPG", "winningGAA"]]
y = data["didTheyComeBack"]
reg = MLPRegressor(hidden_layer_sizes=(10,10), activation = 'logistic', max_iter=1000)
reg.fit(X, y)
# sit = [[time left, GPG, GAA]]
sit1 = [[60, 3, 3]]
sit2 = [[150, 2, 4]]
sit3 = [[300, 2, 4]]
sit4 = [[1200, 4, 3]]
print(reg.predict(sit1))
print(reg.predict(sit2))
print(reg.predict(sit3))
print(reg.predict(sit4))
There is a very clear problem with this model, and that's that it's very close to being proportional to the the time left in the game. Situation 4 above is one in which the losing team is down by 1 and there are 20 minutes left in the game; it's logically impossible for the losing team to have such a high chance of winning if they pull the goalie. Problems show up in other kinds of models as well:
data = pulled_goalies[pulled_goalies["goal_diff"] == 1]
X = data[["time_left", "losingGPG", "winningGAA"]]
y = data["didTheyComeBack"]
reg = SVR()
reg.fit(X, y)
# sit = [[time left, GPG, GAA]]
sit1 = [[60, 3, 3]]
sit2 = [[150, 3, 3]]
sit3 = [[300, 3, 3]]
sit4 = [[1200, 3, 3]]
print(reg.predict(sit1))
print(reg.predict(sit2))
print(reg.predict(sit3))
print(reg.predict(sit4))
This time we trained a support vector machine, and the numbers are practically indistinguishable from each other. Despite all of the data wrangling we did to get to this point, there is very little insight that we can gather from these machine learning models: the situation is simply too complicated to be predicted with the data points that we chose. So if we can't build a machine learning model that tells us when the optimal time to pull a goalie, what insights can we gain from this? Let's look at goal frequency, and how often a goal is scored at even strength versus during a pulled-goalie scenareo.
In this step we will calculate the average goals-per-game in "regular" situations, then calculate the goals-against frequency and goals-for frequency with the goalie pulled.
regular_goals_number = pulled_goalies["losingScore"].sum() + pulled_goalies["winningScore"].sum()
regular_goals_frequency = regular_goals_number/(3600 * pulled_goalies.shape[0] - pulled_goalies["time_left"].sum())/2
print(regular_goals_frequency) # Expected goals scored per second
print(regular_goals_frequency*3600) # Expected goals scored per game
The above block calculates the goals-per-second and goals-per-game that one can expect in all situations where the goalie is not pulled. In any given second, a team has a .0698% chance of scoring. This translates to about 2.51 goals per 60 minutes, or one goal per 23.9 minutes.
goals_for_frequency = pulled_goalies["goalsForLosingTeam"].sum()/pulled_goalies["time_left"].sum()
goals_against_frequency = pulled_goalies["didWinningTeamScore"].sum()/pulled_goalies["time_left"].sum()
print(goals_for_frequency*3600)
print(goals_against_frequency*3600)
The above block calculates the goals-per-game that one can expect 1) if your goalie is pulled and 2) if your opponent has pulled the goalie. If your goalie is pulled, you can expect to score about 3.49 goals per 60 minutes, or one goal per 17.2 minutes. If your opponent's goalie is pulled, you can expect to score about 7.07 goals per 60 minutes, or one goal per 8.49 minutes.
This alligns with what we already suspect about pulling the goalie: it is a negative EV decision, but while your opponent is much more likely to score a goal, your likelihood of scoring also increases.
In this next part, we will simulate 1000 games for every 10-second time interval, and we will determine the percent of times that the losing team comes back from their defecit after pulling their goalie during this time interval.
df = pd.DataFrame([], columns = ["defecit","time_left","even_strength_comeback_percentage","pull_goalie_comeback_percentage"])
escp = []
defecits = []
timeLeft = []
for defecit in [1]: # was [1,2], this took too long to run
for time in range(0,1200,10):
es_comebacks = 0
for i in range(1000): # simulate 1000 games with these conditions
goalsTotal = 0
for sec in range(time):
if random.random() < regular_goals_frequency:
goalsTotal += 1
elif random.random() < regular_goals_frequency:
goalsTotal -= 1
if goalsTotal >= defecit:
es_comebacks += 1
break
escp.append(es_comebacks/1000)
defecits.append(defecit)
timeLeft.append(time)
df["defecit"] = defecits
df["time_left"] = timeLeft
df["even_strength_comeback_percentage"] = escp
pgcp = []
for index, row in df.iterrows():
pg_comebacks = 0
for i in range(1000): # simulate 1000 games
goalsTotal = 0
for sec in range(row["time_left"]):
if random.random() < goals_for_frequency:
goalsTotal += 1
elif random.random() < goals_against_frequency:
goalsTotal -= 1
if goalsTotal >= row["defecit"]:
pg_comebacks += 1
break
pgcp.append(pg_comebacks/1000)
df["pull_goalie_comeback_percentage"] = pgcp
df["diff"] = df["pull_goalie_comeback_percentage"] - df["even_strength_comeback_percentage"]
with pd.option_context('display.max_rows', None, 'display.max_columns', None):
display(df.loc[df["defecit"] == 1])
Shown above is the table of all simulation results for games with a one-goal defecit: it contains the game's time, the percent of the time that the team came back from their defecit without ever pulling the goalie, the percent of the time that the team came back after pulling the goalie at time = "time_left", and the difference between these two percentages. Shown below is a scatter plot of the difference between the percentages versus time left in the game when the goalie was pulled.
temp = df.loc[df["defecit"] == 1]
ax = plt.scatter(temp["time_left"], temp["diff"])
plt.ylabel("Comeback Percent Difference")
plt.xlabel("Time Left In Game (Seconds)")
plt.show()
Coaches should want to give their team the greatest opportunity to succeed, so they should decide to pull the goalie when pulling the goalie gives them the best advantage versus leaving the goalie in. According to this graph, with a one-goal defecit this happens around the 300-second mark: that is the approximate maximum of this graph.
Conclusion:
This is an imperfect analysis, but it indicates something that many hockey enthusiasts would not expect: if teams want to win games, they should pull their goalies sooner. On average, successful teams pull their goalie with around 2 minutes and 30 seconds left in the game, which is a full 30 seconds sooner than unsuccessful teams. When examining times that teams score and get scored on with the goalie pulled, we also found that teams pulling their goalie with up to 5 minutes remaining could be the best possible option. Further analysis needs to be done, but this shows one thing: if hockey teams want to win games, they should pull their goalies earlier.