# Sampling from distribution parametrized by $\log(1-p)$ given $\log (p)$

The context for the problem is that I'm working with a modified genetic algorithm where the fitness score of each chromosome is given by a log-probability $$\log(p) \in (-\infty, 0)$$. Thus, I have a categorical distribution parameterized by (unnormalized) log-probabilities. These scores are large in magnitude, so I want to be able to sample from this distribution without needing to first exponentiate and normalize, as exponentiating would lead to numerical over/underflow errors during program execution. I know that I can sample proportional to $$p$$ without exponentiation by using the Gumbel-max trick.

However, what's giving me trouble is how, given this set of log-probabilities $$\log(p)$$ , I would be able to sample from a distribution parameterized by $$\log(1 - p)$$ . The motivation for this is that I want to "kill off" lowest-scoring chromosomes in the population, hence driving up overall population fitness. It is not possible to directly convert $$\log(p)$$ to $$\log(1 - p)$$ by logarithm rules, and as such I'm struggling to find some kind of transformation that will allow me to do this.

How should I go about thinking of a solution to this? Is this sort of thing even possible?

• What distribution are you sampling from?
– Sycorax
Jan 25 at 20:48
• @Sycorax Thanks for the clarifying question -- it is not a standard distribution, and I've just added some edits to my question for some more context. Jan 27 at 20:42

It is not possible to directly convert $$\log(p)$$ to $$\log(1−p)$$ by logarithm rules...

Au contraire mon ami, it is indeed possible. Since $$p \geqslant 0$$ you have $$p = \exp(\log p)$$ (allowing computation for infinities) so you can express the required quantity as:

\begin{align} \log(1-p) &= \log(1 - \exp(\log p)) \\[6pt] &= \log(1 - \exp(-(-\log p))) \\[6pt] &= \text{log1mexp}(-\log p). \\[6pt] \end{align}

Thus, if you already have the value $$a = \log p$$ you can get the desired quantity as:

$$\log(1-p) = \text{log1mexp}(-a).$$

This transformation can be implemented in R using the log1mexp function in the VGAM package (also see Mächler 2012 for some discussion of computational issues for this function). Incidentally, this function is closely related to the hyperbolic functions in trigonometry.

• Hi Ben, thanks for your answer and the link to the paper. This seems like a sensible approach -- however, the issue is that it exponentiates the log-probabilities during computation. I'm dealing with fairly large log-probabilities (i.e. abs($\log (p) ) > 1000$ ), and I begin hitting numerical overflow errors for abs($\log (p) ) > 750$. Ideally, I wouldn't have to exponentiate at all (as is the case for the Gumbel Max Trick). Jan 28 at 20:54