“It’s Difficult to Make Predictions, Especially About the Future”

We’ve made our predictions about the past (estimating the relative strengths of teams based on past results), now we need to predict the future. I think it also nicely captures that even our predictions about the past are noisy- we can not ever truly know the exact strengths of football teams; the job of analytics is to estimate these are accurately as possible. But any noise in those past predictions will be carried forward and amplified when predicting the future.

Onward to the code, first as always, loading libraries and setting a seed for reproducibility:

library(tidyverse)
library(ggrepel)

set.seed(3459)

We’re then going to load all the stuff we prepped and predicted in the last post. Remember the α parameter below refers to a teams attacking strength (the relative number of goals they are expected to score), and the β parameter refers to the attacking strength (the inverse of the relative number of goals they are expected to concede). Finally, γ refers to the extra advantage of playing at home.

(all these files are on the github repo for this site)

fixtures <- readRDS("../../static/files/dc_2/dc_fixtures.rds")
results <- readRDS("../../static/files/dc_2/dc_results.rds")

model <- readRDS("../../static/files/dc_2/dc_model.rds")

model
## $alpha
##          Arsenal Blackburn_Rovers    Coventry_City   Dover_Athletic 
##        1.1106558        0.6370160        0.3023048       -0.2875353 
##     Enfield_Town    Frimley_Green 
##       -0.3767038       -1.3857376 
## 
## $beta
##          Arsenal Blackburn_Rovers    Coventry_City   Dover_Athletic 
##        0.6457175        0.4289270        0.3647815       -0.1362931 
##     Enfield_Town    Frimley_Green 
##       -0.3852812       -0.9178517 
## 
## $gamma
##    gamma 
## 0.189462

We’ll define a quick function to do our prediction. For a quick explanation of exactly why it’s coded as presented, see the previous post, under the title ‘Tinkering’.

For a given string of a home team and an away team, the function finds the relevant parameters from a third argument (param_list) and calculates the expected goal for each team.

predict_results <- function(home, away, param_list) {
  e_goals_home <- exp(param_list$alpha[home] - param_list$beta[away] + param_list$gamma)
  e_goals_away <- exp(param_list$alpha[away] - param_list$beta[home])
  
  df <- data.frame(home = home, away = away,
                   e_hgoal = as.numeric(e_goals_home), 
                   e_agoal = as.numeric(e_goals_away))
  
  return(df)
}

If we run this for two example teams for example:

#two example teams
home <- "Blackburn_Rovers"
away <- "Arsenal"

prediction <- predict_results(home, away, model) 
prediction
##               home    away  e_hgoal  e_agoal
## 1 Blackburn_Rovers Arsenal 1.198128 1.977293

We can see that it gives Arsenal (the away team) a slightly more optimistic chance than Blackburn. The expected goals for each team of course can be rewritten as the mean, and in our Poisson model refers to λ (lambda)- the mean times an event (goal) happens per a time interval (match). We also set a maximum number of possible goals (7 in this case*) to bound the area under the distribution so we aren’t sampling forever.

*sharp readers might notice that this is actually lower than the lambda for our more extreme cases (e.g. Arsenal at home to Frimley Green), but for realistic matches (even between wildly different professional sides) this is a fair enough assumption.

We then use dpois() to calculate the probability of this Poisson function returning a value (0:7 goals) given it’s lambda value. So if we run this over the prediction we made for Blackburn Rovers vs. Arsenal we get:

#set a limit of where we'll calculate across
max_goals <- 7

#calculate the probability of scoring x goals for either team
blackburn_goal_probs <- lapply(0:max_goals, dpois, lambda = prediction$e_hgoal)
arsenal_goal_probs <- lapply(0:max_goals, dpois, lambda = prediction$e_agoal)

#bind together in a df
df <- data.frame(goals = rep(0:max_goals, 2),
                 team = rep(c(home, away), each = max_goals+1),
                 p = c(unlist(blackburn_goal_probs), unlist(arsenal_goal_probs)))

#plot the p of scoring x goals for either team
p1 <- ggplot(df, aes(x = goals, y = p, fill = team)) +
  geom_density(stat = "identity", alpha = 0.5) +
  scale_fill_manual(values = c("red", "blue")) +
  labs(title = "Predicted goals for Blackburn Rovers and Arsenal",
       y = "probability") +
  theme_minimal()

