Insights Puzzle

Solution: ‘How Many Half-Lives?’

When do negative results from a half-life experiment mean a theory is dead?

This month’s Insights puzzle was inspired by the negative results of an enormous proton decay experiment in Japan, as Natalie Wolchover recently reported. To help puzzle solvers develop an intuition for what a half-life is, I asked you to answer two simple questions in less than 30 seconds each. These two questions were similar in style to the Cognitive Reflection Test, originated by the decision theorist Shane Frederick. They require you to suppress the seemingly obvious but incorrect answer that first pops into mind and reflect more deeply to arrive at the correct, but less obvious, answer.

Question 1:

In a bucket filled with nutrients is a bunch of amoebas. Every day, their population doubles. If it takes 48 days for the bucket to be filled with amoebas, how long did it take for the bucket to become half full?

The answer that pops into the mind instantly is half of 48 — 24 days — because we are very familiar with linear process. However, what we are dealing with here is an exponential process. Since the population doubles every day, it had to have been half of the current value the day before, after 47 days. Exponential processes increase or decrease in proportion to their current value. One of the many mystical properties of the base of natural logarithms e (2.71828….) that makes it so special is that when you raise it to a variable x, the resulting function’s rate of change is exactly equal to its current value — its rise or fall is, in a sense, perfect.

Gerd drew attention to the fact that the illustration should have shown the bucket filled to a different height to make it half full. Since the bucket shown is not cylindrical but has sloping sides, this point is valid and presents a new puzzle, which readers are welcome to attempt (you can assume that the sides of the bucket slope upward at an angle of 100 degrees).

Question 2:

How many half-lives does a pound of radioactive material have?

This question can be answered from a few different perspectives. Mathematically, the answer is infinite, because there is always a finite probability that some atoms will not have decayed, no matter how much time elapses. From a practical point of view, I said it would be a very large number. Ralf B performed the calculation, and described it as follows (with a correction of the units as pointed out by Rational):

The pound contains 454 grams, or 454*6.022e23 atoms (assuming hydrogen). This is equal to about 2^88 atoms, so they will decay to one in 88 half-lives. For materials other than hydrogen, the number will be slightly smaller.

A monthly puzzle celebrating the sudden insights and unexpected twists of scientific problem solving. Your guide is Pradeep Mutalik, a medical research scientist at the Yale Center for Medical Informatics and a lifelong puzzle enthusiast.

Ralf B is correct in that this will be the expected number of half-lives. However, decay curves have long tails, so it is quite possible that there will be more half-lives than that. Consider just 2 atoms. A similar calculation says that you should have just 2 half-lives. But each atom has a 1-in-4 chance of not having decayed after 2 half-lives. So the chances that either one is still around after 2 half-lives are pretty good (7/16 or about 44 percent). It would take 8 half-lives to reduce this probability below 1 percent. For 288 atoms, a calculation using Wolfram Alpha indicates that it would take 95 half-lives to reach a point where the probability that no undecayed atoms exist exceeds 99 percent. So Ralf B’s point is well taken — in practice, the number of half-lives  is quite a bit larger than 2, but not quite as large as my words implied. I guess I hadn’t fully absorbed the lesson of Shane Frederick’s Cognitive Reflection test!

If you want to try the other questions from the Cognitive Reflection Test, remember, just 30 seconds per question:

  1. A bat and a ball cost $1.10 in total. The bat costs $1.00 more than the ball. How much does the ball cost? _____ cents
  2. If it takes 5 machines 5 minutes to make 5 widgets, how long would it take 100 machines to make 100 widgets? _____ minutes

Back to half-lives, after the first two simple reflection questions, I confess that I overcompensated by adding two questions that were far too complex for the amount of detail I provided. These questions therefore admitted to many possible answers, as I realized after looking at the various ways readers filled in the gaps. So here, instead of giving a unique solution when there isn’t one, I will discuss interesting insights, assumptions and techniques that readers brought to bear on Questions 3 and 4, all of which have some degree of validity considering the ambiguity in the questions.

Question 3:

Let us say you manage to obtain 5 atoms of the radioactive isotope of unobtainium, of Avatar fame. After exactly one year, 2 atoms have decayed. You want to figure out the half-life of the substance, and like a true scientist you seek a range that has a 95 percent probability of containing the true value. What is your range for the plausible half-life of unobtainium?

