Stable Complementary Log-Log Probabilities

Why This Link Needs Careful Computation

For the complementary log-log link, the event probability approaches zero gradually as the linear predictor becomes very negative, and approaches one as it becomes very positive. These limits are mathematically straightforward, but direct floating-point calculations can lose information near either boundary. For example, computing an event probability as one minus the non-event probability can return zero when the true event probability is small; computing a log probability after rounding a probability to zero can return negative infinity even when the log probability itself is representable.

This tutorial focuses on calculating the event probability, the non-event probability, and their logarithms. These quantities enter Bernoulli and binomial likelihood calculations, so preserving them matters even when a displayed probability has rounded to zero or one.

Definition: Under the complementary log-log link, the linear predictor $\eta$ and event probability $\mu$, with $0<\mu<1$, satisfy $$ \log\{-\log(1-\mu)\}=\eta. $$ The inverse link gives the event probability $p(\eta)=1-\exp\{-\exp(\eta)\}$ and the non-event probability $q(\eta)=\exp\{-\exp(\eta)\}$. Thus $p(\eta)+q(\eta)=1$.

The defining distinction is between computing a probability and computing its logarithm. A probability may be too close to zero or one to be represented accurately, while its logarithm may still be a useful, finite number. A stable implementation should therefore calculate the quantity actually needed, rather than calculate a rounded probability and then take its logarithm.

Stable Forms for the Probabilities

Let $t=\exp(\eta)$, which is positive in exact arithmetic. The identity $1-\exp(-t)=-\operatorname{expm1}(-t)$ gives a stable formula for the event probability:

p(\eta)=-\operatorname{expm1}\{-\exp(\eta)\}, \qquad q(\eta)=\exp\{-\exp(\eta)\}.

Here $\operatorname{expm1}(x)$ computes $\exp(x)-1$ while retaining accuracy when $x$ is close to zero. When $\eta$ is very negative, $t$ is small and $\exp(-t)$ is close to one. Subtracting that rounded value from one can destroy the small difference that is $p$. By contrast, $-\operatorname{expm1}(-t)$ computes that difference directly.

The non-event probability does not require subtraction: its direct expression is an exponential. In particular, do not compute it as $1-p$ when $p$ is close to one. For large $\eta$, $q$ may be smaller than the smallest representable positive number and round to zero. That is a limit of the floating-point representation, not a reason to compute $q$ by subtracting two nearly equal numbers.

Worked Example: A Very Small Event Probability

Take $\eta=-20$. Then $t=\exp(-20)\approx 2.061153622\times10^{-9}$. The stable event-probability calculation is $$ p=-\operatorname{expm1}(-2.061153622\times10^{-9}) \approx 2.061153620\times10^{-9}. $$

The non-event probability is $$ q=\exp(-2.061153622\times10^{-9}) \approx 0.9999999979388464. $$

In exact arithmetic, $1-q=p$. But subtracting the displayed $q$ from one gives about $2.0611536\times10^{-9}$, with fewer reliable digits than the direct expm1 calculation; at more extreme predictors, rounding $q$ to one can make the subtraction return zero. The stable formula calculates the small event probability directly.

Stable Log Probabilities

For a Bernoulli observation $Y\in\{0,1\}$, the log-likelihood contribution is $\log p$ if $Y=1$, and $\log q$ if $Y=0$. The non-event log probability has a particularly simple exact form:

\log q(\eta)=-\exp(\eta).

It is important not to obtain this value by taking the logarithm of a computed $q$. If $q$ has underflowed to zero, $\log(q)$ is negative infinity, although $-\exp(\eta)$ may still be representable and finite.

For the event, the direct log formula is $\log\{-\operatorname{expm1}(-t)\}$, where $t=\exp(\eta)$. This is appropriate when $t$ can be computed and the event probability is not lost to underflow. At very negative $\eta$, however, $\exp(\eta)$ can itself underflow. In that regime use the factorization

p=t\,h(t),\qquad h(t)=\frac{1-\exp(-t)}{t} =\frac{-\operatorname{expm1}(-t)}{t}, \qquad \log p=\eta+\log h(t).