p1

Because of how maths works, these curves are the same result we would get if we ran rpois() (sampling from the Poisson function) lots of times. We’ll do that quickly because it sets the stage nicely for what will come later.

#sample from the function lots of times for each team
n <- 100000
blackburn_goals_samples <- rpois(n, lambda = prediction$e_hgoal)
arsenal_goals_samples <- rpois(n, lambda = prediction$e_agoal)

df <- data.frame(team = rep(c(home, away), each = n),
                 sampled_goals = c(blackburn_goals_samples, arsenal_goals_samples))

#look its the same plot!
p2 <- ggplot(df, aes(x = sampled_goals, fill = team)) +
  geom_bar(stat = "count", position = "dodge", colour = "black", alpha = 0.5) +
  geom_line(aes(colour = team), stat = "count", size = 3) +
  scale_fill_manual(values = c("red", "blue"), guide = FALSE) +
  scale_colour_manual(values = c("red", "blue"), guide = FALSE) +
  labs(title = "Predicted goals for Blackburn Rovers and Arsenal",
       y = "probability",
       x = "sampled goals") +
  theme_minimal() +
  theme(axis.text.y = element_blank())

p2

Ok great!, in terms of predicting the result, the rightwards shift of the red (Arsenal) curve here is the difference in the teams ability to generate a positive goal differential- it makes it more likely that if we sample event, Arsenal will have scored more goals than Blackburn Rovers at the end of the match. Of course, it’s also obvious that while Arsenal’s curve is right shifted, the bars for Arsenal scoring 0 goals and Blackburn scoring 6 are still sizable enough that it isn’t outside the realm of possibility.

This is a nice way of presenting the chance of each team scoring n goals, but doesn’t really help us in predicting the result of a match given that this relies on the interaction of both these distributions (we need to know how many goals BOTH Arsenal AND Blackburn will score).

To calculate this, we can do an outer product of the probabilities for both teams scoring n goals. We can then plot the probability of each scoreline as a tile plot:

#calculate matrix of possible results and probabilities of those
matrix <- outer(unlist(arsenal_goal_probs), unlist(blackburn_goal_probs)) %>%
  as.data.frame() %>%
  gather() %>%
  #add in scorelines
  mutate(hgoals = rep(0:max_goals, max_goals+1),
         agoals = rep(0:max_goals, each = max_goals+1))

#make the tile plot
p3 <- ggplot(matrix, aes(x = hgoals, y = agoals, fill = value)) +
  geom_tile() +
  geom_text(aes(label = paste(hgoals, agoals, sep = "-"))) +
  scale_fill_gradient2(low = "white", high = "red", guide = FALSE) +
  theme_minimal()

p3

Where we can see that the most common scorelines are low scoring (football is a low scoring game), and slightly biased towards away goals (i.e. Arsenal are more likely to win than lose). The darkest (most likely) tiles being 1-1 or a 2-1 Arsenal win seem very plausible given our calculated λs earlier.

We can then do this for every fixture and build a large graph of the expected results for each using a simple map2_ apply. Because of the huge plot, I’ve restricted it here to a 3x3 matrix of the results for Arsenal, Coventry City, and Enfield Town, but if you click you should be linked to the full image.

#want to predict over the whole fixture space
all_fixtures <- bind_rows(fixtures, results) %>%
  filter(!duplicated(paste(home, away), fromLast = TRUE))

#get the lambda for each team per game
predictions <- map2_df(all_fixtures$home, all_fixtures$away, 
                       predict_results,
                       model)

#calc out probabilities and bind up
all_predictions <- map2_df(
  predictions$e_hgoal, predictions$e_agoal, 
  function(lambda_home, lambda_away, max_goals) {
    hgoal_prob <- dpois(0:max_goals, lambda_home) %>% `names<-`(0:max_goals)
    agoal_prob <- dpois(0:max_goals, lambda_away) %>% `names<-`(0:max_goals)
    outer(hgoal_prob, agoal_prob) %>%
      as.data.frame() %>% 
      gather() %>% 
      rownames_to_column("row") %>%
      mutate(hgoal = as.numeric(row) %% (max_goals+1)-1) %>% 
      mutate(hgoal = case_when(hgoal < 0 ~ max_goals, TRUE ~ hgoal),
             agoal = as.numeric(key)) %>%
      select(sample_hgoal = hgoal, sample_agoal = agoal, prob = value)
}, max_goals) %>%
  cbind(all_fixtures[rep(seq_len(nrow(all_fixtures)), each=(max_goals+1)^2),], .) %>%
  group_by(home, away) %>%
  mutate(prob = prob / sum(prob)) %>%
  ungroup()

