Initial commit with perm_entropy

entropy/__init__.py

@@ -0,0 +1,5 @@
+# Import EntroPy objects
+from .utils import *
+from .entropy import *
+
+__version__ = "0.1.0"

entropy/entropy.py

@@ -0,0 +1,65 @@
+import numpy as np
+from math import factorial
+
+from .utils import _embed
+
+all = ['perm_entropy']
+
+
+def perm_entropy(x, order=3, delay=1, normalize=False):
+    """Permutation Entropy (Bandt and Pompe 2002)
+
+    Parameters
+    ----------
+    x : list or np.array
+        Time series
+    order : int
+        Order of permutation entropy
+    delay : int
+        Time delay
+    normalize : bool
+        If True, divide by log2(factorial(order)) to normalize the entropy
+        between 0 and 1. Otherwise, return the permutation entropy in bit.
+
+    Returns
+    -------
+    pe : float
+        Permutation Entropy
+
+    References
+    ----------
+    .. [1] Massimiliano Zanin et al. Permutation Entropy and Its Main
+        Biomedical and Econophysics Applications: A Review.
+        http://www.mdpi.com/1099-4300/14/8/1553/pdf
+    .. [2] Christoph Bandt and Bernd Pompe. Permutation entropy — a natural
+        complexity measure for time series.
+        http://stubber.math-inf.uni-greifswald.de/pub/full/prep/2001/11.pdf
+
+    Examples
+    --------
+    1. Permutation entropy with order 2
+        >>> x = [4, 7, 9, 10, 6, 11, 3]
+        >>> # Return a value between 0 and log2(factorial(order))
+        >>> print(perm_entropy(x, order=2))
+            0.918
+    2. Normalized permutation entropy with order 3
+        >>> x = [4, 7, 9, 10, 6, 11, 3]
+        >>> # Return a value comprised between 0 and 1.
+        >>> print(perm_entropy(x, order=3, normalize=True))
+            0.589
+    """
+    x = np.array(x)
+    ran_order = range(order)
+    hashmult = np.power(order, ran_order)
+    # Embed x and sort the order of permutations
+    sorted_idx = _embed(x, order=order, delay=delay).argsort(kind='quicksort')
+    # Associate unique integer to each permutations
+    hashval = (np.multiply(sorted_idx, hashmult)).sum(1)
+    # Return the counts
+    _, c = np.unique(hashval, return_counts=True)
+    # Use np.true_divide for Python 2 compatibility
+    p = np.true_divide(c, c.sum())
+    pe = -np.multiply(p, np.log2(p)).sum()
+    if normalize:
+        pe /= np.log2(factorial(order))
+    return pe

entropy/tests/test_entropy.py

@@ -0,0 +1,23 @@
+import numpy as np
+import unittest
+
+from entropy import perm_entropy
+
+np.random.seed(1234567)
+RANDOM_TS = np.random.rand(1000)
+BANDT_PERM = [4, 7, 9, 10, 6, 11, 3]
+
+
+class TestEntropy(unittest.TestCase):
+
+    def test_permutation_entropy(self):
+        self.assertEqual(np.round(perm_entropy(RANDOM_TS, order=3,
+                                               delay=1, normalize=True), 3),
+                         0.999)
+        # Compare with Bandt and Pompe 2002
+        self.assertEqual(np.round(perm_entropy(BANDT_PERM, order=2), 3), 0.918)
+        self.assertEqual(np.round(perm_entropy(BANDT_PERM, order=3), 3), 1.522)
+
+
+if __name__ == '__main__':
+    unittest.main()

entropy/utils.py

@@ -0,0 +1,27 @@
+import numpy as np
+
+all = ['_embed']
+
+
+def _embed(x, order=3, delay=1):
+    """Time-delay embedding.
+
+    Parameters
+    ----------
+    x : 1d-array, shape (n_times)
+        Time series
+    order : int
+        Embedding dimension (order)
+    delay : int
+        Delay.
+
+    Returns
+    -------
+    embedded : ndarray, shape (n_times - (order - 1) * delay, order)
+        Embedded time-series.
+    """
+    N = len(x)
+    Y = np.zeros((order, N - (order - 1) * delay))
+    for i in range(order):
+        Y[i] = x[i * delay:i * delay + Y.shape[1]]
+    return Y.T

notebooks/entropy.ipynb