Question 4:

In the above scenario, for your second year, you obtain 30 atoms of the substance. You need to divide it into three portions — A, B and C — according to the following rules. When you inspect the portions after exactly one more year, you would like A to have exactly 6 undecayed atoms, B to have 7 or 8, and C to have 4 or 5. Remember, you don’t know the half-life exactly. How many atoms must be in A, B and C initially in order to maximize your chances of getting the precise results you want?

The simplest approach to these questions was adopted by Gerd, Luca U. and Tzara Duchamp, who calculated on the basis of the sample measurement alone to obtain an exponential decay constant of ln (5/3), giving a half-life of 1.36 years. For Question 4, the best way to split the groups then turns out to be 10, 13 and 7 atoms.

However, as John Pickering pointed out:

The way the question is phrased (seek a range that has a 95 percent probability of containing the true value) means that what is required is a Bayesian posterior (or credible) interval rather than simply a confidence interval (meaning that if the experiment was repeated multiple times the confidence interval calculated would contain the true value 95 percent of the time).

John’s point is worthy of deep consideration. When we make a statistical observation of a variable like half-life or height by examining a sample population, the measurement by itself cannot tell you the true value of the variable within the population. By random variation, our sample mean might be higher or lower than the true mean. But statisticians can calculate a 95 percent “confidence interval” as John explained: If the experiment were repeated multiple times, the confidence interval calculated would contain the true value 95 percent of the time. This can be done entirely based on the measurements from the samples. Each time you did the experiment, you would get a different confidence interval, but if you did it a hundred times, the true half-life or height would be within your stated intervals about 95 percent of the time.

However, this is subtly different from finding a range that has a 95 percent probability of containing the true value: In this case we cannot just rely on measurements from a sample, but we need some real world knowledge about how the true value is distributed (a Bayesian prior) in order to calculate the probability distribution of which we can exclude tails on the left and right sides representing a total probability of 5 percent. In the absence of clear direction in the questions on these points, readers came up with some ingenious extrapolations of their own. One possible prior assumption could have been to assume an equal chance of the half-life for equal periods of time (there is an equal chance of the half-life in the range 0-1 years as there is over 1-2 years), or that the decay rate is uniform, as Mark P assumed. The problem is that this obviously cannot be true for all the infinite periods of time that exist, or else the probability cannot add up to one. After some arbitrary point in time you have to assume that the chances of the half-life taper off to zero. Where that point is located, at, say, 5, 10 or 20 years, will affect the values of the true half-life and change the answers.

Nevertheless, several readers, such as John Pickering, Mark, Mark P and Joseph Fine, sketched out the general method very well: Figure out the probability distribution and then lop off enough of the left and right tails to leave 95 percent of the probability in the center. Here I liked Mark’s approach of not taking 2.5 percent of the probability distribution from the left and right tails (since they are asymmetric), but rather find a probability density that is equal on the left and right sides and encloses 95 percent of the probability distribution function. Very nice! Unfortunately, because of the differences in their assumptions, for which I take full responsibility, all the above readers obtained different ranges: Jon Pickering got 68 days to 27 years; Mark obtained 0.484 to 6.26 years; Mark P’s range was 0.9 to 6.3 years and Joseph Fine’s 0.755 to 3.106 years. All of you are winners. With such a small sample, a frequentist approach that might have given an accurate answer is untenable.

I find Mark P’s following speculation intriguing:

I kind of think that the obvious division of 10, 12, 8 or 10, 13, 7 is misguided. Those answers assume that the decay rate is exactly 0.4 per year. If it’s much higher or much lower than that — which is reasonable given the tiny sample size — using 10, 12, 8 or 10, 13, 7 is likely to overshoot or undershoot the mark by a wide margin.

That was the kind of situation I was trying to set up. Unfortunately, the different models will all yield different answers. In retrospect, I missed a chance to ask a far clearer question  more relevant to the negative proton decay result. Here is my do-over, for those who may still have the stamina to think about it. The answer is equal parts mathematical and psychological, just like the proton decay result.

Half-Life Do-Over Question:

