Bootstraping repo, and setting up the RPDE notebook.
This commit is contained in:
Executable
BIN
Binary file not shown.
@@ -0,0 +1,464 @@
|
||||
{
|
||||
"cells": [
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 6,
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"data": {
|
||||
"text/html": [
|
||||
" <script type=\"text/javascript\">\n",
|
||||
" window.PlotlyConfig = {MathJaxConfig: 'local'};\n",
|
||||
" if (window.MathJax) {MathJax.Hub.Config({SVG: {font: \"STIX-Web\"}});}\n",
|
||||
" if (typeof require !== 'undefined') {\n",
|
||||
" require.undef(\"plotly\");\n",
|
||||
" requirejs.config({\n",
|
||||
" paths: {\n",
|
||||
" 'plotly': ['https://cdn.plot.ly/plotly-latest.min']\n",
|
||||
" }\n",
|
||||
" });\n",
|
||||
" require(['plotly'], function(Plotly) {\n",
|
||||
" window._Plotly = Plotly;\n",
|
||||
" });\n",
|
||||
" }\n",
|
||||
" </script>\n",
|
||||
" "
|
||||
]
|
||||
},
|
||||
"metadata": {},
|
||||
"output_type": "display_data"
|
||||
}
|
||||
],
|
||||
"source": [
|
||||
"import librosa\n",
|
||||
"import numpy as np\n",
|
||||
"import plotly.express as px\n",
|
||||
"import plotly.graph_objs as go\n",
|
||||
"from plotly.offline import init_notebook_mode, iplot\n",
|
||||
"init_notebook_mode(connected=True)\n",
|
||||
"from numpy.linalg import norm\n",
|
||||
"\n",
|
||||
"from tqdm import tqdm\n"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"source": [
|
||||
"### First things first\n",
|
||||
"\n",
|
||||
"_Loading the audio data, and embedding the audio time series_\n",
|
||||
" Our algorithm operates on time series\n",
|
||||
"(here, it's audio, what a coincidence!).\n"
|
||||
],
|
||||
"metadata": {
|
||||
"collapsed": false,
|
||||
"pycharm": {
|
||||
"name": "#%% md\n"
|
||||
}
|
||||
}
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"# loads the sound data into a float32\n",
|
||||
"data, rate = librosa.core.load(\"pho_a.wav\")\n",
|
||||
"\n",
|
||||
"px.line(data, x=range(len(data)), y=data)\n"
|
||||
],
|
||||
"metadata": {
|
||||
"collapsed": false,
|
||||
"pycharm": {
|
||||
"name": "#%%\n"
|
||||
}
|
||||
}
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 7,
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"def embed_time_series(data: np.ndarray, dim: int, tau: int):\n",
|
||||
" \"\"\"This creates an special embedding for the input data. The\n",
|
||||
" output shape of this function is (N, D) with\n",
|
||||
" N = len(data) - (dim - 1) * tau\n",
|
||||
" D = dim\"\"\"\n",
|
||||
" embed_points_nb = data.shape[0] - (dim - 1) * tau\n",
|
||||
" embed_mask = np.arange(embed_points_nb)[:, np.newaxis] + np.arange(dim)\n",
|
||||
" tau_offsets = np.arange(dim) * (tau - 1)\n",
|
||||
" embed_mask = embed_mask + tau_offsets\n",
|
||||
" return data[embed_mask]"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 51,
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"# this creates an embedding of a slice of the time series\n",
|
||||
"ts = embed_time_series(data[:2000], dim=4, tau=22)\n",
|
||||
"print(ts.shape)"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"trace = go.Scatter3d(\n",
|
||||
" x=ts[:,0],\n",
|
||||
" y=ts[:,1],\n",
|
||||
" z=ts[:,2],\n",
|
||||
" name=\"Pho_a\",\n",
|
||||
" mode='markers',\n",
|
||||
" marker=dict(\n",
|
||||
" color=\"blue\",\n",
|
||||
" size=2,\n",
|
||||
" symbol='circle',\n",
|
||||
" line=dict(\n",
|
||||
" color='rgb(204, 204, 204)',\n",
|
||||
" width=1\n",
|
||||
" ),\n",
|
||||
" opacity=0.7\n",
|
||||
" )\n",
|
||||
" )\n",
|
||||
"layout = go.Layout(margin=dict(l=0,r=0,b=0,t=0),\n",
|
||||
" scene = dict(\n",
|
||||
" xaxis={\"title\":\"f1\"},\n",
|
||||
" yaxis={\"title\":\"f2\"},\n",
|
||||
" zaxis={\"title\":\"f3\"})\n",
|
||||
" )\n",
|
||||
"fig = go.Figure(data=[trace], layout=layout)\n",
|
||||
"iplot(fig, filename=\"/tmp/pho_a_trajectory.png\")"
|
||||
],
|
||||
"metadata": {
|
||||
"collapsed": false,
|
||||
"pycharm": {
|
||||
"name": "#%%\n"
|
||||
}
|
||||
}
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"source": [
|
||||
"### This is where the fun begins\n",
|
||||
"\n",
|
||||
"Now, let's implement the RPDE algorithm.\n",
|
||||
"I've provided 4 different implementations:\n",
|
||||
"\n",
|
||||
"1. The first one is a plain python implem of the recurrence diagram\n",
|
||||
"2. The second one is numba implem without any parallelization\n",
|
||||
"2. The third one is a naïve parallelized numba implem. It's \"almost\" working,\n",
|
||||
"as it outputs a value close to the actual value, but some race conditions make it\n",
|
||||
"undeterministic.\n",
|
||||
"3. The third one is a fixed parallelized version of the previous numba implementation,\n",
|
||||
"fully adapted to work in parallel"
|
||||
],
|
||||
"metadata": {
|
||||
"collapsed": false,
|
||||
"pycharm": {
|
||||
"name": "#%% md\n"
|
||||
}
|
||||
}
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"def plain_python_recurrence_histogram(ts: np.ndarray, epsilon: float, t_max: int):\n",
|
||||
" distances_histogram = np.zeros(len(ts))\n",
|
||||
" epsilon = 0.25\n",
|
||||
" for i in tqdm(np.arange(len(ts))):\n",
|
||||
" # finding the first \"out of ball\" index\n",
|
||||
" first_out = len(ts) # security\n",
|
||||
" for j in np.arange(i + 1, len(ts)):\n",
|
||||
" d = norm(ts[i] - ts[j])\n",
|
||||
" if d > epsilon:\n",
|
||||
" first_out = j\n",
|
||||
" break\n",
|
||||
"\n",
|
||||
" # finding the first \"back to the ball\" index\n",
|
||||
" for j in np.arange(first_out + 1, len(ts)):\n",
|
||||
" d = norm(ts[i] - ts[j])\n",
|
||||
" if d < epsilon:\n",
|
||||
" distances_histogram[j - i] += 1\n",
|
||||
" break\n",
|
||||
" return distances_histogram"
|
||||
],
|
||||
"metadata": {
|
||||
"collapsed": false,
|
||||
"pycharm": {
|
||||
"name": "#%%\n"
|
||||
}
|
||||
}
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 44,
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"from numba import jit, float32, prange, int32"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 42,
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"@jit(int32[:](float32[:,:], float32, int32), nopython=True)\n",
|
||||
"def recurrence_histogram(ts: np.ndarray, epsilon: float, t_max: int):\n",
|
||||
" distances_histogram = np.zeros(len(ts), dtype=np.int32)\n",
|
||||
" for i in np.arange(len(ts)):\n",
|
||||
" # finding the first \"out of ball\" index\n",
|
||||
" first_out = len(ts) # security\n",
|
||||
" for j in np.arange(i + 1, len(ts)):\n",
|
||||
" if t_max > 0 and j - i > t_max:\n",
|
||||
" break\n",
|
||||
" d = norm(ts[i] - ts[j])\n",
|
||||
" if d > epsilon:\n",
|
||||
" first_out = j\n",
|
||||
" break\n",
|
||||
"\n",
|
||||
" # finding the first \"back to the ball\" index\n",
|
||||
" for j in np.arange(first_out + 1, len(ts)):\n",
|
||||
" if 0 < t_max < j - i:\n",
|
||||
" break\n",
|
||||
" d = norm(ts[i] - ts[j])\n",
|
||||
" if d < epsilon:\n",
|
||||
" distances_histogram[j - i] += 1\n",
|
||||
" break\n",
|
||||
" return distances_histogram"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"# Here, we naively just switched out all \"np.arange\" occurences to numba.prange.\n",
|
||||
"@jit(int32[:](float32[:,:], float32, int32), nopython=True)\n",
|
||||
"def wrong_parallel_recurrence_histogram(ts: np.ndarray, epsilon: float, t_max: int):\n",
|
||||
" distances_histogram = np.zeros(len(ts), dtype=np.int32)\n",
|
||||
" for i in prange(len(ts)):\n",
|
||||
" # finding the first \"out of ball\" index\n",
|
||||
" first_out = len(ts) # security\n",
|
||||
" for j in prange(i + 1, len(ts)):\n",
|
||||
" if t_max > 0 and j - i > t_max:\n",
|
||||
" break\n",
|
||||
" d = norm(ts[i] - ts[j])\n",
|
||||
" if d > epsilon:\n",
|
||||
" first_out = j\n",
|
||||
" break\n",
|
||||
"\n",
|
||||
" # finding the first \"back to the ball\" index\n",
|
||||
" for j in prange(first_out + 1, len(ts)):\n",
|
||||
" if 0 < t_max < j - i:\n",
|
||||
" break\n",
|
||||
" d = norm(ts[i] - ts[j])\n",
|
||||
" if d < epsilon:\n",
|
||||
" distances_histogram[j - i] += 1\n",
|
||||
" break\n",
|
||||
" return distances_histogram"
|
||||
],
|
||||
"metadata": {
|
||||
"collapsed": false,
|
||||
"pycharm": {
|
||||
"name": "#%%\n"
|
||||
}
|
||||
}
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 46,
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"@jit(int32[:](float32[:,:], float32, int32), parallel=True, nopython=True)\n",
|
||||
"def parallel_recurrence_histogram(ts: np.ndarray, epsilon: float, t_max: int):\n",
|
||||
" return_distances = np.zeros(len(ts), dtype=np.int32)\n",
|
||||
" for i in prange(len(ts)):\n",
|
||||
" # finding the first \"out of ball\" index\n",
|
||||
" first_out = len(ts) # security\n",
|
||||
" for j in np.arange(i + 1, len(ts)):\n",
|
||||
" if t_max > 0 and j - i > t_max:\n",
|
||||
" break\n",
|
||||
" d = norm(ts[i] - ts[j])\n",
|
||||
" if d > epsilon:\n",
|
||||
" first_out = j\n",
|
||||
" break\n",
|
||||
"\n",
|
||||
" # finding the first \"back to the ball\" index\n",
|
||||
" for j in np.arange(first_out + 1, len(ts)):\n",
|
||||
" if t_max > 0 and j - i > t_max:\n",
|
||||
" break\n",
|
||||
" d = norm(ts[i] - ts[j])\n",
|
||||
" if d < epsilon:\n",
|
||||
" return_distances[i] = j - i\n",
|
||||
" break\n",
|
||||
" \n",
|
||||
" # building histogram, can't be parallel\n",
|
||||
" distance_histogram = np.zeros(len(ts), dtype=np.int32)\n",
|
||||
" for i in np.arange(len(ts)):\n",
|
||||
" if return_distances[i] != 0:\n",
|
||||
" distance_histogram[return_distances[i]] += 1\n",
|
||||
" return distance_histogram"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"source": [
|
||||
"### Running stuff\n",
|
||||
"\n",
|
||||
"Let's run the good ol' python implementation, just to take a look at what does\n",
|
||||
"the algorithm's output look like.\n",
|
||||
"\n",
|
||||
"Also, let's check that the various numba implementations's outputs match the plain\n",
|
||||
"python \"ground truth\"."
|
||||
],
|
||||
"metadata": {
|
||||
"collapsed": false,
|
||||
"pycharm": {
|
||||
"name": "#%% md\n"
|
||||
}
|
||||
}
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 18,
|
||||
"metadata": {},
|
||||
"outputs": [
|
||||
{
|
||||
"name": "stdout",
|
||||
"output_type": "stream",
|
||||
"text": [
|
||||
"CPU times: user 1min 30s, sys: 67.9 ms, total: 1min 30s\n",
|
||||
"Wall time: 55.6 s\n"
|
||||
]
|
||||
},
|
||||
{
|
||||
"data": {
|
||||
"text/plain": [
|
||||
"array([0, 0, 0, ..., 0, 0, 0], dtype=int32)"
|
||||
]
|
||||
},
|
||||
"execution_count": 18,
|
||||
"metadata": {},
|
||||
"output_type": "execute_result"
|
||||
}
|
||||
],
|
||||
"source": [
|
||||
"distances_histogram = plain_python_recurrence_histogram(ts, 0.12, 10000)\n",
|
||||
"px.line(distances_histogram, x=range(len(distances_histogram)), y=distances_histogram)"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "markdown",
|
||||
"source": [
|
||||
"### The ~~hour~~ seconds of truth"
|
||||
],
|
||||
"metadata": {
|
||||
"collapsed": false
|
||||
}
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"%time plain_python_recurrence_histogram(ts, 0.12, 10000)"
|
||||
],
|
||||
"metadata": {
|
||||
"collapsed": false,
|
||||
"pycharm": {
|
||||
"name": "#%%\n"
|
||||
}
|
||||
}
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"%time numba_recurrence_histogram(ts, 0.12, 10000)"
|
||||
],
|
||||
"metadata": {
|
||||
"collapsed": false,
|
||||
"pycharm": {
|
||||
"name": "#%%\n"
|
||||
}
|
||||
}
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": null,
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"%time parallel_recurrence_histogram(ts, 0.12, 10000)"
|
||||
],
|
||||
"metadata": {
|
||||
"collapsed": false,
|
||||
"pycharm": {
|
||||
"name": "#%%\n"
|
||||
}
|
||||
}
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 52,
|
||||
"metadata": {
|
||||
"pycharm": {
|
||||
"name": "#%%\n"
|
||||
}
|
||||
},
|
||||
"outputs": [
|
||||
{
|
||||
"name": "stdout",
|
||||
"output_type": "stream",
|
||||
"text": [
|
||||
"CPU times: user 16min 10s, sys: 3.01 s, total: 16min 13s\n",
|
||||
"Wall time: 6min 59s\n"
|
||||
]
|
||||
},
|
||||
{
|
||||
"data": {
|
||||
"text/plain": [
|
||||
"array([0, 0, 0, ..., 0, 0, 0], dtype=int32)"
|
||||
]
|
||||
},
|
||||
"execution_count": 52,
|
||||
"metadata": {},
|
||||
"output_type": "execute_result"
|
||||
}
|
||||
],
|
||||
"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.7.6"
|
||||
}
|
||||
},
|
||||
"nbformat": 4,
|
||||
"nbformat_minor": 4
|
||||
}
|
||||
Reference in New Issue
Block a user