#plot again
p3 <- all_predictions %>%
  #filter only a few out to scale plot 
  filter(home %in% c("Arsenal", "Coventry_City", "Enfield_Town"),
         away %in% c("Arsenal", "Coventry_City", "Enfield_Town")) %>%
  ggplot(aes(x = sample_hgoal, y = sample_agoal, fill = prob)) +
  geom_tile() +
  geom_point(aes(x = hgoal, y = agoal), 
             colour = "blue", size = 5, alpha = 0.5 / max_goals^2) +
  geom_text(aes(label = paste(sample_hgoal, sample_agoal, sep = "-")), size = 2.3) +
  scale_fill_gradient2(low = "white", high = "red", guide = FALSE) +
  labs(
    title = "predictions for final score across all fixtures",
    y = "away goals",
    x = "home goals") +
  theme_minimal() +
  facet_grid(away~home, scales = "free")

p3

For the whole matrix, click here

So what?

These graphs are nice, but whats important is what they show: we have a way to quantify how likely any result is in a match between two given teams. Why is this useful

• Firstly, we can use the output of this to build betting models. Given the odds on final scores for any match, we can hedge effectively by betting on (e.g.) the five overwhelmingly most likely results.
• Secondly, we can simulate leagues. This is perhaps especially of interest given the context of writing this post. I’m going to focus on this application because I don’t bet on football, and also because it’s hard to get a nice database of odds at the moment given the aforementioned situation.

The Verve - Monte Carlo

We can do this using a technique called Monte Carlo simulation. There are lots of good explanation of the technique on the internet, but it basically boils down to this:

nested_probabilities <- all_predictions %>%
  filter(is.na(hgoal)) %>%
  select(-hgoal, -agoal) %>%
  nest(probabilities = c(sample_hgoal, sample_agoal, prob))

nested_probabilities$probabilities[[1]] %>%
  rename("Coventry City" = sample_hgoal, "Arsenal" = sample_agoal) %>%
  arrange(-prob) %>%
  #show first 15 rows
  .[1:15,]
## # A tibble: 15 x 3
##    `Coventry City` Arsenal   prob
##              <dbl>   <dbl>  <dbl>
##  1               0       2 0.115 
##  2               0       1 0.109 
##  3               1       2 0.0983
##  4               1       1 0.0933
##  5               0       3 0.0806
##  6               1       3 0.0691
##  7               0       0 0.0516
##  8               1       0 0.0442
##  9               0       4 0.0425
## 10               2       2 0.0422
## 11               2       1 0.0400
## 12               1       4 0.0364
## 13               2       3 0.0296
## 14               2       0 0.0190
## 15               0       5 0.0179

sum(nested_probabilities$probabilities[[1]]$prob)
## [1] 1

nested_probabilities$probabilities[[1]] %>%
  rename("Coventry_City" = sample_hgoal, "Arsenal" = sample_agoal) %>%
  sample_n(1, weight = prob)
## # A tibble: 1 x 3
##   Coventry_City Arsenal   prob
##           <dbl>   <dbl>  <dbl>
## 1             1       3 0.0691

nested_probabilities %>%
    mutate(sampled_result = map(probabilities, sample_n, 1, weight = prob)) %>%
    select(-probabilities) %>%
    unnest(cols = c(sampled_result))
## # A tibble: 6 x 6
##   home             away             gameweek sample_hgoal sample_agoal   prob
##   <chr>            <chr>               <int>        <dbl>        <dbl>  <dbl>
## 1 Coventry_City    Arsenal                 9            0            5 0.0179
## 2 Blackburn_Rovers Dover_Athletic          9            1            1 0.0575
## 3 Frimley_Green    Enfield_Town            9            0            4 0.0418
## 4 Arsenal          Blackburn_Rovers       10            2            1 0.0966
## 5 Coventry_City    Frimley_Green          10            3            0 0.170 
## 6 Dover_Athletic   Enfield_Town           10            2            1 0.0839

