Haaland Goal Forecast Evalution

Last fall, I put together an algorithm to forecast how many goals to expect from Erling Haaland in the Premier League 2022-23 season. The median then was 41, with 90% of simulations coming between 32 and 51. Now that the season is done and we know his actual final total (36), we can go back and see how the model performed over time. The blue line is the median forecast, dark grey contains 90% of sims and light grey contains 95%. The dashed line is his final tally of 36.

library(hrbrthemes)
library(ggridges)
library(plotly)
library(tidyverse)
library(worldfootballR)

set.seed(82722)

haaland_summary_2022 <- fb_player_match_logs("https://fbref.com/en/players/1f44ac21/Erling-Haaland", season_end_year = 2022, stat_type = 'summary')
haaland_summary_2023 <- fb_player_match_logs("https://fbref.com/en/players/1f44ac21/Erling-Haaland", season_end_year = 2023, stat_type = 'summary')

haaland_game_logs <- haaland_summary_2022 %>%
  rbind(haaland_summary_2023) %>%
  as_tibble()

injury_days <- tibble(Date = as.Date('2022-10-29'),
                      Comp = 'Premier League',
                      Round = 'Matchweek 14',
                      Opponent = 'Leicester City',
                      Gls = NA,
                      xG_Expected = NA) %>%
  bind_rows(tibble(Date = as.Date('2023-04-01'),
                   Comp = 'Premier League',
                   Round = 'Matchweek 29',
                   Opponent = 'Liverpool',
                   Gls = NA,
                   xG_Expected = NA))

eh_2023_pl_goals <- haaland_summary_2023 %>%
  filter(Comp == 'Premier League') %>%
  summarise(goals = sum(Gls_Performance, na.rm = TRUE)) %>%
  as.numeric()

haaland_prepped <- haaland_game_logs %>%
  filter(Comp %in% c('Bundesliga', 'Premier League')) %>%
  rename(Gls = Gls_Performance) %>%
  select(Date, Comp, Round, Opponent, Gls, xG_Expected) %>%
  mutate(Date = as.Date(Date)) %>%
  rbind(injury_days) %>%
  arrange(Date) %>%
  mutate(rownum = row_number(),
         running_Gls_avg = cummean(ifelse(is.na(Gls), 0, Gls)))
  
haaland_gls_vctrs <- NULL

for(i in 1:max(haaland_prepped$rownum)){
  
  tmp_vector <- haaland_prepped %>%
    filter(rownum <= i & is.na(Gls) == FALSE) %>%
    select(Gls) %>%
    as.vector()
  
  tmp_tibble <- i %>%
    as_tibble %>%
    mutate(gls_vector = tmp_vector) %>%
    rename(rownum = value)
  
  haaland_gls_vctrs <- haaland_gls_vctrs %>%
    rbind(tmp_tibble)
}

haaland_prepped_final <- haaland_prepped %>%
  group_by(Comp, rownum, Gls) %>%
  nest() %>%
  inner_join(haaland_gls_vctrs) %>%
  filter(Comp == 'Premier League') %>%
  ungroup() %>%
  mutate(total_gls = cumsum(ifelse(is.na(Gls), 0, Gls)))

mc_sim <- function(goal_vector, n_games, goals_scored){
  
  sample_sim <- NULL
  
  for(i in 1:10000){
    sim <- sum(sample(goal_vector, size = n_games, replace = T)) + goals_scored
    sample_sim[i] <- sim  
  }
  
  sample_sim
  #quantile(sample_sim, c(.05, .5, .95))
}

poisson_sim <- function(goal_vector, n_games, goals_scored){
  p_sim <- NULL
  
  for(i in 1:10000){
    sim <- sum(rpois(n_games, mean(goal_vector))) + goals_scored
    p_sim[i] <- sim  
  }
  
  p_sim
}

