Adds weighted permutation entropy (#19)

pyentrp/entropy.py

@@ -288,6 +288,76 @@ def multiscale_permutation_entropy(time_series, m, delay, scale):
     return mspe
 
 
+def weighted_permutation_entropy(time_series, order=2, delay=1, normalize=False):
+        """Calculate the weighted permuation entropy. Weighted permutation
+        entropy captures the information in the amplitude of a signal where
+        standard permutation entropy only measures the information in the
+        ordinal pattern, "motif."
+
+        Parameters
+        ----------
+        time_series : list or np.array
+            Time series
+        order : int
+            Order of permutation entropy
+        delay : int
+            Time delay
+        normalize : bool
+            If True, divide by log2(factorial(m)) to normalize the entropy
+            between 0 and 1. Otherwise, return the permutation entropy in bit.
+
+        Returns
+        -------
+        wpe : float
+            Weighted Permutation Entropy
+
+        References
+        ----------
+        .. [1] Bilal Fadlallah, Badong Chen, Andreas Keil, and José Príncipe
+               Phys. Rev. E 87, 022911 – Published 20 February 2013
+
+        Notes
+        -----
+        Last updated (March 2021) by Samuel Dotson (samgdotson@gmail.com)
+
+        Examples
+        --------
+        1. Weighted permutation entropy with order 2
+
+            >>> x = [4, 7, 9, 10, 6, 11, 3]
+            >>> # Return a value between 0 and log2(factorial(order))
+            >>> print(permutation_entropy(x, order=2))
+                0.912
+
+        2. Normalized weighted permutation entropy with order 3
+
+            >>> x = [4, 7, 9, 10, 6, 11, 3]
+            >>> # Return a value comprised between 0 and 1.
+            >>> print(permutation_entropy(x, order=3, normalize=True))
+                0.547
+        """
+        x = _embed(time_series, order=order, delay=delay)
+
+        weights = np.var(x, axis=1)
+        sorted_idx = x.argsort(kind='quicksort', axis=1)
+        motifs, c = np.unique(sorted_idx, return_counts=True, axis=0)
+        pw = np.zeros(len(motifs))
+
+        # TODO hashmap
+        for i, j in zip(weights, sorted_idx):
+            idx = int(np.where((j==motifs).sum(1)==order)[0])
+            pw[idx] += i
+
+        pw /= weights.sum()
+
+        b = np.log2(pw)
+        wpe = -np.dot(pw, b)
+        if normalize:
+            wpe /= np.log2(factorial(order))
+        return wpe
+
+
+
 # TODO add tests
 def composite_multiscale_entropy(time_series, sample_length, scale, tolerance=None):
     """Calculate the Composite Multiscale Entropy of the given time series.

tests/test_entropy.py

@@ -44,6 +44,12 @@ def test_permutationEntropy(self):
         # Assert that a fully random vector has an entropy of 0.99999...
         self.assertEqual(np.round(ent.permutation_entropy(RANDOM_TIME_SERIES, order=3, delay=1, normalize=True), 3), 0.999)
 
+    def test_weightedPermuationEntropy(self):
+        self.assertEqual(np.round(ent.weighted_permutation_entropy(PERM_ENTROPY_BANDT, order=2, delay=1), 3), 0.913)
+        self.assertEqual(np.round(ent.weighted_permutation_entropy(PERM_ENTROPY_BANDT, order=3, delay=1), 3), 1.414)
+        # Assert that a fully random vector has an entropy of 0.99999...
+        self.assertEqual(np.round(ent.weighted_permutation_entropy(RANDOM_TIME_SERIES, order=3, delay=1, normalize=True), 3), 0.999)
+
     def test_multiScalePermutationEntropy(self):
         np.testing.assert_array_equal(np.round(ent.multiscale_permutation_entropy(TS_SAMPLE_ENTROPY, 3, 5, 2), 4),
                                       np.array([2.4699, 2.5649]))