rerun(10, nested_probabilities %>%
    filter(home == "Coventry_City" & away == "Arsenal") %>%
    mutate(sampled_result = map(probabilities, sample_n, 1, weight = prob)) %>%
    select(-probabilities) %>%
    unnest(cols = c(sampled_result))
) %>%
  bind_rows() %>%
  arrange(-prob)
## # A tibble: 10 x 6
##    home          away    gameweek sample_hgoal sample_agoal   prob
##    <chr>         <chr>      <int>        <dbl>        <dbl>  <dbl>
##  1 Coventry_City Arsenal        9            0            2 0.115 
##  2 Coventry_City Arsenal        9            1            2 0.0983
##  3 Coventry_City Arsenal        9            1            1 0.0933
##  4 Coventry_City Arsenal        9            0            3 0.0806
##  5 Coventry_City Arsenal        9            1            3 0.0691
##  6 Coventry_City Arsenal        9            1            0 0.0442
##  7 Coventry_City Arsenal        9            0            4 0.0425
##  8 Coventry_City Arsenal        9            0            4 0.0425
##  9 Coventry_City Arsenal        9            0            4 0.0425
## 10 Coventry_City Arsenal        9            1            5 0.0154

n <- 100
fixture_sims <- rerun(n, nested_probabilities %>%
    mutate(sampled_result = map(probabilities, sample_n, 1, weight = prob)) %>%
    select(-probabilities) %>%
    unnest(cols = c(sampled_result)) %>%
  select(-gameweek, -prob) %>%
  pivot_longer(c(home, away), names_to = "location", values_to = "team") %>%
  mutate(points = case_when(
    location == "home" & sample_hgoal > sample_agoal ~ 3,
    location == "away" & sample_agoal > sample_hgoal ~ 3,
    sample_hgoal == sample_agoal ~ 1,
    TRUE ~ 0
  )) %>%
  mutate(gd = case_when(
    location == "home" ~ sample_hgoal - sample_agoal,
    location == "away" ~ sample_agoal - sample_hgoal
  )))

fixture_sims[1]
## [[1]]
## # A tibble: 12 x 6
##    sample_hgoal sample_agoal location team             points    gd
##           <dbl>        <dbl> <chr>    <chr>             <dbl> <dbl>
##  1            0            0 home     Coventry_City         1     0
##  2            0            0 away     Arsenal               1     0
##  3            4            0 home     Blackburn_Rovers      3     4
##  4            4            0 away     Dover_Athletic        0    -4
##  5            0            0 home     Frimley_Green         1     0
##  6            0            0 away     Enfield_Town          1     0
##  7            3            0 home     Arsenal               3     3
##  8            3            0 away     Blackburn_Rovers      0    -3
##  9            6            1 home     Coventry_City         3     5
## 10            6            1 away     Frimley_Green         0    -5
## 11            1            1 home     Dover_Athletic        1     0
## 12            1            1 away     Enfield_Town          1     0

fixture_sims %>%
  bind_rows() %>%
  group_by(team) %>%
  summarise(av_points = sum(points)/n,
            av_gd = sum(gd) / n)
## # A tibble: 6 x 3
##   team             av_points av_gd
##   <chr>                <dbl> <dbl>
## 1 Arsenal               4.19  2.44
## 2 Blackburn_Rovers      3.16  0.7 
## 3 Coventry_City         3.61  2.42
## 4 Dover_Athletic        2.26 -1.26
## 5 Enfield_Town          2.95  0.23
## 6 Frimley_Green         0.6  -4.53

table <- results %>%
  pivot_longer(c(home, away), names_to = "location", values_to = "team") %>%
  mutate(points = case_when(
    location == "home" & hgoal > agoal ~ 3,
    location == "away" & agoal > hgoal ~ 3,
    hgoal == agoal ~ 1,
    TRUE ~ 0
  )) %>%
  mutate(gd = case_when(
    location == "home" ~ hgoal - agoal,
    location == "away" ~ agoal - hgoal
  )) %>%
  group_by(team) %>%
  summarise(points = sum(points),
            gd = sum(gd))