haaland_simmed <- haaland_prepped_final %>%
  mutate(games_left = 38 - row_number()) %>%
  #rowwise() %>%
  mutate(simmed_data = pmap(list(gls_vector, games_left, total_gls), mc_sim)) %>%
  mutate(p_simmed_data = pmap(list(gls_vector, games_left, total_gls), poisson_sim)) %>%
  unnest(data) %>%
  rowwise() %>%
  mutate(median_sim = median(simmed_data),
         lower_025 = quantile(simmed_data, .025),
         lower_05 = quantile(simmed_data, .05),
         upper_95 = quantile(simmed_data, .95),
         upper_975 = quantile(simmed_data, .975)) %>%
  mutate(Matchweek = gsub('Matchweek ', '', Round)) %>%
  ungroup() %>%
  mutate(Matchweek = fct_reorder(Matchweek, -rownum))

haaland_sim_trend_chart <- haaland_simmed %>%
  ggplot(aes(x = as.factor(Date), 
             y = median_sim,
             text = paste('Date:', Date, 
                          "<br>Goals:", Gls,
                          "<br>Goals to Date:", total_gls,
                          "<br>Forecasted Total:", median_sim,
                          "<br><br>90% Range:", lower_05, "to", upper_95,
                          "<br>95% Range:", lower_025, "to", upper_975)
             ))+
  geom_ribbon(aes(ymin = lower_025, ymax = upper_975), 
              fill = 'grey80',
              alpha = 0.5,
              group = 1)+
  geom_ribbon(aes(ymin = lower_05, ymax = upper_95), 
              fill = 'grey70', 
              alpha = 0.5,
              group = 1)+
  geom_line(group = 1, color = '#6CABDD')+
  geom_hline(yintercept = eh_2023_pl_goals, linetype = 'dashed')+
  #geom_text(aes(as.Date('2022-08-14'),32, label = "Salah's Record", vjust = 6), size = 3)+
  ggtitle('Erling Haaland 2022-23 Projected Goal Total by Week')+
  theme_ipsum()+
  xlab('Date')+
  theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1),
        panel.grid.major.x = element_blank(), 
        panel.grid.minor.x = element_blank(),
        panel.grid.major.y = element_line(color = "grey", 
                                          linewidth = 0.5,
                                          linetype = 2))

ggplotly(haaland_sim_trend_chart, tooltip = 'text') %>% 
  layout(hoverlabel = list(bgcolor = "#6CABDD"))

The model did a pretty decent job at the start - week one it forecasted 37, so only one off. But several early hat tricks got it overly excited, to the point it forecasted 48 on Oct. 2 and our wider error bars did not include 36. If got more reasonable over time, which we’d expect as it learns more about Haaland’s goals per game distribution and we get more actual games than forecasted.

If I were to to this again, I’d take the time to pull the goal distributions for strikers at the same age on similar quality teams. This would add more 0s and 1s, which would temper expectations a bit.* Of course, I built this model for one player, so perhaps testing on the rest of the league would yield some more insights into what can be improved.

*Similar to how we’d adjust batting averages in baseball (see this fantastic walkthrough), though obviously the nuts and bolts would be different.

Using the current model, we can take a look at potential outcomes for Haaland’s next PL campaign:

pl23_gls_vector <- haaland_summary_2023 %>%
  select(Gls_Performance) %>%
  filter(!is.na(Gls_Performance)) %>%
  as.vector() 

p24_sim <- NULL

for(i in 1:10000){
  sim <- rpois(38, mean(pl23_gls_vector$Gls_Performance))
  p24_sim[i] <- sum(sim)
}

summary(p24_sim)

   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
  17.00   34.00   38.00   38.63   43.00   62.00 

p24_sim %>%
  as_tibble() %>%
  ggplot(aes(x = value))+
  geom_histogram(binwidth = 1)+
  theme_ipsum()

This seems reasonable - 39 isn’t a massive jump from this season, and presumably he should get better as he continues to improve and City figures out how to get the ball to him more often. The min and max values are pretty silly - 15 is far too low, and 62 would break the game - but saying “we think he’ll get between 34 and 43 goals next season” is pretty hard to argue with.