I’ve been prescribed to take 1.5 pills of a certain medication every day for 10 days, so I have a bottle with 15 pills. Each morning, I take two pills out of the bottle at random.

On the first morning, these are guaranteed to be two full pills. I consume one of them, split the other in half using a precision blade, consume half of that second pill, and place the remaining half back into the bottle.

On subsequent mornings when I take out two pills, there are three possibilities:

- I get two full pills. As on the first morning, I split one and place the unused half back into the bottle.
- I get one full pill and one half-pill, both of which I consume.
- I get two half-pills. In this case, I take out another pill at random. If it’s a half-pill, then I consume all three halves. But if it’s a full pill, I split it and place the unused half back in the bottle.

Assume that each pill — whether it is a full pill or a half-pill — is equally likely to be taken out of the bottle.

On the 10th day, I again take out two pills and consume them. In a rush, I immediately throw the bottle in the trash before bothering to check whether I had just consumed full pills or half-pills. What’s the probability that I took the full dosage, meaning I don’t have to dig through the trash for a remaining half-pill?

Problems of this sort may seem tedious at first glance, but there is a simple and systematic way to approach them. We can encode the pills in the bottle using *polynomials*. Specifically, we’ll write $x^n y^m$ when there are $n$ full pills and $m$ half-pills. Since there is an element of randomness with how the pills are picked, we will use a weighted sum of such polynomials to describe a mixture or possible states. For example, we could write:

\[

\tfrac{2}{3} x y + \tfrac{1}{3} y^3

\]to describe a scenario where there is a $\frac{2}{3}$ probability that the bottle contains one full pill and one half-pill, and a $\frac{1}{3}$ probability that the bottle contains three half-pills. Note that the coefficients sum to $1$, so exactly one of the two scenarios must occur.

So what happens when we take our daily dose? We can consider each scenario in turn.

- With probability $\frac{\binom{n}{2}}{\binom{n+m}{2}}$, we will pick two full pills. In this case, we cut one in half and return half to the bottle. This means the bottle now contains $x^{n-2}y^{m+1}$.
- With probability $\frac{nm}{\binom{n+m}{2}}$, we will pick one full pill and one half-pill. In this case, we have our exact dose and the bottle now contains $x^{n-1}y^{m-1}$.
- With probability $\frac{\binom{m}{2}}{\binom{n+m}{2}}$, we will pick two half-pills. In this case, we must go back into the bottle and pick one more pill. There are two sub-cases here:
- With probability $\frac{n}{n+m-2}$, we pick one of the full pills. So we split it in half, and return half to the bottle. The jar has now lost a total of one full pill and one half-pill, so it contains $x^{n-1}y^{m-1}$.
- With probability $\frac{m-2}{n+m-2}$, we pick one of the half-pills. So we have our exact dose and the bottle now contains $x^n y^{m-3}$.

Combining all the possible cases, we find that every day, the polynomial describing the contents of the bottle change according to the mapping

\begin{multline}

x^n y^m \mapsto \tfrac{n (n-1)}{(m+n) (m+n-1)}x^{n-2} y^{m+1} + \tfrac{m n (3 m+2 n-5)}{(m+n) (m+n-1) (m+n-2)} x^{n-1} y^{m-1} \\

+ \tfrac{m (m-1) (m-2)}{(m+n) (m+n-1) (m+n-2)} x^n y^{m-3}

\end{multline}

There are two nice things about this formulation. First, edge cases are taken care of automatically. For example, if the bottle contains $n=2$ full pills and $m=2$ half-pills, the third term ($x^n y^{m-3}$) in the expansion is impossible, since the bottle does not contain three half-pills. But the associated coefficient is zero, so this term cancels automatically. Similarly, if $m=0$ (the bottle only contains full pills), two of the coefficients vanish and the third becomes $1$, so we obtain $x^n \mapsto x^{n-2} y$.

Second, when the bottle contains a mixture of possible states (described by a polynomial whose coefficients sum to 1), we can apply the above transformation to each term of the polynomial separately, and the final result will again be a mixture of possible states. This follows because the sum of the coefficients of the mapping above sum to 1 as well.

To compute the mixture of states after 10 days, I wrote a Mathematica script that starts at $x^{15}$ (15 full pills) and applies the above transformation recursively. Here is the Mathematica script:

P = x^15; (* initial polynomial *)
Print[P];
For[iter = 0, iter < 9, iter++,
Q = 0;
ML = MonomialList[P];
For[i = 1, i <= Length[ML], i++,
term = ML[[i]];
coef = term /. {x -> 1, y -> 1};
n = Exponent[ term, x];
m = Exponent[term, y];
Q = Q + coef ((m(m-1)(m-2))/((n+m)(n+m-1)(n+m-2)) x^n y^(m-3)
+ (n m(2n+3m-5))/((n+m)(n+m-1)(n+m-2)) x^(n-1) y^(m-1)
+ (n(n-1))/((n+m)(n+m-1)) x^(n-2) y^(m+1));
];
P = Q // Expand;
Print[P];
]

and here is the output from the script:

\begin{gather}

x^{15} \\

x^{13} y\\

\tfrac{1}{7}x^{12}+\tfrac{6 }{7}x^{11} y^2\\

\tfrac{36 }{91}x^{10} y+\tfrac{55 }{91}x^9 y^3\\

\tfrac{23 }{308}x^9+\tfrac{2385}{4004} x^8 y^2+\tfrac{30}{91} x^7 y^4\\

\tfrac{4 }{13}x^7 y+\tfrac{81}{143} x^6 y^3+\tfrac{18}{143} x^5 y^5\\