predicted_finishes <- map_df(fixture_sims, function(simulated_fixtures, table) {
  simulated_fixtures %>%
    select(team, points, gd) %>%
    bind_rows(., table) %>%
    group_by(team) %>%
    summarise(points = sum(points),
              gd = sum(gd)) %>%
    arrange(-points, -gd) %>%
    mutate(predicted_finish = 1:n())
}, table) %>%
  group_by(team, predicted_finish) %>%
  summarise(perc = n() / n)

predicted_finishes
## # A tibble: 10 x 3
## # Groups:   team [6]
##    team             predicted_finish  perc
##    <chr>                       <int> <dbl>
##  1 Arsenal                         1  0.82
##  2 Arsenal                         2  0.18
##  3 Blackburn_Rovers                1  0.18
##  4 Blackburn_Rovers                2  0.82
##  5 Coventry_City                   3  0.97
##  6 Coventry_City                   4  0.03
##  7 Dover_Athletic                  3  0.03
##  8 Dover_Athletic                  4  0.97
##  9 Enfield_Town                    5  1   
## 10 Frimley_Green                   6  1

p4 <- ggplot(predicted_finishes, aes(x = predicted_finish, y = perc, fill = team)) +
  geom_bar(stat = "identity", colour = "black") +
  scale_fill_manual(values = c("red", "blue", "skyblue", "white", "dodgerblue4", "blue")) +
  labs(
    title = "Predicted finish position of teams",
    subtitle = "with two gameweeks left to play",
    y = "fraction of finishes",
    x = "final position"
  ) +
  theme_minimal() +
  facet_wrap(~team)

p4

The Real Thing

We’re now at the stage where we can start to look at real data. One of the motivating forces which drew me back to this putative blog series was the current football situation- with season ending with games left to play.

We can use the knowledge we’ve built up over these last posts to see what we expect to happen in these unplayed games, if they cannot be completed.

To make code more concise, I’ve used Ben Torvaney’s code in his regista package (he also has some nice usage blogs similar to this post at his blog). The underlying maths is exactly the same as in my previous post though with a few different design choices. If we run the simulations using the code from the previous post we should get exactly the same answer.

The code following is also extremely similar to the final chunks of one of my previous posts in analysing the current Liverpool team’s achievements.

library(rvest)
library(regista)

First we need to download the data on the current English Premier League season. Once we have this we can split it into played matches (where we 100% know the result) and unplayed matches which we need to predict the result of. For the basis of the team strength estimates I’ve used the xg created and allowed per game, as I believe these give a better estimate of team strength (indeed Ben Torvaney has a nice post on using even the shot-by-shot xg to produce Dixon-Coles models).

#download the match data from 2019/2020
fixtures_2020 <- "https://fbref.com/en/comps/9/schedule/Premier-League-Fixtures" %>%
  read_html() %>%
  html_nodes("#sched_ks_3232_1") %>%
  html_table() %>%
  as.data.frame() %>%
  separate(Score, into = c("hgoal", "agoal"), sep = "–") %>%
  #only care about goals and expected goals
  select(home = Home, away = Away, home_xg = xG, away_xg = xG.1, hgoal, agoal) %>%
  filter(home != "") %>%
  mutate(home = factor(home), away = factor(away)) %>%
  #round expected goals to nearest integer
  mutate_at(c("home_xg", "away_xg", "hgoal", "agoal"), .funs = funs(round(as.numeric(.))))

#matches with a known result
#used for modelling
played_matches <- fixtures_2020 %>%
  filter(!is.na(home_xg))

#matches with an unknown result
#used for simulation
unplayed_matches <- fixtures_2020 %>%
  filter(is.na(home_xg)) %>%
  select_if(negate(is.numeric))

#fit the dixon coles model
#use xg per game, not 'actual' goals
fit_2020 <- dixoncoles(home_xg, away_xg, home, away, data = played_matches)

To get a look at what the team parameters in a real-life league look like we can extract them from the model and plot them:

#extract Dixon-Coles team strenth parameters
pars_2020 <- fit_2020$par %>%
  .[grepl("def_|off_", names(.))] %>%
  matrix(., ncol = 2) %>%
  as.data.frame() %>%
  rename(attack = V1, defence = V2)
