# Implementation of Kullback–Leibler divergence in Mathematica: dealing with negative infinity?

I’m trying to implement a function to calculate Kullback-Liebler Divergence but I am running into a problem with complex infinities in my intermediate result. Here is how I have implemented it thus far:

pmfA={1/6,1/6,1/6,1/6,1/6,1/6};
pmfB={1/4,1/4,1/8,1/8,1/8,1/8};

Function[{p, q}, p*Log[p] /Log[q]][pmfA,pmfB]
Total[%]

For some “nice” sample data (containing no zeros), we get the nice result:

1.00526

However, if we use some “not so nice” sample data (containing zeros for some components of the PMF), it won’t work, because of a negative infinity for Log[0].

eg. Try it out with…

pmfA={1/6,1/6,1/6,1/6,1/6,1/6};
pmfB={1/4,1/4,1/2,0,0,0};

For the purpose of computing the total, I think it would be acceptable to treat the instances as “0” but I am unsure of how to implement that.

=================

1

p Log[p]/Log[q] /. -\[Infinity] -> 0 or DeleteCases[p Log[p]/Log[q], -\[Infinity]] are two possible implementations.
– C. E.
Jan 20 at 21:07

1

Your definition is incorrect. The ratio p/q should be inside the logarithm. @Pickett You would also have to guard against Indeterminate when machine-precision numbers are used, I think. But anyway, I think the point is moot – there’s only one logarithm.
– Jens
Jan 20 at 23:30

=================

1

=================

My interpretation of the question is that, after fixing the error in the definition of the divergence, the mathematical problem that is being asked about arises when the numerator under the Log goes to zero – not (as is written in the question) the denominator. Here is how you can then get the correct result in this case:

pmfA = {1/6, 1/6, 1/6, 1/6, 1/6, 1/6};
pmfB = {1/4, 1/4, 1/2, 0, 0, 0};

Function[{p, q},
Limit[p*Log[(p + ϵ)/(q + ϵ)], ϵ ->
0]][pmfB, pmfA]
Total[%]

(* ==> {1/4 Log[3/2], 1/4 Log[3/2], Log[3]/2, 0, 0, 0} *)

(* ==> 1/2 Log[3/2] + Log[3]/2 *)

So I switched the roles of pmfB, pmfA to get to the mathematically tractable case where the logarithm goes to negative infinity while being multiplied by 0. This is calculated by introducing a formal Limit.