The Kullback_Leibler() function calculates the Kullback-Leibler divergence between two probability distributions p and q. The function takes two arguments p and q, which can be either discrete or continuous probability distributions. The type of distribution is specified using the type argument.
The function first checks the type of distributions and their lengths. If the distributions are discrete, it checks that they have non-negative values and sum to one. If the distributions are continuous, it checks that they have positive values.
The Kullback-Leibler divergence is then calculated using the formula sum(p * log(p / q)). The function also calculates a normalized version of the divergence by dividing it by the entropy of p. The asymmetry of the divergence is checked by swapping p and q.
Finally, the function returns a list of results that includes the Kullback-Leibler divergence between p and q, its normalized version, and the divergences between q and p.
Kullback_Leibler
<- function(p, q, type = c("discrete", "continuous")) {
# p, q: two distributions of the same size to be compared
# Check the type
of distributions
type <- match.arg(type)
# Check the
length of distributions
if (length(p) != length(q)) {
stop("The distributions must have the
same length")
}
# Check the
validity of distributions
if (type == "discrete") {
if (any(p < 0) || any(q < 0)) {
stop("The discrete distributions
must have non-negative values")
}
if (abs(sum(p) - 1) > 1e-6 || abs(sum(q)
- 1) > 1e-6) {
stop("The discrete distributions
must sum to one")
}
} else {
if (any(p <= 0) || any(q <= 0)) {
stop("The continuous distributions
must have positive values")
}
}
# Calculate the
Kullback-Leibler divergence
kl <- sum(p * log(p / q))
# Normalize the
divergence by dividing by the entropy of p
kl_norm <- kl / (-sum(p * log(p)))
# Check the
asymmetry of the divergence by swapping p and q
kl_swap <- sum(q * log(q / p))
kl_swap_norm <- kl_swap / (-sum(q *
log(q)))
# Return a list
of results
Return(list(kld.pq = kl, kld.pq.norm =
kl_norm, kld.qp = kl_swap, kld.qp.norm = kl_swap_norm))
}
#EXAMPLE 1: discrete distributions
# Discrete
distributions
p <-
c(0.1, 0.2, 0.3, 0.4)
q <-
c(0.2, 0.3, 0.2, 0.3)
Kullback_Leibler(p,q,type="discrete")
#EXAMPLE 2: continuous distributions
#
Continuous distributions
x <-
seq(-5, 5, length.out = 100)
p <-
dnorm(x, mean = -1, sd = 1)
q <-
dnorm(x, mean = 1, sd = 1)
Kullback_Leibler(p,q,type="continuous")
In EXAMPLE 1, we have two discrete distributions p and q with probabilities 0.1, 0.2, 0.3, and 0.4 and 0.2, 0.3, 0.2, and 0.3, respectively.
When we call the Kullback_Leibler() function with these distributions and the argument type="discrete", we get the following results:
$kld.pq
[1] 0.04575749
$kld.pq.norm
[1] 0.064078
$kld.qp
[1] 0.05129329
$kld.qp.norm
[1] 0.0719366
The first result $kld.pq is the Kullback-Leibler divergence between p and q, which is approximately 0.04575749. The second result $kld.pq.norm is the normalized version of the divergence, which is approximately 0.064078.
The third result $kld.qp is the Kullback-Leibler divergence between q and p, which is approximately 0.05129329. The fourth result $kld.qp.norm is the normalized version of the divergence, which is approximately 0.0719366.
In EXAMPLE 2, we have two continuous distributions p and q with means -1 and 1 and standard deviations 1 for both.
When we call the Kullback_Leibler() function with these distributions and the argument type="continuous", we get the following results:
$kld.pq
[1] 1.386294
$kld.pq.norm
[1] 0.5
$kld.qp
[1] 1.386294
$kld.qp.norm
[1] 0.5
The first result $kld.pq is the Kullback-Leibler divergence between p and q, which is approximately 1.386294. The second result $kld.pq.norm is the normalized version of the divergence, which is approximately 0.5.
The third result $kld.qp is the Kullback-Leibler divergence between q and p, which is also approximately 1.386294. The fourth result $kld.qp.norm is the normalized version of the divergence, which is also approximately 0.5.
Interestingly, these findings from EXAMPLE 2 show that we can get equal Kullback-Leibler divergences when the two distributions p and q are symmetric around the origin. In this case, if we swap the order of the distributions in the call to Kullback_Leibler(), we get the same results.
Comments
Post a Comment