Let us say you manage to obtain 5 atoms of the radioactive isotope of unobtainium, of Avatar fame. On theoretical grounds you have determined that its half-life is exactly one year. After exactly one year, no atoms have decayed. After how many years of obtaining this negative result will it be less than 5 percent likely that your theory is true? When will you be ready to admit that your theory is dead?

Thank you to all who attempted these ambiguous questions. I’m happy that some of you still had fun with it and contributed insights worth sharing.

The Quanta T-shirt for this month is still up for grabs. Do submit your contributions here.

Updated on February 10, 2017

The comments below have raised interesting issues:

  • Binomial vs. Poisson: Ashish used the binomial distribution in his calculation for the do-over question, whereas Joseph Fine used the Poisson distribution. Which one applies here? The Poisson distribution is an easy approximation of the binomial distribution that’s useful when the time interval of interest is a small fraction of the half-life — instances in which binomial calculations quickly become intractable. For this reason the Poisson distribution is the one used in serious scientific calculations such as the proton decay experiment. However, in the case of our question, the decay time is very short — on the order of the observation time — and therefore the binomial distribution, which is more accurate, should be used.
  • Truth and Probability: Ethaniel made a point about tightening the language regarding probability and the truth of a given theory, which is well taken. As he pointed out, the right question to answer is, “If you assume that your theory is true, how long does it take for the no-decay observation to fall under the 5 percent likelihood threshold?” The answer for the do-over question is 0.864 years, as determined by Ashish, and for the die-hards who want the likelihood to fall under the 5-sigma threshold (1 in 3.5 million, much used in physics experiments such as the proton decay calculation), it is 4.35 years, as determined by Ethaniel.
  • Letting Go of a Theory: Will some die-hards cling to their pet theories even longer than that? Mark explored the psychological nature of the question “When will you be ready to admit that your theory is dead?” by carrying out a number of elegant investigations that show how the strength of the initial belief affects how long you are willing to wait before losing hope. There are two famous anecdotes regarding this point.

The first concerns Albert Einstein, who in 1906 was confronted with an experimental result that seemed to go against his special theory of relativity while supporting some rival theories. According to his biographer Banesh Hoffmann, Einstein looked at the rival theories with “aesthetic disapproval” and confidently suggested that the experimenter might be mistaken. As we now know, Einstein’s self-belief was vindicated.

The second is a famous quote by Max Planck: “Science advances one funeral at a time.” Even in science, some die-hards never give up their pet theories. Sadly, unlike the other occasion, this was Einstein’s fate with quantum mechanics. It’s the new generation that embraces and advances new results.

For showing this mathematically, and for his earlier neat method capturing 95 percent of the probability distribution mentioned above, the Quanta T-shirt for this Insights puzzle goes to Mark. Congratulations!

 

View Reader Comments (13)

Leave a Comment

