Probability of winning a game in tennis?
Here is my attempt, I would appreciate feedback:
Let $p$ be the probability that player A wins a single point.
(1) Player A can win after 40-0, with probability $p^4$
(2) Player A can win after 40-15 with probability ${4\choose 1}\times p^4(1-p)$
(3) Player A can win after 40-30 with probability ${5\choose 2}\times p^4(1-p)^2$
(4) Player A can reach deuce, with probability ${6\choose 3}\times p^3(1-p)^3$, and once deuce is reached, player A will win with a probability of $\frac{p^2}{1-2p(1-p)}$ (http://www.austinrochford.com/posts/2013-04-25-probability-and-deuces-in-tennis.html)
Therefore the probability player A will win the game is given by,
$$P(Win) = p^4 + {4\choose 1}\cdot p^4(1-p) + {5\choose2}\cdot p^4(1-p)^2 + {6\choose 3}\cdot \frac{p^5(1-p)^3}{1-2p(1-p)}$$
Numerically, this problem is easy using Markov chains. It is a little bit tricky to create the transition matrix. I will paste some R code I created just for fun.
transition_matrix <- function(p_win) {
half_states <- c(0, 15, 30, 40)
state <- function(own, opp) paste0(half_states[own], '/', half_states[opp])
states <- c(state(expand.grid(1:4, 1:4)[, 1], expand.grid(1:4, 1:4)[, 2]), 'won', 'lost')
tm <- array(0, dim = c(length(states), length(states)))
colnames(tm) <- states
rownames(tm) <- states
rex <- '^([0-9]+)/([0-9]+)$'
for (s in states) {
if (grepl(rex, s)) {
own <- which(half_states == gsub(rex, '\\1', s))
opp <- which(half_states == gsub(rex, '\\2', s))
stopifnot(length(own) == 1, length(opp) == 1)
if (own < 4) {
j_win <- which(states == state(own + 1, opp))
} else {
if (opp < 4) {
j_win <- which(states == 'won')
} else {
j_win <- which(states == '40/30')
}
}
if (opp < 4) {
j_lost <- which(states == state(own, opp + 1))
} else {
if (own < 4) {
j_lost <- which(states == 'lost')
} else {
j_lost <- which(states == '30/40')
}
}
stopifnot(length(j_win) == 1, length(j_lost) == 1)
i <- which(rownames(tm) == s)
tm[i, j_win] <- p_win
tm[i, j_lost] <- 1 - p_win
} else {
if (s == 'won') tm['won', 'won'] <- 1 else tm['lost', 'lost'] <- 1
}
}
tm
}
For example, transition_matrix(0.3) is:
0/0 15/0 30/0 40/0 0/15 15/15 30/15 40/15 0/30 15/30 30/30 40/30 0/40 15/40 30/40 40/40 won lost
0/0 0 0.3 0.0 0.0 0.7 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
15/0 0 0.0 0.3 0.0 0.0 0.7 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
30/0 0 0.0 0.0 0.3 0.0 0.0 0.7 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
40/0 0 0.0 0.0 0.0 0.0 0.0 0.0 0.7 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.3 0.0
0/15 0 0.0 0.0 0.0 0.0 0.3 0.0 0.0 0.7 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
15/15 0 0.0 0.0 0.0 0.0 0.0 0.3 0.0 0.0 0.7 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
30/15 0 0.0 0.0 0.0 0.0 0.0 0.0 0.3 0.0 0.0 0.7 0.0 0.0 0.0 0.0 0.0 0.0 0.0
40/15 0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.7 0.0 0.0 0.0 0.0 0.3 0.0
0/30 0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.3 0.0 0.0 0.7 0.0 0.0 0.0 0.0 0.0
15/30 0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.3 0.0 0.0 0.7 0.0 0.0 0.0 0.0
30/30 0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.3 0.0 0.0 0.7 0.0 0.0 0.0
40/30 0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.7 0.3 0.0
0/40 0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.3 0.0 0.0 0.0 0.7
15/40 0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.3 0.0 0.0 0.7
30/40 0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.3 0.0 0.7
40/40 0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.3 0.0 0.0 0.7 0.0 0.0 0.0
won 0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0
lost 0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0
>
Now, what you are looking for is the result after playing infinite points. For practical purposes, multiplying the matrix by itself 20 times is like playing 2^20 == 1,048,576 points. You can expect all paths to be either 'won' or 'lost' after that many points.
for (time in 1:20) tm <- tm %*% tm
And what you get is the probability of winning for each possible initial state, not only 0/0.
> tm[1:16, c('won', 'lost')]
won lost
0/0 0.099211034 0.9007890
15/0 0.210981724 0.7890183
30/0 0.412168966 0.5878310
40/0 0.710224138 0.2897759
0/15 0.051309310 0.9486907
15/15 0.124758621 0.8752414
30/15 0.284431034 0.7155690
40/15 0.586034483 0.4139655
0/30 0.019831034 0.9801690
15/30 0.056327586 0.9436724
30/30 0.155172414 0.8448276
40/30 0.408620690 0.5913793
0/40 0.004189655 0.9958103
15/40 0.013965517 0.9860345
30/40 0.046551724 0.9534483
40/40 0.155172414 0.8448276
>
Additionally, you can solve the inverse problem using optimize()
sq_error <- function (p_win, p_win_game) {
tm <- transition_matrix(p_win)
for (time in 1:20) tm <- tm %*% tm
(tm['0/0', 'won'] - p_win_game)^2
}
p_win <- optimize(function(p) sq_error(p, p_win_game = 0.3), c(0, 1))$minimum
> p_win
[1] 0.4166342