pars_2020$team <- unique(gsub("def_*|off_*", "", names(fit_2020$par)))[1:20]

#plot as before
p5 <- pars_2020 %>%
  mutate(defence = 1 - defence) %>%
  ggplot(aes(x = attack, y = defence, colour = attack + defence, label = team)) +
  geom_point(size = 3, alpha = 0.7) +
  geom_text_repel() +
  scale_colour_continuous(guide = FALSE) +
  labs(title = "Dixon-Coles parameters for the 2019/2020 EPL",
       x = "attacking strength",
       y = "defensive strength") +
  theme_minimal()

p5

It might surprise some that Manchester City are predicted to be better than Liverpool by this model, but it shouldn’t given the underlying numbers for both teams. Liverpool have run very hot and Manchester City have run very cold this season.

Finally, we can then calculate the current Premier League table, and simulate remaining games to predict where teams will finish the season if the remainder of games were to be played. I’ve chosen 1000 sims just for sake of processing time, but you can scale up and down as desired.

#calculate the current EPL table
current_epl_table <- played_matches %>%
  select(home, away, hgoal, agoal) %>%
  pivot_longer(c(home, away), names_to = "location", values_to = "team") %>%
  mutate(points = case_when(
    location == "home" & hgoal > agoal ~ 3,
    location == "away" & agoal > hgoal ~ 3,
    hgoal == agoal ~ 1,
    TRUE ~ 0
  )) %>%
  mutate(gd = case_when(
    location == "home" ~ hgoal - agoal,
    location == "away" ~ agoal - hgoal
  )) %>%
  group_by(team) %>%
  summarise(points = sum(points),
            gd = sum(gd))

#the number of sims to run
n <- 10000

#simulate remaining matches
fixture_sims_2020 <- rerun(
  n,
  augment.dixoncoles(fit_2020, unplayed_matches, type.predict = "scorelines") %>%                     
    mutate(sampled_result = map(.scorelines, sample_n, 1, weight = prob)) %>%
    select(-.scorelines) %>%
    unnest(cols = c(sampled_result)) %>%
  pivot_longer(c(home, away), names_to = "location", values_to = "team") %>%
  mutate(points = case_when(
    location == "home" & hgoal > agoal ~ 3,
    location == "away" & agoal > hgoal ~ 3,
    hgoal == agoal ~ 1,
    TRUE ~ 0
  )) %>%
  mutate(gd = case_when(
    location == "home" ~ hgoal - agoal,
    location == "away" ~ agoal - hgoal
  )) %>%
    select(team, points, gd))

#calculate final EPL tables
predicted_finishes_2020 <- map_df(fixture_sims_2020, function(sim_fixtures, table) {
  sim_fixtures %>%
    select(team, points, gd) %>%
    bind_rows(., table) %>%
    group_by(team) %>%
    summarise(points = sum(points),
              gd = sum(gd)) %>%
    arrange(-points, -gd) %>%
    mutate(predicted_finish = 1:n())
}, current_epl_table) %>%
  group_by(team, predicted_finish) %>%
  summarise(perc = n() / n) %>%
  group_by(team) %>%
  mutate(mean_finish = mean(predicted_finish)) %>%
  arrange(mean_finish) %>%
  ungroup() %>%
  mutate(team = factor(team, levels = unique(team)))

If we then plot these predicted finishes (ordered by the chance of their highest finish position), we can get an idea of where we expect teams to end the season:

#list of team colours
team_cols <- c("red", "skyblue", "darkblue", "darkblue", "darkred",
               "orange", "red", "white", "red", "blue", "maroon", 
               "blue", "white", "red", "dodgerblue", "yellow", 
               "maroon", "red", "maroon", "yellow")

#plot the finishing position by chance based on these simualtions
p6 <- ggplot(predicted_finishes_2020, 
             aes(x = predicted_finish, y = perc, fill = team)) +
  geom_bar(stat = "identity", colour = "black") +
  scale_fill_manual(values = team_cols, guide = FALSE) +
  labs(
    title = "Predicted finish position of teams",
    subtitle = "for incomplete 2019/2020 EPL season",
    y = "fraction of finishes",
    x = "final position"
  ) +
  theme_minimal() +
  facet_wrap(~team)

p6

