# Written June 22nd, 2006 # Author: Markus Triska triska@metalevel.at # Public domain code. # Adjust the matrix m of success probabilities to run your own # simulations. m <- matrix(rep(0.5, 32^2), ncol=32) m[1,] <- 0.9 m[,1] <- 0.1 ###################################################################### win <- function(pair) { npair <- as.numeric(unlist(strsplit(pair, " "))) p1 <- npair[1] p2 <- npair[2] prob1 <- m[p1,p2] w <- sample(c(p1,p2), 1, prob=c(prob1, 1-prob1)) return(w) } winners <- function(group, n) { o <- outer(group, group, "paste") pairs <- o[upper.tri(o)] wins <- unlist(lapply(pairs, win)) t <- table(wins) val <- (as.numeric(names(sort(t, decreasing=T)[1:n]))) # print(val) return(val) } worldcup <- function() { # stage 1 (group) g1 <- winners(1:4, 2) g2 <- winners(5:8, 2) g3 <- winners(9:12, 2) g4 <- winners(13:16, 2) g5 <- winners(17:20, 2) g6 <- winners(21:24, 2) g7 <- winners(25:28, 2) g8 <- winners(29:32, 2) # stage 2 (knockout) k1 <- winners(c(g1[1], g2[2]), 1) k2 <- winners(c(g1[2], g2[1]), 1) k3 <- winners(c(g3[1], g4[2]), 1) k4 <- winners(c(g3[2], g4[1]), 1) k5 <- winners(c(g5[1], g6[2]), 1) k6 <- winners(c(g5[2], g6[1]), 1) k7 <- winners(c(g7[1], g8[2]), 1) k8 <- winners(c(g7[2], g8[1]), 1) # quarter finals q1 <- winners(c(k1, k2), 1) q2 <- winners(c(k3, k4), 1) q3 <- winners(c(k5, k6), 1) q4 <- winners(c(k7, k8), 1) # semi finals s1 <- winners(c(q1, q2), 1) s2 <- winners(c(q3, q4), 1) # final f <- winners(c(s1, s2), 1) return(f) } simulations <- function(n) { s <- c() for (i in 1:n) { print(paste("simulating ", i, "/", n)) s <- c(worldcup(), s) } return(s) }