Reader CommentsLeave a Comment

  • DoOver: For 5 atoms the probability of any of them not decaying in their half-life is 1/32. So, the time in which they should decay with 95% confidence is .864 yrs. Once we have observed that none of the atoms decayed the process "restarts". There is no bearing of the atoms having not decayed in a given time on their future. So after any observation of all atoms we have .864 yrs to observe a decay before our confidence in the theory will drop to less than 5%.

  • Given that the half life of unobtainium is exactly 1 year, the original 5 atoms should have decreased to 5/4ths of an atom after 2 years. But none of the atoms have been observed to decay!

    One approach is to find the value of the decay rate 'lambda' used in the Poisson distribution.

    Using the formula ln(2) * ln(5/lambda) = 1 (half-life in years);

    The value of lambda is found to be 1.18145044 decays/year.

    Check on Half-life: ln(2) * ln(5/1.18145044) = 1.000000001 (approximately).

    Plugging that into the Poisson formula:
    http://keisan.casio.com/exec/system/1180573179

    The probability of x = 0 decays in 1 year = 0.30683337177 (more or less).

    And the probability of having less than 5% success of observing at least 1 atom decay is reached after 2.53564 years (half lives).

    So the theory/hypothesis that the half life is 1 year should be rejected after 3 years, if not before.

    Joseph Fine

  • "The pound contains 0.454 grams, or 0.454*6.022e23 atoms (assuming hydrogen). This is equal to about 2^88 atoms, so they will decay to one in 88 half-lives."

    No, a pound is 454 grammes, or .454 Kg.

    The answer is about 87.8 half lives but you need to watch your units!

  • > After how many years of obtaining this negative result will it be less than 5 percent likely that your theory is true?

    The p-value is the probability to observe a specific data (or something more extreme), given that the theory is true.
    It is thus NOT the probability that the theory is true, given that we have observed a specific data (here, no decay, which is already extremal).
    P(A|B)≠P(B|A), the most famous example of that obvious statement being the “false positive paradox”: The probability to be really sick given that one has been positively tested is not the probability to be positively tested given that you are really sick.
    So p-values (as used by Ashish above) do not give the likelihood of the theory as asked in the do-over question (“5 percent likely that your theory is true”).

    I may be wrong, but I think that it’s generally impossible to calculate the likelihood of a theory/hypothesis, at least in the frequentist approach.
    Following the frequentist approach I’ll imagine 1 billion of parallel universes and say that the “1-year half-life theory” is true in 5 % of them (i.e., in 50 million universes: So by the way no decay after 1 year will be observed in 1,562,500 universes on average) and that the theory is false in the other 95 %.
    The real problem is: How false it is?
    Among the 950 million universes where the theory is wrong, I can imagine:
    • a few million universes where the right theory is in fact “half-life = infinite”,
    • a few million universes where the right theory is in fact “half-life = 1000 years”,
    • a few million universes where the right theory is in fact “half-life = 1.001 years”,
    • a few million universes where the right theory is in fact “half-life = 0.999 years”,
    • a few million universes where the right theory is in fact “half-life = 1 day”,
    • a few million universes where the right theory is in fact “half-life = 1 second”,
    • etc.
    Whatever the right theories I chose for those 950 million parallel universes, the statement “the 1-year half-life theory has a likelihood of 5 %” stands, but obviously I will interpret very differently the no-decay observation if I chose all the 950 million universes to have infinite half-life or if I chose them all to have 1-second half-life.
    It where the Bayesian approach enters the game, by giving a prior state of those 950 million parallel universes.

    A way to answer the do-over question would be to assume that the classical misunderstanding of p-values apply, and that the right question to answer is “if you assume that your theory is true, how long does it take for the no-decay observation to fall under the 5 % likelihood threshold?”
    That question has been already answered by Ashish above: log(.05)/log(1/32) years ≈ 0.8644 years.
    Then for the second part of the question (“When will you be ready to admit that your theory is dead?”) I’ll take the standard one-sided 5-sigma threshold, giving: log(0.0000005733/2)/log(1/32) years ≈ 4.346 years.

    Another way to answer the do-over question would be to indeed search the likelihood of the theory, but I don’t see how to do that, even with a Bayesian approach…

  • The formula should have been:
    ln(2)/ln(5-lambda) = 1.0 (half-life in years)
    So (5-lambda) = 1.18145044
    lambda = 5 – 1.18145044 = 3.8185456 decays/yr.

    Plugging that into the Poisson calculator:
    http://keisan.casio.com/exec/system/1180573179
    The probability of 0 decays in 1 year = 0.021959(22893) and the theory is 'disproved' after 1 year and has less than 5% likelihood at t = ln(.05)/ln(.021959) = 0.784517 years.

    Joseph Fine

  • "… means that what is required is a Bayesian posterior (or credible) interval rather than simply a confidence interval ". Not completely true. The description given is precisely the definition of a "Tolerance Interval". Do we have two different right answers?

  • I'm disagree with the answer to the question 1.

    For me the good answer is ….NEVER !
    Why ? Cause amoebas are consuming nutriments ! And in fact the global volume never change… More amoebas less nutriments….

  • If the half life is theorized as exactly 1 year, lambda = ln(2) = 0.693147181.
    And the probability of zero decays in the first year(half-life( is 0.5.

    So the theory that half-life of 1 year can be "disproved" at the 95% level when (0.5)^x = .05 or at x = ln(0.05)/ln(0.5) = 4.321928 years (more or less).

    That should be enough iterations.

    Sorry for the previous errors.

    Joseph

  • For this problem I used a Bayesian approach again. It turns out that the Bayesian prior you use makes all the difference, so this probably explains Pradeep's comment about this being half psychological. Apologies in advance for the length of this comment…

    Define:

    t = half-life in years
    x = probability an atom does not decay in one year = 2^(-1/t)

    Therefore t = -1 / log2(x). The theoretical value is t0 = 1 year, so x0 = 1/2.

    Given n nuclei and probability x, then the probability that k nuclei do not decay in one year is given by the binomial distribution:

    P_binomial(k|n,x) = (n choose k) x^k (1-x)^(n-k)

    Bayes's Theorem gives us what we need, the PDF of x after a measurement that k nuclei do not decay:

    P(x|n,k) = P_binomial(k|n,x) P(x) / integral(P_binomial(k|n,x) P(x) dx, x=0..1)

    For clarity, let's call the PDF for x after m years P_m:

    P_m(x|n,k) = P_binomial(k|n,x) P_(m-1)(x) / integral(P_binomial(k|n,x) P_(m-1)(x) dx, x=0..1)

    The prior probability we use before starting our measurements is P_0(x), and this choice decides the whole rest of the question.

    If we ignore the theory and let the measurements speak for themselves (I am an experimentalist myself, so herein lies my bias!) then we should start with a flat or uninformative prior:

    P_(0, flat)(x) = 1

    In that case we get the following posterior probability after the first year, where we measure k=n=5 surviving nuclei:

    P_(1, flat)(x) = P_binomial(n|n,x) P_(0, flat)(x) / integral(P_binomial(n|n,x) P_(0, flat)(x) dx, x=0..1)
    = x^n / integral(x^n dx, x=0..1)
    = (n+1) x^n
    = 6 x^5

    I thought it would be appropriate to reject the hypothesis if the posterior probability of x being <= 1/2 is less than 5%, in other words that the integral(P(x), x=0..1/2) < 0.05. Another way of saying this is that the probability that the theory "x<=0.5" is true, given the data and the initial prior, is less than 5%. For the flat prior this integral is (1/2)^6 = 0.016 so we would reject the theory after one year since this number is < 0.05. This is what Ashish and Ethaniel found.

    Then I looked at some other priors. It sounds like we are very sure about our theory, so how do we incorporate that information? After doing a little google research, it looks like the typical way to parameterize priors for a binomial distribution is by using the beta distribution (https://en.wikipedia.org/wiki/Beta_distribution). I found a great explanation for how this works here:

    http://stats.stackexchange.com/questions/58564/help-me-understand-bayesian-prior-and-posterior-distributions/58792#58792

    Our initial prior must have a mean of 1/2, so we need to set alpha=beta. If we start with a prior like this:

    P_(0, alpha)(x) = Beta(alpha, alpha)(x)

    then after measuring k successes and n-k failures, the nice thing about the beta distribution is that the Bayesian prior for the next year becomes:

    P_(1, alpha)(x) = Beta(alpha+k, alpha+n-k)(x)

    and similarly for subsequent measurements.

    Since we continue measuring k=n=5, this becomes:

    P_(1, alpha)(x) = Beta(alpha+5, alpha)(x)
    P_(2, alpha)(x) = Beta(alpha+10, alpha)(x)
    P_(3, alpha)(x) = Beta(alpha+15, alpha)(x)

    Et cetera. As we increase alpha, the initial prior Beta(alpha, alpha) gets narrower and narrower around x=1/2. The standard deviation of Beta(alpha, alpha) is alpha/((2*alpha) sqrt(2*alpha+1)) = 1/(2*sqrt(2*alpha+1)). For alpha=1, the prior is flat. For alpha=2, the standard deviation is 1/(2*sqrt(5))=0.22. For alpha=3, the standard deviation is 1/(2*sqrt(7))=0.19, etc. Here are some plots of these distributions:

    https://en.wikipedia.org/wiki/Beta_distribution#/media/File:Beta_distribution_pdf.svg

    https://www.wolframalpha.com/input/?i=Plot%5BPDF%5BBetaDistribution%5B2,+2%5D,+x%5D,+PDF%5BBetaDistribution%5B3,+3%5D,+x%5D,+PDF%5BBetaDistribution%5B4,+4%5D,+x%5D,+x%3D0..1%5D

    The more confident we are in our theory, the narrower we should choose the prior to be. It seems that the choice of alpha is roughly like imagining we have already made a measurement of 2*alpha atoms and we measured 50% (alpha) of them decaying in a year.

    Going back to the flat prior from earlier, we could do the same calculation using beta distributions becase Beta(1,1) is the flat prior. So:

    P_(0, a=1)(x) = Beta(1,1)(x)

    After measuring 5 successes and 0 failures (decays), the posterior is

    P_(1, a=1)(x) = Beta(6,1)(x) = 6! / (5! 0!) x^5 = 6 x^5

    Which is what we already found. So that math checks out.

    Likewise starting with a prior

    P_(0, alpha)(x) = Beta(alpha, alpha)(x)

    After m years our posterior will be

    P_(m, alpha)(x) = Beta(alpha + 5m, alpha)(x)

    And the probability that x is below 1/2 (the p-value) comes from integrating this from 0 to 1/2.

    I wrote a Python script to perform this integral for different values of alpha and to determine the number of years until the p-value drops below 0.05. Here are the results:

    Beta distribution parameters: alpha = beta = 1
    Survival probability (x) +/- 1-sigma prior interval: 0.211–0.789
    Half-life (t) +/- 1-sigma prior interval: 0.446–2.920 years
    Minimum years until p-value < 0.05: 1 years

    Beta distribution parameters: alpha = beta = 3
    Survival probability (x) +/- 1-sigma prior interval: 0.311–0.689
    Half-life (t) +/- 1-sigma prior interval: 0.593–1.861 years
    Minimum years until p-value < 0.05: 2 years

    Beta distribution parameters: alpha = beta = 10
    Survival probability (x) +/- 1-sigma prior interval: 0.391–0.609
    Half-life (t) +/- 1-sigma prior interval: 0.738–1.398 years
    Minimum years until p-value < 0.05: 2 years

    Beta distribution parameters: alpha = beta = 30
    Survival probability (x) +/- 1-sigma prior interval: 0.436–0.564
    Half-life (t) +/- 1-sigma prior interval: 0.835–1.210 years
    Minimum years until p-value < 0.05: 3 years

    Beta distribution parameters: alpha = beta = 100
    Survival probability (x) +/- 1-sigma prior interval: 0.465–0.535
    Half-life (t) +/- 1-sigma prior interval: 0.905–1.109 years
    Minimum years until p-value < 0.05: 5 years

    Beta distribution parameters: alpha = beta = 1000
    Survival probability (x) +/- 1-sigma prior interval: 0.489–0.511
    Half-life (t) +/- 1-sigma prior interval: 0.968–1.033 years
    Minimum years until p-value < 0.05: 16 years

    Beta distribution parameters: alpha = beta = 10000
    Survival probability (x) +/- 1-sigma prior interval: 0.496–0.504
    Half-life (t) +/- 1-sigma prior interval: 0.990–1.010 years
    Minimum years until p-value < 0.05: 47 years

    All of this is to say quantitatively that depending on how much you believe the theory to really be true, you can wait as many years as you would like! As I said earlier, I am an experimentalist, and I am also impatient, so I would be satisfied with rejecting the theory after only 1 year of measurements. But someone very devoted to the theory might have to wait a generation until they are convinced…

  • Sorry to be dense, but isn't .454 * 6e23 closer to 2^78 not 2^88 ?

    log2(6*.454) + 23*log2(10) = 77.85

    (or just estimate it, we know that 10^3 is roughly 2^10, so 10^24 is roughly 2^80)

  • @Rational and Slipper.Mystery

    Correction
    As Rational pointed out, a pound is 454 grams or 0.454 kg. Thanks for pointing this out.

    The units in the quote by Ralf B are wrong, but the calculation of 2^88 atoms is correct. A more exact value, as Rational points out, is 2^87.8.

    We will correct this in the text.

  • For anybody interested in exponential growth/decay problems: You might find interesting a series of exercises I developed for exploring the exponential function both mathematically (with varying degrees of technicality, from experimental mathematics to algebraic proof) and empirically (using real world data on population, resource use etc.) To my knowledge there is very little literature that covers this ground (correct me if I'm wrong). Anybody who wants to use these materials for teaching is welcome to do so! (I'm not in teaching myself anymore).

    http://www.slideshare.net/amenning/exponential-growth-casestudies

Comments are closed.