So great news for Liverpool fans who the model believes have a 100% chance of finishing in first place. Leicester also might be happy with a nailed on 3rd place, with Chelsea or Manchester United probably rounding out the top four, and Wolves joining the loser of the two in the Europa League.

#get the predictions for the 2019/2020 champion
winner <- predicted_finishes_2020 %>%
  filter(predicted_finish < 2)%>%
  mutate(prediction = "Champion chance")

winner
## # A tibble: 2 x 5
##   team            predicted_finish   perc mean_finish prediction     
##   <fct>                      <int>  <dbl>       <dbl> <chr>          
## 1 Liverpool                      1 1.00           1.5 Champion chance
## 2 Manchester City                1 0.0001         2.5 Champion chance
#get prediction for those who qualify for champions league
#and for europa league
champs_league <- predicted_finishes_2020 %>%
  filter(predicted_finish < 5) %>%
  group_by(team) %>%
  summarise(perc = sum(perc)) %>%
  arrange(-perc) %>%
  mutate(prediction = "Champions League chance")

champs_league
## # A tibble: 10 x 3
##    team              perc prediction             
##    <fct>            <dbl> <chr>                  
##  1 Liverpool       1      Champions League chance
##  2 Manchester City 1      Champions League chance
##  3 Leicester City  0.933  Champions League chance
##  4 Chelsea         0.479  Champions League chance
##  5 Manchester Utd  0.46   Champions League chance
##  6 Wolves          0.106  Champions League chance
##  7 Sheffield Utd   0.0155 Champions League chance
##  8 Tottenham       0.004  Champions League chance
##  9 Arsenal         0.0018 Champions League chance
## 10 Everton         0.0005 Champions League chance
europa_league  <- predicted_finishes_2020 %>%
  filter(predicted_finish < 7) %>%
  group_by(team) %>%
  summarise(perc = sum(perc)) %>%
  arrange(-perc) %>%
  mutate(prediction = "(at least) Europa League chance")

europa_league
## # A tibble: 13 x 3
##    team              perc prediction                     
##    <fct>            <dbl> <chr>                          
##  1 Liverpool       1      (at least) Europa League chance
##  2 Manchester City 1      (at least) Europa League chance
##  3 Leicester City  0.999  (at least) Europa League chance
##  4 Manchester Utd  0.954  (at least) Europa League chance
##  5 Chelsea         0.954  (at least) Europa League chance
##  6 Wolves          0.729  (at least) Europa League chance
##  7 Sheffield Utd   0.196  (at least) Europa League chance
##  8 Tottenham       0.096  (at least) Europa League chance
##  9 Arsenal         0.0479 (at least) Europa League chance
## 10 Everton         0.0139 (at least) Europa League chance
## 11 Burnley         0.0089 (at least) Europa League chance
## 12 Crystal Palace  0.0009 (at least) Europa League chance
## 13 Southampton     0.0008 (at least) Europa League chance

(obviously this model does not account for any ramifications of Manchester City’s European ban)

At the foot of the table, the model is fairly bullish on Norwich being relegated, with Aston Villa probably joining them, and then probably West Ham rounding out the relegation spots.

#get predictions for those who would be relegated
relegated  <- predicted_finishes_2020 %>%
  filter(predicted_finish > 17) %>%
  group_by(team) %>%
  summarise(perc = sum(perc)) %>%
  arrange(-perc) %>%
  mutate(prediction = "Relegation chance")

relegated
## # A tibble: 8 x 3
##   team             perc prediction       
##   <fct>           <dbl> <chr>            
## 1 Norwich City  0.934   Relegation chance
## 2 Aston Villa   0.700   Relegation chance
## 3 Bournemouth   0.507   Relegation chance
## 4 West Ham      0.402   Relegation chance
## 5 Watford       0.270   Relegation chance
## 6 Brighton      0.171   Relegation chance
## 7 Newcastle Utd 0.0126  Relegation chance
## 8 Southampton   0.00270 Relegation chance

Final Remarks

I want to make it clear at the end of this post that this probably isn’t the most sophisticated model for predicting football matches (more to come in a part 3, maybe this time within less than a year), but does a pretty good job!

In any case though, I don’t think that running these sims is a good way to end the season- in truth, there’s probably no good way. This post is more about how to use this technique than whether to use it.

Best, stay safe as always!