In $$ p(Y=C|\mathbf{x}) = \frac{p(\mathbf{x}|Y=C)p(Y=C)}{~\sum_{k=1}^{|C|}{}p(\mathbf{x}|Y=C_k)p(Y=C_k)} $$
both the denominator and the numerator can become very small, typically because the $p(x_i \vert C_k)$ can be close to 0 and we multiply many of them with each other. To prevent underflows, one can simply take the log of the numerator, but one needs to use the log-sum-exp trick for the denominator.
More specifically, in order to prevent underflows:
If we only care about knowing which class $(\hat{y})$ the input $(\mathbf{x}=x_1, \dots, x_n)$ most likely belongs to with the maximum a posteriori (MAP) decision rule, we don't have to apply the log-sum-exp trick, since we don't have to compute the denominator in that case. For the numerator one can simply take the log to prevent underflows: $log \left( p(\mathbf{x}|Y=C)p(Y=C) \right) $. More specifically:
$$\hat{y} = \underset{k \in \{1, \dots, |C|\}}{\operatorname{argmax}}p(C_k \vert x_1, \dots, x_n) = \underset{k \in \{1, \dots, |C|\}}{\operatorname{argmax}} \ p(C_k) \displaystyle\prod_{i=1}^n p(x_i \vert C_k)$$
which becomes after taking the log:
$$ \begin{align} \hat{y} &= \underset{k \in \{1, \dots, |C|\}}{\operatorname{argmax}} \log \left( p(C_k \vert x_1, \dots, x_n) \right)\\ &= \underset{k \in \{1, \dots, |C|\}}{\operatorname{argmax}} \log \left( \ p(C_k) \displaystyle\prod_{i=1}^n p(x_i \vert C_k) \right) \\ &= \underset{k \in \{1, \dots, |C|\}}{\operatorname{argmax}} \left( \log \left( p(C_k) \right) + \ \displaystyle\sum_{i=1}^n \log \left(p(x_i \vert C_k) \right) \right) \end{align}$$
If we want to compute the class probability $p(Y=C|\mathbf{x})$, we will need to compute the denominator:
$$ \begin{align} \log \left( p(Y=C|\mathbf{x}) \right) &= \log \left( \frac{p(\mathbf{x}|Y=C)p(Y=C)}{~\sum_{k=1}^{|C|}{}p(\mathbf{x}|Y=C_k)p(Y=C_k)} \right)\\ &= \log \left( \underbrace{p(\mathbf{x}|Y=C)p(Y=C)}_{\text{numerator}} \right) - \log \left( \underbrace{~\sum_{k=1}^{|C|}{}p(\mathbf{x}|Y=C_k)p(Y=C_k)}_{\text{denominator}} \right)\\ \end{align} $$
The element $\log \left( ~\sum_{k=1}^{|C|}{}p(\mathbf{x}|Y=C_k)p(Y=C_k) \right)\\ $ may underflow because $p(x_i \vert C_k)$ can be very small: it is the same issue as in the numerator, but this time we have a summation inside the logarithm, which prevents us from transforming the $p(x_i \vert C_k)$ (can be close to 0) into $\log \left(p(x_i \vert C_k) \right)$ (negative and not close to 0 anymore, since $0 \leq p(x_i \vert C_k) \leq 1$). To circumvent this issue, we can use the fact that $p(x_i \vert C_k) = \exp \left( {\log \left(p(x_i \vert C_k) \right)} \right)$ to obtain:
$$\log \left( ~\sum_{k=1}^{|C|}{}p(\mathbf{x}|Y=C_k)p(Y=C_k) \right) =\log \left( ~\sum_{k=1}^{|C|}{} \exp \left( \log \left( p(\mathbf{x}|Y=C_k)p(Y=C_k) \right) \right) \right)$$
At that point, a new issue arises: $\log \left( p(\mathbf{x}|Y=C_k)p(Y=C_k) \right)$ may be quite negative, which implies that $ \exp \left( \log \left( p(\mathbf{x}|Y=C_k)p(Y=C_k) \right) \right) $ may become very close to 0, i.e. underflow. This is where we use the log-sum-exp trick:
$$\log \sum_k e^{a_k} = \log \sum_k e^{a_k}e^{A-A} = A+ \log\sum_k e^{a_k -A}$$
with:
We can see that introducing the variable $A$ avoids underflows. E.g. with $k=2, a_1 = - 245, a_2 = - 255$, we have:
Using the log-sum-exp trick we avoid the underflow, with $A=\max ( -245, -255 )=-245$: $\begin{align}\log \sum_k e^{a_k} &= \log \sum_k e^{a_k}e^{A-A} \\&= A+ \log\sum_k e^{a_k -A}\\ &= -245+ \log\sum_k e^{a_k +245}\\&= -245+ \log \left(e^{-245 +245}+e^{-255 +245}\right) \\&=-245+ \log \left(e^{0}+e^{-10}\right) \end{align}$
We avoided the underflow since $e^{-10}$ is much farther away from 0 than $3.96143\times 10^{- 107}$ or $1.798486 \times 10^{-111}$.
Suppose we want to identify which of two databases is more likely to have generated a phrase (for example, which novel is this phrase more likely to have come from). We could assume independence of the words conditional on the database (Naive Bayes assumption).
Now look up the second link you have posted. There $a$ would represent the joint probability of observing the sentence given a database and the $e^{b_{t}}$s would represent the probability of observing each of the words in the sentence.
We can see from this answer that the smallest number in Python(just take it for example) is 5e-324 due to the IEEE754, and the hardware cause applies to other languages as well.
In [2]: np.nextafter(0, 1)
Out[2]: 5e-324
And any float smaller than that would lead to 0.
In [3]: np.nextafter(0, 1)/2
Out[3]: 0.0
And let's see the function of Naive Bayes with discrete features and two classes as you required:
$$ p(S=1|w_1, ... w_n) = \frac{p(S=1) \prod_{i=1}^n p(\mathbf{w_i}|S=1)}{~\sum_{s=\{0, 1\}}p(S=s)\prod_{i=1}^n p(\mathbf{w_i}|S=s)} $$
Let me instantiate that function by a simple NLP task bellow.
We decide to detect if the coming email is spam($S=1$) or not spam($S=0$) and we have a word vocabulary of 5,000 in size($n=5,000$) and the only concern is if a word($w_i$) occurs($p(w_i|S=1)$) in the email or not($1-p(w_i|S=1)$) for simplicity(Bernoulli naive Bayes).
In [1]: import numpy as np
In [2]: from sklearn.naive_bayes import BernoulliNB
# let's train our model with 200 samples
In [3]: X = np.random.randint(2, size=(200, 5000))
In [4]: y = np.random.randint(2, size=(200, 1)).ravel()
In [5]: clf = BernoulliNB()
In [6]: model = clf.fit(X, y)
We can see that $p(S=s)\prod_{i=1}^n p(\mathbf{w_i}|S=s)$ would be very small because of the probabilities(both $p(w_i|S=1)$ and $1-p(w_i|S=1)$ would be between 0 and 1) in $\prod_i^{5000}$, and hence we are sure that the product would be smaller than $5e^{-324}$ and we just get $0/0$.
In [7]: (np.nextafter(0, 1)*2) / (np.nextafter(0, 1)*2)
Out[7]: 1.0
In [8]: (np.nextafter(0, 1)/2) / (np.nextafter(0, 1)/2)
/home/lerner/anaconda3/bin/ipython3:1: RuntimeWarning: invalid value encountered in double_scalars
#!/home/lerner/anaconda3/bin/python
Out[8]: nan
In [9]: l_cpt = model.feature_log_prob_
In [10]: x = np.random.randint(2, size=(1, 5000))
In [11]: cls_lp = model.class_log_prior_
In [12]: probs = np.where(x, np.exp(l_cpt[1]), 1-np.exp(l_cpt[1]))
In [13]: np.exp(cls_lp[1]) * np.prod(probs)
Out[14]: 0.0
Then the problem arises: How can we calculate the probability of the email is a spam $p(S=1|w_1, ... w_n)$? Or how can we calculate the numerator and the denominator?
We can see the official implementation in sklearn:
jll = self._joint_log_likelihood(X)
# normalize by P(x) = P(f_1, ..., f_n)
log_prob_x = logsumexp(jll, axis=1)
return jll - np.atleast_2d(log_prob_x).T
For the numerator it converted the product of probabilities to the sum of log likelihood and for the denominator it used the logsumexp in scipy which is:
out = log(sum(exp(a - a_max), axis=0))
out += a_max
Because we cannot add two joint probabilities by adding its joint log likelihood, and we should go out from the log space to the probability space. But we cannot add the two true probabilities because they are too small and we should scale them and do the addition: $\sum_{s=\{0,1\}} e^{jll_s - max\_jll}$ and put the result back into the log space $\log\sum_{s=\{0,1\}} e^{jll_s - max\_jll}$ then rescale it back: $max\_jll+ \log\sum_{s=\{0,1\}} e^{jll_s - max\_jll}$ in log space by adding the $max\_jll$.
And here is the derivation:
$\begin{align} \log \sum_{s=\{0,1\}} e^{jll_s} & = \log \sum_{s=\{0,1\}} e^{jll_s}e^{max\_jll-max\_jll} \\& = \log e ^{max\_jll}+ \log\sum_{s=\{0,1\}} e^{jll_s - max\_jll} \\& = max\_jll+ \log\sum_{s=\{0,1\}} e^{jll_s - max\_jll} \end{align}$
where $max\_jll$ is the $a\_max$ in the code.
Once we get both the numerator and the denominator in log space we can get the log conditional probability($\log p(S=1|w_1, ... w_n)$) by subtracting the denominator from the numerator:
return jll - np.atleast_2d(log_prob_x).T
Hope that helps.
Reference:
1. Bernoulli Naive Bayes Classifier
2. Spam Filtering with Naive Bayes – Which Naive Bayes?
