Calculating the median in the St. Petersburg paradox

Question

I am studying a recreational probability problem (which from the comments here I discovered it has a name and long history). One way to address the paradox created by the problem is to study the median value instead of the expected value. I want to calculate the median value exactly (not only find bounds or asymptotic values). I have found a certain approach and I am stuck in a specific step. I present my analysis and I would like some help on that specific step.

[Note: Other solutions to the general problem are welcome (however after the revelation of the long history I found a lot of material) but what I really want is to know the answer to the sub-problem that my approach raises.]

The problem

We have the following game: I toss a coin as many times is needed to get tails. Then I count the number of consecutive heads that preceded (call it h) and I give you $2^h$ dollars. How much are you willing to pay to play such a game? In other words, what is the maximum buy-in for that game you are willing to pay? Note also, that we can play this game any amount of finite times (each time with you paying the buy-in).

A naive answer

One straightforward way to answer this is to calculate the expected value of one game. This should be the upper limit for the buy-in. The expected value is the infinite sum of the return of each case times the probability of each case. More specifically $$\sum_{i=0}^\infty (2^{i-1}\cdot\frac{1}{2^i}) = \sum_{i=0}^\infty \frac{1}{2} = \infty$$ This might seem counter-intuitive but it is true: Whatever constant and finite amount you bet per game, you are expected to win on the long run! Why is this so counter-intuitive though? Would you be willing to play this in practice with say 1000 dollars per game? The answer is no, because you would need an immensely large amount of games to actually win. So if we care about a more practical measure, the expected value is of no help. What we need is the median (or any other percentile value). If we know the median return for N games, we can at least know that if the buy-in is $\frac{median}{N}$, half of the possible cases you will lose and for half you will win. We will not know how much we will win or lose (we do have an upper bound on the losses though) but at least we know the chances to win or lose for a finite N number of games.

Finding the median

So how do you calculate the median return from N games (or more generally any ith percentile)?

If we play only one game (N=1) then it is trivial. The median is 1. For N=2 it starts getting more complicated. With probability 0.25 we'll get back 1+1, with 0.125 1+2, with 0.125 2+1. These 3 cases already bring us to a total of 0.5, so the median is 3 (and so the maximum bet is 1.5 per game). For any N, how do we enumerate all the cases and find the 50% point (or any i% point)? I realized that this is (partly) an ordering problem. We do not want just to enumerate random cases, we have to order them, starting from the case with the smallest possible return, then getting the one(s) with the next smallest return and so on. As we are doing this ordering we are adding the probabilities of these cases. When we reach 50% (or i%) we stop. The return value for that case is our median value (ith percentile value). The ordering is where I am stuck.

Sub-problem formulation

We can depict the possible space of returns with a matrix where the N columns are the N games and the infinite rows are the return for each game: $$\begin{array}{c} \text{row 1} \\ \text{row 2} \\ \text{row 3} \\ \vdots \\ \text{row i} \\ \vdots \end{array} \;\;\;\; \overbrace{\begin{array}{cccc} 1 & 1 & \cdots & 1 \\ 2 & 2 & \cdots & 2 \\ 4 & 4 & \cdots & 4 \\ \vdots & \vdots & \ddots & \vdots \\ 2^{i-1} & 2^{i-1} & \cdots & 2^{i-1} \\ \vdots & \vdots & & \vdots \end{array}}^N$$

A series of N games consists of picking values for each column (i.e., picking a game outcome for each game). The smallest possible total return is when all game outcomes are 1. So total return = N. The next possible one is when we get one outcome from the second row (total return N+1). The next smallest total return is N+2 (2 game outcomes from the second row). Notice though that for total return N+3 we have two "configurations": 1) cases where we have N-3 outcomes from the first row and 3 from the second row, OR 2) cases where we have N-1 outcomes from the 1st row and 1 outcome from the 3rd row! So ordering is not such an easy process.

Configurations vs. cases

Notice how I talked about "configurations" instead of individual cases. An individual case is a sequence of game outcomes (which are completely described by the game returns). For example a case of 4 games could be (1, 1, 16, 8) for a total return of 26. A configuration on the other hand is a more general construct which specifies how many outcomes we have from each row. A configuration completely determines the total return, but not the individual order that the outcomes happened. For example, the case given above is part of the configuration "2 outcomes from row 1, 1 outcome from row 4, 1 outcome from row 5". Cases (1,16,1,8) and (8,1,1,16) belong to the same configuration. From a configuration I can calculate how many distinct cases it has and what is the probability of each case. For example, for the configuration " $N_i$ outcomes from row i, $N_j$ from row j, $N_k$ from row k" we have:

The number of distinct cases is ${N\choose {N_i}}\cdot{{N-N_i}\choose{N_j}}\cdot{{N-N_i-N_j}\choose{N_k}}$

The probability for each of these cases is $2^{-(i\cdot N_i + j\cdot N_j + k\cdot N_k)}$

The total return value for any of these cases is $N_i \cdot 2^{i-1}+N_j \cdot 2^{j-1}+N_k \cdot 2^{k-1}$

The example above shows a configuration with 3 rows, just to get a taste of the complexity of the problem. I can generalise the formulas to find distinct cases, their probabilities and their total returns for any given configuration. The problem is ordering the configurations. Can we find an algorithm that orders and lists the configurations based on their total return value? Let's describe each configuration as a series of pairs {(x,i), (y,j), ...} where the first number of a pair denotes the row number and the second number of a pair denotes how many outcomes do we have from that row. For example, {(1,4), (3,1), (4,2)} means that we get 4 outcomes from row 1, 1 outcome from row 3, and 2 outcomes from row 4. This also means that we played 4 + 1 + 2 = 7 games. I manually computed the first terms of the ordered configurations list, for N games. I give the configuration(s) on the left and the total return on the right. Note that some total returns have more than one configurations that produce them.

