information/
joint.rs

1/// # Joint Entropy
2/// Generalization of Entropy to multiple dimensions
3///
4/// <https://en.wikipedia.org/wiki/Joint_entropy>
5///
6/// Joint entropy is calculated for multiple variables as:
7/// ```math
8/// H(Xi ... Xn) = - Σi ... Σn P( xi, ..., xn ) * ln[ P(xi, ..., xn) ]
9/// ```
10///
11/// # Usage
12/// ```
13/// use ndarray::{Array1, Array2, Array3, Array4};
14/// use ndarray_rand::{RandomExt, rand_distr::Uniform};
15/// use information::joint_entropy;
16///
17/// // 1D Entropy
18/// let c_x = Array1::random(10, Uniform::new(0.1, 0.8));
19/// let p_x = &c_x / c_x.sum();
20/// let h = joint_entropy!(&p_x);
21/// assert!(h >= 0.0);
22///
23/// // 2D Entropy
24/// let c_xy = Array2::random((2, 10), Uniform::new(0.1, 0.8));
25/// let p_xy = &c_xy / c_xy.sum();
26/// let h = joint_entropy!(&p_xy);
27/// assert!(h >= 0.0);
28///
29/// // 3D Entropy
30/// let c_xy = Array3::random((2, 2, 10), Uniform::new(0.1, 0.8));
31/// let p_xy = &c_xy / c_xy.sum();
32/// let h = joint_entropy!(&p_xy);
33/// assert!(h >= 0.0);
34///
35/// // 4D Entropy
36/// let c_xy = Array4::random((2, 2, 2, 10), Uniform::new(0.1, 0.8));
37/// let p_xy = &c_xy / c_xy.sum();
38/// let h = joint_entropy!(&p_xy);
39/// assert!(h >= 0.0);
40/// ```
41///
42#[macro_export]
43macro_rules! joint_entropy {
44    ($prob:expr) => {
45        $prob.iter().fold(0.0, |acc, p| {
46            if *p == 0.0 {
47                acc
48            } else {
49                acc - (p * (*p as f64).ln())
50            }
51        })
52    };
53}
54
55#[cfg(test)]
56mod testing {
57
58    use crate::{
59        entropy::entropy,
60        joint_entropy,
61        prob::{prob1d, prob2d},
62    };
63    use approx::assert_relative_eq;
64    use ndarray::{array, Array1, Array2, Array3, Array4};
65    use ndarray_rand::{rand_distr::Uniform, RandomExt};
66
67    const N_ITER: usize = 1000;
68    const ARRAY_SIZE: usize = 100;
69
70    #[test]
71    fn test_joint_macro() {
72        for _ in 0..N_ITER {
73            // 1D Entropy
74            let c_x = Array1::random(ARRAY_SIZE, Uniform::new(0.1, 0.8));
75            let p_x = &c_x / c_x.sum();
76            let h = joint_entropy!(&p_x);
77            assert!(h >= 0.0);
78            assert_relative_eq!(h, entropy(&p_x));
79
80            // 2D Entropy
81            let c_xy = Array2::random((2, ARRAY_SIZE), Uniform::new(0.1, 0.8));
82            let p_xy = &c_xy / c_xy.sum();
83            let h = joint_entropy!(&p_xy);
84            assert!(h >= 0.0);
85
86            // 3D Entropy
87            let c_xy = Array3::random((2, 2, ARRAY_SIZE), Uniform::new(0.1, 0.8));
88            let p_xy = &c_xy / c_xy.sum();
89            let h = joint_entropy!(&p_xy);
90            assert!(h >= 0.0);
91
92            // 4D Entropy
93            let c_xy = Array4::random((2, 2, 2, ARRAY_SIZE), Uniform::new(0.1, 0.8));
94            let p_xy = &c_xy / c_xy.sum();
95            let h = joint_entropy!(&p_xy);
96            assert!(h >= 0.0);
97        }
98    }
99
100    #[test]
101    fn test_joint() {
102        let px = array![[0.5, 0.0], [0.25, 0.25]];
103        let hx = joint_entropy!(&px);
104        assert_eq!(hx, 1.0397207708399179);
105    }
106
107    #[test]
108    /// https://en.wikipedia.org/wiki/Joint_entropy#Nonnegativity
109    fn test_nonnegative() {
110        for _ in 0..N_ITER {
111            let c_xy = Array2::random((2, ARRAY_SIZE), Uniform::new(0.1, 0.8));
112            let p_xy = &c_xy / c_xy.sum();
113            let h = joint_entropy!(&p_xy);
114            assert!(h >= 0.0);
115        }
116    }
117
118    #[test]
119    /// https://en.wikipedia.org/wiki/Joint_entropy#Nonnegativity
120    fn test_nonnegative_3d() {
121        for _ in 0..N_ITER {
122            let c_xyz = Array3::random((2, 2, ARRAY_SIZE), Uniform::new(0.1, 0.8));
123            let p_xyz = &c_xyz / c_xyz.sum();
124            let h = joint_entropy!(&p_xyz);
125            assert!(h >= 0.0);
126        }
127    }
128
129    #[test]
130    /// https://en.wikipedia.org/wiki/Joint_entropy#Greater_than_individual_entropies
131    fn test_gte_individual_entropies() {
132        for _ in 0..N_ITER {
133            let c_x = Array1::random(ARRAY_SIZE, Uniform::new(0, 3));
134            let c_y = Array1::random(ARRAY_SIZE, Uniform::new(0, 3));
135
136            let p_xy = prob2d(&c_x, &c_y, 4, 4).unwrap();
137            let p_x = prob1d(&c_x, 4).unwrap();
138            let p_y = prob1d(&c_y, 4).unwrap();
139
140            let h_xy = joint_entropy!(&p_xy);
141            let h_x = entropy(&p_x);
142            let h_y = entropy(&p_y);
143
144            assert!(h_xy >= h_x);
145            assert!(h_xy >= h_y);
146        }
147    }
148
149    #[test]
150    /// https://en.wikipedia.org/wiki/Joint_entropy#Less_than_or_equal_to_the_sum_of_individual_entropies
151    fn test_lte_individual_entropies() {
152        for _ in 0..N_ITER {
153            let c_x = Array1::random(ARRAY_SIZE, Uniform::new(0, 3));
154            let c_y = Array1::random(ARRAY_SIZE, Uniform::new(0, 3));
155
156            let p_xy = prob2d(&c_x, &c_y, 4, 4).unwrap();
157            let p_x = prob1d(&c_x, 4).unwrap();
158            let p_y = prob1d(&c_y, 4).unwrap();
159
160            let h_xy = joint_entropy!(&p_xy);
161            let h_x = entropy(&p_x);
162            let h_y = entropy(&p_y);
163
164            assert!(h_xy <= h_x + h_y);
165        }
166    }
167}