@@ -0,0 +1,234 @@
+{
+ "cells": [
+  {
+   "cell_type": "code",
+   "execution_count": 2,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "import itertools\n",
+    "import numpy as np\n",
+    "from math import factorial\n",
+    "\n",
+    "def _embed(x, order=3, delay=1):\n",
+    "    N = len(x)\n",
+    "    Y = np.zeros((order, N - (order - 1) * delay))\n",
+    "    for i in range(order):\n",
+    "        Y[i] = x[i * delay : i * delay + Y.shape[1]]\n",
+    "    return Y.T\n",
+    "\n",
+    "def permutation_entropy(x, order=3, delay=1, normalize=False):\n",
+    "    \"\"\"Permutation Entropy (Bandt and Pompe 2002)\n",
+    "\n",
+    "    Parameters\n",
+    "    ----------\n",
+    "    x : list or np.array\n",
+    "        Time series\n",
+    "    order : int\n",
+    "        Order of permutation entropy\n",
+    "    delay : int\n",
+    "        Time delay\n",
+    "    normalize : bool\n",
+    "        If True, divide by log2(factorial(order)) to normalize the entropy\n",
+    "        between 0 and 1. Otherwise, return the permutation entropy in bit.\n",
+    "\n",
+    "    Returns\n",
+    "    -------\n",
+    "    pe : float\n",
+    "        Permutation Entropy\n",
+    "\n",
+    "    References\n",
+    "    ----------\n",
+    "    .. [1] Massimiliano Zanin et al. Permutation Entropy and Its Main\n",
+    "        Biomedical and Econophysics Applications: A Review.\n",
+    "        http://www.mdpi.com/1099-4300/14/8/1553/pdf\n",
+    "    .. [2] Christoph Bandt and Bernd Pompe. Permutation entropy — a natural\n",
+    "        complexity measure for time series.\n",
+    "        http://stubber.math-inf.uni-greifswald.de/pub/full/prep/2001/11.pdf\n",
+    "\n",
+    "    Examples\n",
+    "    --------\n",
+    "    1. Permutation entropy with order 2\n",
+    "        >>> x = [4, 7, 9, 10, 6, 11, 3]\n",
+    "        >>> # Return a value between 0 and log2(factorial(order))\n",
+    "        >>> print(permutation_entropy(x, order=2))\n",
+    "            0.918\n",
+    "    2. Normalized permutation entropy with order 3\n",
+    "        >>> x = [4, 7, 9, 10, 6, 11, 3]\n",
+    "        >>> # Return a value comprised between 0 and 1.\n",
+    "        >>> print(permutation_entropy(x, order=3, normalize=True))\n",
+    "            0.589\n",
+    "    \"\"\"\n",
+    "    x = np.array(x)\n",
+    "    ran_order = range(order)\n",
+    "    hashmult = np.power(order, ran_order)\n",
+    "    # Embed x and sort the order of permutations\n",
+    "    sorted_idx = _embed(x, order=order, delay=delay).argsort(kind='quicksort')\n",
+    "    # Associate unique integer to each permutations\n",
+    "    hashval = (np.multiply(sorted_idx, hashmult)).sum(1)\n",
+    "    # Return the counts\n",
+    "    _, c = np.unique(hashval, return_counts=True)\n",
+    "    # Use np.true_divide for Python 2 compatibility\n",
+    "    p = np.true_divide(c, c.sum())\n",
+    "    pe = -np.multiply(p, np.log2(p)).sum()\n",
+    "    if normalize:\n",
+    "        pe /= np.log2(factorial(order))\n",
+    "    return pe"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 3,
+   "metadata": {},
+   "outputs": [
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "0.9182958340544896\n",
+      "0.5887621559162939\n"
+     ]
+    }
+   ],
+   "source": [
+    "x = [4, 7, 9, 10, 6, 11, 3]\n",
+    "print(permutation_entropy(x, order=2))\n",
+    "print(permutation_entropy(x, order=3, normalize=True))"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 4,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "def hash_term(perm):\n",
+    "    \"\"\"Associate unique integer to a permutation.\n",
+    "    \"\"\"\n",
+    "    deg = len(perm)\n",
+    "    return sum([perm[k]*deg**k for k in range(deg)])\n",
+    "\n",
+    "\n",
+    "def permutation_entropy_orig(x, order=3, delay=1, normalize=False):\n",
+    "    \"\"\"Permutation Entropy (Bandt and Pompe 2002)\n",
+    "\n",
+    "    Parameters\n",
+    "    ----------\n",
+    "    x : list or np.array\n",
+    "        Time series\n",
+    "    order : int\n",
+    "        Order of permutation entropy\n",
+    "    delay : int\n",
+    "        Time delay\n",
+    "    normalize : bool\n",
+    "        If True, divide by log2(factorial(order)) to normalize the entropy\n",
+    "        between 0 and 1. Otherwise, return the permutation entropy in bit.\n",
+    "\n",
+    "    Returns\n",
+    "    -------\n",
+    "    pe : float\n",
+    "        Permutation Entropy\n",
+    "\n",
+    "    References\n",
+    "    ----------\n",
+    "    .. [1] Massimiliano Zanin et al. Permutation Entropy and Its Main\n",
+    "        Biomedical and Econophysics Applications: A Review.\n",
+    "        http://www.mdpi.com/1099-4300/14/8/1553/pdf\n",
+    "\n",
+    "    .. [2] Christoph Bandt and Bernd Pompe. Permutation entropy — a natural\n",
+    "        complexity measure for time series.\n",
+    "        http://stubber.math-inf.uni-greifswald.de/pub/full/prep/2001/11.pdf\n",
+    "\n",
+    "    Examples\n",
+    "    --------\n",
+    "    1. Permutation entropy with order 2\n",
+    "\n",
+    "        >>> x = [4, 7, 9, 10, 6, 11, 3]\n",
+    "        >>> # Return a value between 0 and log2(factorial(order))\n",
+    "        >>> print(permutation_entropy(x, order=2))\n",
+    "            0.918\n",
+    "\n",
+    "    2. Normalized permutation entropy with order 3\n",
+    "\n",
+    "        >>> x = [4, 7, 9, 10, 6, 11, 3]\n",
+    "        >>> # Return a value comprised between 0 and 1.\n",
+    "        >>> print(permutation_entropy(x, order=3, normalize=True))\n",
+    "            0.589\n",
+    "    \"\"\"\n",
+    "    x = np.array(x)\n",
+    "    permutations = list(itertools.permutations(range(order)))\n",
+    "    hashlist = [hash_term(perm) for perm in permutations]\n",
+    "    c = np.zeros(len(permutations), dtype=int)\n",
+    "    ran = np.arange(len(x) - delay * (order - 1))\n",
+    "    step = ran + order * delay\n",
+    "\n",
+    "    for i in ran:\n",
+    "        sorted_idx = x[i:step[i]:delay].argsort(kind='quicksort')\n",
+    "        c[np.nonzero(hashlist == hash_term(sorted_idx))[0][0]] += 1\n",
+    "\n",
+    "    c = c[np.nonzero(c)]\n",
+    "    # Use np.true_divide for Python 2 compatibility\n",
+    "    p = np.true_divide(c, c.sum())\n",
+    "    pe = -np.multiply(p, np.log2(p)).sum()\n",
+    "    if normalize:\n",
+    "        pe /= np.log2(factorial(order))\n",
+    "    return pe"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 7,
+   "metadata": {},
+   "outputs": [
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "0.9999781776726874\n",
+      "0.9999781776726874\n",
+      "2.76 ms ± 135 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)\n",
+      "288 ms ± 20.2 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)\n"
+     ]
+    }
+   ],
+   "source": [
+    "x = np.random.rand(30000)\n",
+    "# x = [4, 7, 9, 10, 6, 11, 3]\n",
+    "\n",
+    "print(permutation_entropy(x, order=2, normalize=True))\n",
+    "print(permutation_entropy_orig(x, order=2, normalize=True))\n",
+    "\n",
+    "%timeit permutation_entropy(x, order=3)\n",
+    "%timeit permutation_entropy_orig(x, order=3)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": []
+  }
+ ],
+ "metadata": {
+  "kernelspec": {
+   "display_name": "Python 3",
+   "language": "python",
+   "name": "python3"
+  },
+  "language_info": {
+   "codemirror_mode": {
+    "name": "ipython",
+    "version": 3
+   },
+   "file_extension": ".py",
+   "mimetype": "text/x-python",
+   "name": "python",
+   "nbconvert_exporter": "python",
+   "pygments_lexer": "ipython3",
+   "version": "3.6.5"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 2
+}

requirements.txt

@@ -0,0 +1,2 @@
+numpy
+scipy

setup.cfg

@@ -0,0 +1,19 @@
+[aliases]
+test=pytest
+
+# [tool:pytest]
+# addopts = --showlocals --durations=10 --cov --cov-report=
+# markers =
+#     slow: mark a test as slow.
+
+[flake8]
+ignore = N806, N803
+exclude =
+    .git,
+    __pycache__,
+    docs,
+    external,
+    __init__.py,
+    examples,
+    setup.py,
+statistics=True