$\begin{array}{ll} \text{Configurations} & \text{Total return} \\ \{(1,N)\} & N \\ \{(1,N-1),\; (2,1)\} & N+1 \\ \{(1,N-2),\; (2,2)\} & N+2 \\ \{(1,N-3),\; (2,3)\},\;\; \{(1,N-1),\; (3,1)\} & N+3 \\ \{(1,N-4),\; (2,4)\},\;\; \{(1,N-2),\; (2,1),\; (3,1)\} & N+4 \\ \{(1,N-5),\; (2,5)\},\;\; \{(1,N-3),\; (2,2),\; (3,1)\} & N+5 \\ \{(1,N-6),\; (2,6)\},\;\; \{(1,N-4),\; (2,3),\; (3,1)\},\;\; \{(1,N-2),\; (3,2)\} & N+6 \\ \end{array}$

If I can produce this order algorithmically then I will be able to calculate the median (or ith percentile) for any N.

I would also appreciate any help in formulating the problem in more accepted/mainstream terms. I believe that the formulation is valid and clear(?), but if we use a formulation from an established subfield maybe it will point to the solution too. Thanks!

Answer 1

It is easy enough to simulate this to get a reasonable estimate of the median values per game. For example using R

StPPmedian <- function(sampsize, numgames) {
   median( rowMeans( matrix(2^floor(-log2(runif(numgames*sampsize))), 
                      nrow=sampsize, ncol=numgames) ) )              }

you can get something like the following estimates of the median for 100,000 simulations each for various numbers of games:

> StPPmedian(100000,5)
[1] 2.4
> StPPmedian(100000,10)
[1] 2.9
> StPPmedian(100000,20)
[1] 3.4
> StPPmedian(100000,50)
[1] 4.08
> StPPmedian(100000,100)
[1] 4.59
> StPPmedian(100000,200)
[1] 5.115
> StPPmedian(100000,500)
[1] 5.78

and drawing these estimates of the median against the logarithm of the number of games certainly suggests some kind of logarithmic relationship is plausible, possibly with the median value per game close to $\log_4(N) + O(1)$ as also suggested in the following simulations

> StPPmedian(100000,4)
[1] 2.25
> StPPmedian(100000,16)
[1] 3.3125
> StPPmedian(100000,64)
[1] 4.296875
> StPPmedian(100000,256)
[1] 5.289062

Added December 2015:

In the comments I said that empirically the median appeared to be $\log_4(N)+O(1)$ and that to one decimal place the $O(1)$ term seems to be about $1.3$ for large $N$. A.S. suggested that the value of $O(1)$ depends on along which sub-sequence you take the limit. The following chart may demonstrate this: there seem to be visual patterns, combined with the noise resulting from simulation.

maxn <- 512 ; n <- 1:maxn ; meds <- rep(NA,maxn)
for (i in n){ meds[i] <- StPPmedian(1000000,i) } 
plot(meds-log(n,4), ylim=c(1.2,1.4))

Answer 2

Feller was able to determine the median asymptotically, namely $N \log_2 N$, or $\log_2 N$ per game. He actually proved an even stronger concentration bound. This result can also be found in his "Introduction to Probability".

Economists are still coming up with new interpretations of this problem, known as the St. Petersburg Paradox, after the city in which Euler, who invented it, was working (the name might be due to Feller). See, for example, this article, which also suggests the same value.

Answer 3

You can estimate the median by considering events with probability more than $1/2$. For example, let $X_t$ be the event "exactly one game reached stage $t$". The probability that this happens is $$ p_t = N 2^{-t} (1-2^{-t})^{N-1}. $$ For $t \approx \log_2 N$, we have $p_t \approx 1/2$. This means that with probability more than $1/2$, the total outcome is at least $2^t \approx N$. So the median is at least $\approx 2N$.

By considering the event "no game reached stage $t$", we get (in the same way) that the median is at most $N^2$. This estimate is probably very bad.

Both these values are not normalized - to get the value per round, divide by $N$.

Probably by devising better events you can get closer to the actual median. Your enumeration approach, on the other hand, seems to only work for small $N$, even if you did find a way to enumerate the probabilities in increasing order of total outcome. What I'm suggesting is very similar to this enumeration, only you don't have to be entirely accurate; the more accurate you are, the closer to the median you'll get.

Answer 4

I had to write an ugly recursive algorithm in C to calculate the exact median for this game. Please find the first 30 exact median values in the graph at page 62 philpapers.org/archive/ERGTEO.pdf

Answer 5

The usual resolution to this paradox is to consider that whoever is paying you off cannot pay off an infinite amount of money. So you terminate your sum at the point you have all the bank's money. If that is \$1B, they can only pay 30 rounds, so the fair value is \$30. As it is logarithmic, it doesn't vary much as you change the estimate of what they can pay.

Added: You actually don't care about the order of results, just what they are. So in $N$ games, your chance of getting $N$ is just $2^{-N}$, your chance of getting $N+1$ is $N2^{-N}$, your chance of $N+2$ is $(N+N(N-1)/2)2^{-N}$ and so on. The generating function will keep track of this. If the cutoff is $2^m$, the generating function is $\left(\frac{x^{2^m}}{2^{m+1}} + \sum_{i=0}^m \frac{x^{2^i}}{2^{i+1}}\right)^N$ where the coefficient of $x^p$ gives the chance of winning $p$