As $t$ tends to zero, $h(t)$ tends to one, so $\log h(t)$ tends to zero and $\log p$ is close to $\eta$. A useful expansion, obtained from the Taylor series of the exponential and logarithm, is $$ \log h(t)=-\frac{t}{2}+\frac{t^2}{24}+O(t^4) \quad\text{as }t\to0. $$ For sufficiently small $t$, this correction can be evaluated as a short series using $\operatorname{log1p}(u)$, which computes $\log(1+u)$ accurately for small $u$. For example, $h(t)=1+u$ with $u=-t/2+t^2/6+O(t^3)$; evaluating $\operatorname{log1p}(u)$ avoids taking the logarithm of a rounded value of $h(t)$ that may equal one. In an implementation, the series and its switching threshold should be chosen for the precision and error tolerance in use.

Worked Example: A Log Event Probability Beyond Probability Range

Take $\eta=-800$. In ordinary double-precision arithmetic, $\exp(-800)$ underflows to zero. Applying the probability formula to that computed zero would give $p=0$, and then $\log(p)$ would be negative infinity.

Mathematically, $t=\exp(-800)$ is positive but extraordinarily small. From the factorization above, $$ \log p=-800+\log h(t). $$ Since $h(t)=1-t/2+O(t^2)$, the correction $\log h(t)$ is approximately $-t/2$, far smaller than the precision available when adding it to $-800$. Therefore $\log p$ rounds to $-800$, a finite and useful value, even though $p$ itself is too small to represent.

This result does not assert that the probability equals $\exp(-800)$ in floating-point arithmetic. It shows why a log-probability routine should use the predictor and the small-$t$ limit rather than take the logarithm of an underflowed probability.

Derivations and Boundary Behavior

Result: For finite real $\eta$, the complementary log-log inverse link gives $0<p(\eta)<1$ and $0<q(\eta)<1$. Moreover, as $\eta\to-\infty$, $\log p(\eta)-\eta\to0$; as $\eta\to+\infty$, $\log q(\eta)=-\exp(\eta)$.

Proof. For finite $\eta$, $t=\exp(\eta)>0$. Hence $0<\exp(-t)<1$, so $q=\exp(-t)$ is strictly between zero and one. Since $p=1-q$, $p$ is also strictly between zero and one.

For the negative-predictor limit, write $p=t h(t)$, where $h(t)=(1-\exp(-t))/t$. The derivative of $\exp(-t)$ at zero is $-1$, so $(1-\exp(-t))/t\to1$ as $t\to0$. Continuity of the logarithm at one then gives $\log h(t)\to0$. Since $\log t=\eta$, $\log p-\eta=\log h(t)\to0$. The last identity, $\log q=-\exp(\eta)$, follows exactly from $q=\exp\{-\exp(\eta)\}$; no approximation is involved. $\square$

Result: For $t>0$, $p=1-\exp(-t)=-\operatorname{expm1}(-t)$, and $q=\exp(-t)$. For $t$ near zero, $p=t+O(t^2)$, while $\log p=\log t+O(t)$.

Derivation. The expm1 identity follows from $\operatorname{expm1}(-t)=\exp(-t)-1$, so negating both sides gives the stated event probability. For the small-$t$ behavior, Taylor expansion gives $\exp(-t)=1-t+t^2/2+O(t^3)$. Therefore $p=t-t^2/2+O(t^3)=t+O(t^2)$. Equivalently, $p=t h(t)$ with $h(t)\to1$, and taking logarithms yields $\log p=\log t+\log h(t)=\log t+O(t)$. These identities explain why a small event probability is approximately $\exp(\eta)$, while its log probability is approximately $\eta$. $\square$

For the other boundary, as $\eta$ grows, $q=\exp\{-\exp(\eta)\}$ decreases rapidly toward zero. The event probability approaches one, and its logarithm approaches zero. When $q$ is small, one may calculate $\log p$ as $\operatorname{log1p}(-q)$, which avoids loss from computing $\log(1-q)$ by a less suitable route. When $q$ is so small that it rounds to zero, $\log p$ rounds to zero, as it should at ordinary floating-point precision.

Worked Example: A Large Predictor and a Finite Log Non-Event Probability

