{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Independent Component Analysis Lab\n",
    "\n",
    "In this notebook, we'll use Independent Component Analysis to retrieve original signals from three observations each of which contains a different mix of the original signals. This is the same problem explained in the ICA video.\n",
    "\n",
    "## Dataset\n",
    "Let's begin by looking at the dataset we have. We have three WAVE files, each of which is a mix, as we've mentioned. If you haven't worked with audio files in python before, that's okay, they basically boil down to being lists of floats.\n",
    "\n",
    "\n",
    "Let's begin by loading our first audio file, **[ICA mix 1.wav](ICA mix 1.wav)** [click to listen to the file]:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import wave\n",
    "\n",
    "# Read the wave file\n",
    "mix_1_wave = wave.open('ICA mix 1.wav','r')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's peak at the parameters of the wave file to learn more about it"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "_wave_params(nchannels=1, sampwidth=2, framerate=44100, nframes=264515, comptype='NONE', compname='not compressed')"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "mix_1_wave.getparams()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "So this file has only channel (so it's mono sound). It has a frame rate of 44100, which means each second of sound is represented by 44100 integers (integers because the file is in the common PCM 16-bit format). The file has a total of 264515 integers/frames, which means its length in seconds is:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "5.998072562358277"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "264515/44100"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's extract the frames of the wave file, which will be a part of the dataset we'll run ICA against:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/opt/conda/lib/python3.6/site-packages/ipykernel_launcher.py:3: DeprecationWarning: Numeric-style type codes are deprecated and will result in an error in the future.\n",
      "  This is separate from the ipykernel package so we can avoid doing imports until\n"
     ]
    }
   ],
   "source": [
    "# Extract Raw Audio from Wav File\n",
    "signal_1_raw = mix_1_wave.readframes(-1)\n",
    "signal_1 = np.fromstring(signal_1_raw, 'Int16')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "signal_1 is now a list of ints representing the sound contained in the first file."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "('length: ',\n",
       " 264515,\n",
       " 'first 100 elements: ',\n",
       " array([  879,  1268,  1460,  1756,  1943,  2216,  2407,  2668,  2866,\n",
       "         3106,  3308,  3546,  3752,  3981,  4175,  4395,  4588,  4790,\n",
       "         4966,  5146,  5292,  5436,  5550,  5643,  5717,  5759,  5790,\n",
       "         5798,  5789,  5756,  5713,  5649,  5576,  5478,  5381,  5267,\n",
       "         5146,  4999,  4856,  4682,  4502,  4308,  4097,  3875,  3637,\n",
       "         3380,  3107,  2825,  2514,  2194,  1847,  1472,  1087,   671,\n",
       "          227,  -219,  -691, -1176, -1666, -2167, -2669, -3179, -3668,\n",
       "        -4170, -4643, -5116, -5559, -5985, -6380, -6765, -7105, -7422,\n",
       "        -7706, -7955, -8163, -8339, -8470, -8557, -8600, -8618, -8585,\n",
       "        -8524, -8425, -8298, -8129, -7947, -7720, -7475, -7205, -6916,\n",
       "        -6606, -6266, -5922, -5556, -5165, -4774, -4353, -3922, -3476, -3021], dtype=int16))"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "'length: ', len(signal_1) , 'first 100 elements: ',signal_1[:100]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "If we plot this array as a line graph, we'll get the familiar wave form representation:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "[encoded data, 21588 characters]\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f544220ada0>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "\n",
    "fs = mix_1_wave.getframerate()\n",
    "timing = np.linspace(0, len(signal_1)/fs, num=len(signal_1))\n",
    "\n",
    "\n",
    "plt.figure(figsize=(12,2))\n",
    "plt.title('Recording 1')\n",
    "plt.plot(timing,signal_1, c=\"#3ABFE7\")\n",
    "plt.ylim(-35000, 35000)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In the same way, we can now load the other two wave files, **[ICA mix 2.wav](ICA mix 2.wav)** and **[ICA mix 3.wav](ICA mix 3.wav)**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/opt/conda/lib/python3.6/site-packages/ipykernel_launcher.py:6: DeprecationWarning: Numeric-style type codes are deprecated and will result in an error in the future.\n",
      "  \n",
      "/opt/conda/lib/python3.6/site-packages/ipykernel_launcher.py:13: DeprecationWarning: Numeric-style type codes are deprecated and will result in an error in the future.\n",
      "  del sys.path[0]\n"
     ]
    },
    {
     "data": {
      "image/png": "[encoded data, 19760 characters]\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f54421ef390>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "[encoded data, 17016 characters]\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f54421a1f98>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "\n",
    "mix_2_wave = wave.open('ICA mix 2.wav','r')\n",
    "\n",
    "#Extract Raw Audio from Wav File\n",
    "signal_raw_2 = mix_2_wave.readframes(-1)\n",
    "signal_2 = np.fromstring(signal_raw_2, 'Int16')\n",
    "\n",
    "\n",
    "mix_3_wave = wave.open('ICA mix 3.wav','r')\n",
    "\n",
    "#Extract Raw Audio from Wav File\n",
    "signal_raw_3 = mix_3_wave.readframes(-1)\n",
    "signal_3 = np.fromstring(signal_raw_3, 'Int16')\n",
    "\n",
    "\n",
    "plt.figure(figsize=(12,2))\n",
    "plt.title('Recording 2')\n",
    "plt.plot(timing,signal_2, c=\"#3ABFE7\")\n",
    "plt.ylim(-35000, 35000)\n",
    "plt.show()\n",
    "\n",
    "\n",
    "plt.figure(figsize=(12,2))\n",
    "plt.title('Recording 3')\n",
    "plt.plot(timing,signal_3, c=\"#3ABFE7\")\n",
    "plt.ylim(-35000, 35000)\n",
    "plt.show()\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now that we've read all three files, we're ready to [zip](https://docs.python.org/3/library/functions.html#zip) them to create our dataset.\n",
    "\n",
    "* Create dataset ```X``` by zipping signal_1, signal_2, and signal_3 into a single list"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[(879, 157, 879),\n",
       " (1268, 899, 2584),\n",
       " (1460, 1610, 4081),\n",
       " (1756, 2325, 5682),\n",
       " (1943, 3022, 7147),\n",
       " (2216, 3718, 8689),\n",
       " (2407, 4386, 10114),\n",
       " (2668, 5056, 11583),\n",
       " (2866, 5690, 12954),\n",
       " (3106, 6318, 14327)]"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "X = list(zip(signal_1, signal_2, signal_3))\n",
    "\n",
    "# Let's peak at what X looks like\n",
    "X[:10]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We are now ready to run ICA to try to retrieve the original signals.\n",
    "\n",
    "* Import sklearn's [FastICA](http://scikit-learn.org/stable/modules/generated/sklearn.decomposition.FastICA.html) module\n",
    "* Initialize FastICA look for three components\n",
    "* Run the FastICA algorithm using fit_transform on dataset X"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [],
   "source": [
    "# TODO: Import FastICA\n",
    "from sklearn.decomposition import FastICA\n",
    "\n",
    "# TODO: Initialize FastICA with n_components=3\n",
    "ica = FastICA(n_components=3)\n",
    "\n",
    "# TODO: Run the FastICA algorithm using fit_transform on dataset X\n",
    "ica_result = ica.fit_transform(X)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "```ica_result``` now contains the result of FastICA, which we hope are the original signals. It's in the shape:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(264515, 3)"
      ]
     },
     "execution_count": 19,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "ica_result.shape"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's split into separate signals and look at them"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [],
   "source": [
    "result_signal_1 = ica_result[:,0]\n",
    "result_signal_2 = ica_result[:,1]\n",
    "result_signal_3 = ica_result[:,2]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's plot to see how the wave forms look"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "[encoded data, 21660 characters]\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f54325fc588>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "[encoded data, 12476 characters]\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f54325a6940>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "[encoded data, 24080 characters]\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f54325d94e0>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Plot Independent Component #1\n",
    "plt.figure(figsize=(12,2))\n",
    "plt.title('Independent Component #1')\n",
    "plt.plot(result_signal_1, c=\"#df8efd\")\n",
    "plt.ylim(-0.010, 0.010)\n",
    "plt.show()\n",
    "\n",
    "# Plot Independent Component #2\n",
    "plt.figure(figsize=(12,2))\n",
    "plt.title('Independent Component #2')\n",
    "plt.plot(result_signal_2, c=\"#87de72\")\n",
    "plt.ylim(-0.010, 0.010)\n",
    "plt.show()\n",
    "\n",
    "# Plot Independent Component #3\n",
    "plt.figure(figsize=(12,2))\n",
    "plt.title('Independent Component #3')\n",
    "plt.plot(result_signal_3, c=\"#f65e97\")\n",
    "plt.ylim(-0.010, 0.010)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Do some of these look like musical wave forms? \n",
    "\n",
    "The best way to confirm the result is to listen to resulting files. So let's save as wave files and verify. But before we do that, we'll have to:\n",
    "* convert them to integer (so we can save as PCM 16-bit Wave files), otherwise only some media players would be able to play them and others won't\n",
    "* Map the values to the appropriate range for int16 audio. That range is between -32768 and +32767. A basic mapping can be done by multiplying by 32767.\n",
    "* The sounds will be a little faint, we can increase the volume by multiplying by a value like 100"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [],
   "source": [
    "from scipy.io import wavfile\n",
    "\n",
    "# Convert to int, map the appropriate range, and increase the volume a little bit\n",
    "result_signal_1_int = np.int16(result_signal_1*32767*100)\n",
    "result_signal_2_int = np.int16(result_signal_2*32767*100)\n",
    "result_signal_3_int = np.int16(result_signal_3*32767*100)\n",
    "\n",
    "\n",
    "# Write wave files\n",
    "wavfile.write(\"result_signal_1.wav\", fs, result_signal_1_int)\n",
    "wavfile.write(\"result_signal_2.wav\", fs, result_signal_2_int)\n",
    "wavfile.write(\"result_signal_3.wav\", fs, result_signal_3_int)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The resulting files we have now are: [note: make sure to lower the volume on your speakers first, just in case some problem caused the file to sound like static]\n",
    "* [result_signal_1.wav](result_signal_1.wav)\n",
    "* [result_signal_2.wav](result_signal_2.wav)\n",
    "* [result_signal_3.wav](result_signal_3.wav)\n",
    "\n",
    "\n",
    "\n",
    "Music:\n",
    "* Piano - The Carnival of the Animals - XIII. The Swan (Solo piano version). Performer: Markus Staab\n",
    "* Cello - Cello Suite no. 3 in C, BWV 1009 - I. Prelude. Performer: European Archive"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "anaconda-cloud": {},
  "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.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}
