monott

### NHL PLAYOFF SIM

library(dplyr)
library(PlayerRatings)
library(magrittr)
library(ggplot2)
library(curleylab) #only using this for ggplot theme - not necessary to install if don't want. just remove 'curleytheme' in plot

#copy and paste initial ELO ratings from my nhlelo shiny app
#   https://jalapic.shinyapps.io/nhlelo/#

#End of Regular Season Ratings
ratings<-
  read.table(text="
             1 = NYR - 2310
             2 = STL - 2271
             3 = ANA - 2268
             4 = MIN - 2261
             5 = TBL - 2261
             6 = OTT - 2249
             7 = CBJ - 2247
             8 = WAS - 2246
             9 = VAN - 2246
             10 = MON - 2245
             11 = WIN - 2236
             12 = DAL - 2233
             13 = CHI - 2233
             14 = CGY - 2230
             15 = COL - 2219
             16 = NAS - 2210
             17 = NYI - 2207
             18 = LAK - 2203
             19 = DET - 2199
             20 = BOS - 2198
             21 = FLA - 2189
             22 = SJS - 2181
             23 = PIT - 2175
             24 = PHI - 2137
             25 = NJD - 2129
             26 = CAR - 2125
             27 = EDM - 2092
             28 = TOR - 2088
             29 = BUF - 2057
             30 = ARZ - 2056
             ", header=F, stringsAsFactors=F)

ratings %<>% select(V3,V5)
colnames(ratings)<-c("team", "initial")

# set up bracket
bracket <- data.frame(
  team=c("MON", "OTT", "TBL", "DET", "NYR", "PIT", "WAS", "NYI", "STL","MIN", "NAS", "CHI", "ANA", "WIN", "VAN", "CGY"),
  r1 = rep(LETTERS[1:8], each=2),
  r2 = rep(LETTERS[1:4], each=4),
  r3 = rep(LETTERS[1:2], each=8),
  r4 = rep(LETTERS[1:1], each=16),
  stringsAsFactors=F
)

bracket <- merge(bracket, ratings, by="team") %>% arrange(r1,r2,r3,r4)

## get 1st round games
r1sim <- bracket %$% split(., r1)

# existing results
teams <- r1sim[[1]][,1]
game1win <- "MON"
game2win <- "MON"
game1lose <- setdiff(teams,game1win)
game2lose <- setdiff(teams,game2win)

elodf.init <- data.frame(event=1:2, winner=c(game1win, game2win),loser=c(game1lose, game2lose),score=1, stringsAsFactors=F)



# initialize dataframe for Montreal-Ottawa
df0 <- r1sim[[1]] %>% select(Player=team, Rating=initial) %>% mutate(Games=0, Win=0, Draw=0, Loss=0, Lag=0)
df1 <- elo(elodf.init, kfac=15, status=df0)

# run sim  once to check it works

j<-1
tmp<-NULL
tmp[[1]]<-df1
  
while(TRUE) {
  
  tmp[[j+1]] <- getelo(tmp[[j]]$ratings)
  
  if( sum(tail(tmp,1)[[1]]$ratings[,4]==4)==1 )
    
    break
  
  j = j + 1
  
}



result <- tmp[[j+1]]
result





##### Run 5,000 times....

#make function
getresult<-function(df1) {
j<-1
tmp<-NULL
tmp[[1]]<-df1

while(TRUE) {
  
  tmp[[j+1]] <- getelo(tmp[[j]]$ratings)
  
  if( sum(tail(tmp,1)[[1]]$ratings[,4]==4)==1 )
    
    break
  
  j = j + 1
  
}



result <- tmp[[j+1]]$ratings
winrow<-which.max(result[,4])
loserow<-which.min(result[,4])
ddf <- data.frame(winner=result[winrow,1], loser=result[loserow,1], games=result[winrow,3]) 
return(ddf)
}

getresult(df1)

nperms <- 5000
monott <- NULL
for(i in 1:nperms) { monott[[i]] <- getresult(df1) }

monott.df <- do.call("rbind", monott)


#plotting

monott.df %>% 
  mutate(caption = paste(winner, "in", games)) %>%
  group_by(caption) %>%
  summarize(probability = n()/nperms) %>%
  ungroup() %$%
  ggplot(., aes(caption, probability)) + 
  xlab("") + ggtitle("Montreal - Ottawa Series Win Probabilites") +
    geom_bar(position="dodge", stat="identity") +
    curleytheme()


monott.df %>% 
  mutate(caption = paste(winner, "in", games)) %>%
  group_by(caption) %>%
  summarize(probability = n()/nperms)