Take $\eta=7$, so $t=\exp(7)\approx1096.633158$. The non-event probability is $$ q=\exp(-1096.633158), $$ which underflows to zero in standard double precision. The event probability calculated by the expm1 expression rounds to one, because the difference of $q$ from one is much smaller than the spacing of representable numbers near one.

The log non-event probability is different: $$ \log q=-\exp(7)\approx-1096.633158. $$ This is a finite representable number. Thus a non-event observation at this predictor has a very large negative, but finite, log-likelihood contribution. Computing $\log(q)$ after $q$ has underflowed would instead produce negative infinity and incorrectly discard the available finite value.

For an event observation, $\log p$ is so close to zero that it rounds to zero at this precision. A stable calculation can use $\operatorname{log1p}(-q)$ when $q$ is representable; if $q$ has rounded to zero, returning zero for $\log p$ is consistent with the floating-point resolution.

A Practical Calculation Strategy

A stable routine should distinguish the range where an intermediate is useful from the range where it overflows or underflows. The exact switching thresholds depend on the floating-point type and software library, so they should not be treated as universal constants. A sound general strategy is:

Keep the predictor.
Use $\eta$ directly for limiting branches; do not first calculate a probability that may round to a boundary.

For the event probability, use expm1.
When $t=\exp(\eta)$ is safely computable, evaluate $-\operatorname{expm1}(-t)$, not $1-\exp(-t)$.

For the event log probability, protect the small-$t$ case.
Use $\log\{-\operatorname{expm1}(-t)\}$ when appropriate; for very small $t$, use $\eta+\log h(t)$ with a stable small-$t$ evaluation, or its justified series limit.

For the non-event log probability, avoid taking a log of a probability.
Evaluate $-\exp(\eta)$ directly, guarding the exponential if its magnitude exceeds the representable range.

Interpret rounded boundaries correctly.
A computed probability of zero or one may be the unavoidable rounded result; retain a finite log probability whenever its direct formula is representable.

For example, when $\eta$ is so positive that $\exp(\eta)$ exceeds the largest finite floating-point number, the mathematical value of $\log q=-\exp(\eta)$ is below the smallest representable finite value. Returning negative infinity in that case reflects a genuine range limitation for that log probability. It is not the same issue as obtaining negative infinity by taking $\log(0)$ after avoidable underflow of $q$.

Common Mistake: Subtracting or Logging Rounded Probabilities

A frequent mistake is to calculate $p=1-q$ and then use that result for every subsequent calculation. The reasoning is that $p+q=1$, so either probability can be recovered from the other. That identity is exact, but floating-point subtraction is not equally accurate in all ranges. When $q$ is close to one, the subtraction $1-q$ loses relative accuracy in the small event probability. When $p$ rounds to one, $1-p$ can lose the non-event probability entirely. Taking a logarithm afterward cannot recover information already lost.

The correction is to calculate each target quantity from its own stable expression: use $-\operatorname{expm1}(-\exp(\eta))$ for the event probability, $-\exp(\eta)$ for the log non-event probability, and a small-$t$ factorization for the log event probability. Use $\operatorname{log1p}$ when a logarithm has the form $\log(1+u)$ with small $u$; do not assume that an algebraically equivalent expression has the same numerical accuracy.

Key takeaway: The complementary log-log probabilities have simple exact formulas, but stable computation depends on choosing the right form for the quantity and predictor range. Use expm1 for a small event probability, log1p or a justified small-argument expansion for its logarithm, and the direct log formula for the non-event probability.

Check Your Understanding

Use the formulas and reasoning above to answer the following questions.

For $\eta=-12$, which expression should be used to compute the event probability, and why is $1-\exp\{-\exp(\eta)\}$ less desirable?
Explain why a non-event probability may round to zero while its log probability remains finite.
Derive the approximation $\log p(\eta)\approx\eta$ as $\eta\to-\infty$ from the factorization $p=t h(t)$.
For a large positive predictor, why is calculating $q$ as $1-p$ potentially unreliable?
When might $\operatorname{log1p}(-q)$ be a suitable calculation for $\log p$, and what should be considered if $q$ has already underflowed?