\tfrac{335}{4004} x^6+\tfrac{87}{154} x^5 y^2+\tfrac{185}{572} x^4 y^4+\tfrac{4}{143} x^3 y^6\\

\tfrac{22573}{56056} x^4 y+\tfrac{42605}{84084} x^3 y^3+\tfrac{101}{1144} x^2 y^5+\tfrac{1}{429}x y^7\\

\tfrac{44783}{240240} x^3+\tfrac{362603}{560560} x^2 y^2+\tfrac{27185}{168168} x y^4 + \tfrac{61}{12012} y^6\\

\tfrac{80529}{101920} x y+\tfrac{21391}{101920} y^3\\

\end{gather}

Therefore, on the $10^\text{th}$ day, there is a $\tfrac{21391}{101920} \approx 20.988\%$ chance the bottle contains three half-pills, and a $\tfrac{80529}{101920} \approx 79.012\%$ chance the bottle contains one full pill and one half-pill. This is the probability that the correct dosage was taken, which is what is asked in the original problem statement.

### Large-pill limit

We saw that on the final day, the probability that the bottle contains one full pill and one half-pill is about $79\%$. A natural question to ask is: what happens to this probability as we vary the number of starting pills? If we start with only 3 pills, then we *always* have one full pill and one half-pill on the next day. So the probability is $100\%$. Here is a plot of what happens as we vary the number of initial pills:

It appears that as we increase the number of initial pills, the probability of having one full pill and one half-pill on the last day decreases. I don’t have a great explanation for why this happens; if you made it this far and have any ideas, write a comment; I would love to hear it!

### Generating functions?

I am well aware that the polynomials I used look a lot like generating functions. I looked for a nice way to express the pill-taking transformation as a differential operator acting on the generating function, but I was thwarted by the denominators. The numerators are easy, e.g.

\[

\tfrac{\mathrm{d}^3}{\mathrm{d}y^3}(x^n y^m) = m(m-1)(m-2) x^n y^{m-3}

\]but I couldn’t find a nice way to deal with the denominators. For the case above, the denominator should be $(n+m)(n+m-1)(n+m-2)$. If anybody found a way to solve this using the full generating function machinery, I would love to hear about it!

Interesting. Thank you for this.

With respect to why the probability of ending up with all half-pills increases with the number of starting pills: I don’t have a rigorously mathematical explanation, but my sense is that you’re more likely to end up with a lot of half-pills if you ‘start’ with a lot of half pills. In the Riddler scenario, you start 10 days out with 15 whole pills, and 0 half-pills. If, instead, you started 12 days out with 18 whole pills, then when you got to 10 days out (the Riddler starting point), you’d be pretty likely to have one or two half-pills in the jar, rather than all whole pills, and therefore (I think) you’d be slightly more likely to finish with all half pills on the last day.

Again, I haven’t crunched any number on this, it’s just my intuition.

Nice solution! I also solved the problem using Mathematica, but using a different technique, and I reproduced your results. I was able to do an exact calculation using 7500 pills and the probability of ending up with a whole pill plus a half pill was 43.1446%. I tired applying Wynn’s epsilon method, for whatever that’s worth (probably not much here), and got a probability for an infinite number of pills of 28.5%. Relying entirely on intuition, I suspect the probability goes very slowly to zero as the number of pills goes to infinity.

Bravo, Laurent! I’ve played around quite a bit with the question of whether the probability of a full pill and a half pill converges to something as the number of starting pills increases. I haven’t been able to rigorously answer the question, but I do suspect it doesn’t converge and instead very slowly goes to zero. I don’t have any good intuition as to why.

But I should add that Joe’s suggestion makes perfect sense to me. Increasing the number of starting pills means there will be more half-pills around later, which increases likelihood of all half-pills at the end.

Let W be the number of whole pills, H the number of half pills, and T=W+H the total. Let w=W/T. Consider the situation W>>1 and H>>1. The expectation value for the one-day change in W is then ΔW = -2w^2 – (1-w^2-(1-w)^3). Call the right-hand side f(w), and assume we can replace w by its expectation value. We also have Δ(2W+H)=Δ(W+T)=-3 exactly. Manipulating all this eventually yields g(w)Δw = ΔT/T, where g=(3-f)/(f-w(3-f)). Treating this as a differential equation and integrating from a start of w=1, T=T0 to an end of w=wf, T=1, we find (for very large T0) that wf = 3/[ln(T0) + c], where c=3+π/√3. This doesn’t seem to fit the numbers very well, and of course falls apart at least at the beginning (where H is not large) and end (where neither W nor H is large), and maybe in the middle if treating the expectation value as exact is misleading. But I still feel like it captures the essence of what’s happening.

I was mistakenly identifying wf (the expectation value) with pf (the probability). Since the final state is either W=H=1, w=1/2, with probability pf, or W=0,H=3,w=0, with probability 1-pf, we have wf = (1/2)pf + (0)(1-pf) = pf/2. Correcting this factor of two makes the numbers work much better (in particular Emily Boyajian @Emily8191 found pf=0.31 for T0=10^6). So my claim is that pf ~ 6/log(T0) in the limit of large T0.

Nice solution. This is definitely more efficient than banging out the probabilities in excel…

I figure that the reason for the probability of the 1W1H state becoming less likely for larger initial populations ultimately comes down to this:

You’ll never put a whole pill back into the bottle, so the only way to reach the 1W1H state on day n is to avoid drawing a specific whole pill for the entire (n-1) days beforehand.