mlp_small_2500.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Preamble"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "import os\n",
    "import numpy as np                                       # fast vectors and matrices\n",
    "import matplotlib.pyplot as plt                          # plotting\n",
    "from scipy.fftpack import fft\n",
    "\n",
    "from intervaltree import Interval,IntervalTree\n",
    "\n",
    "from time import time\n",
    "\n",
    "import tensorflow as tf\n",
    "\n",
    "from sklearn.metrics import precision_recall_curve\n",
    "from sklearn.metrics import average_precision_score\n",
    "\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Constants"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "d = 2048        # input dimensions\n",
    "m = 128         # number of notes\n",
    "fs = 44100      # samples/second\n",
    "features = 0\n",
    "labels = 1\n",
    "\n",
    "restore_weights = True\n",
    "folder = 'mlp_small_2500/'\n",
    "\n",
    "musicnet = os.environ['MUSICNET']"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Load MusicNet"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "train_data = dict(np.load(open(musicnet,'rb')))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# split our the test set\n",
    "test_data = dict()\n",
    "for id in (2303,2382,1819): # test set\n",
    "    test_data[str(id)] = train_data.pop(str(id))\n",
    "    \n",
    "train_ids = train_data.keys()\n",
    "test_ids = test_data.keys()\n",
    "    \n",
    "print len(train_data)\n",
    "print len(test_data)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Create the test set"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# create the test set\n",
    "Xtest = np.empty([3*7500,d])\n",
    "Ytest = np.zeros([3*7500,m])\n",
    "for i in range(len(test_ids)):\n",
    "    for j in range(7500):\n",
    "        index = fs+j*512 # start from one second to give us some wiggle room for larger segments\n",
    "        Xtest[7500*i + j] = test_data[test_ids[i]][features][index:index+d]\n",
    "        Xtest[7500*i + j] /= np.linalg.norm(Xtest[7500*i + j]) + 10e-6\n",
    "        \n",
    "        # label stuff that's on in the center of the window\n",
    "        for label in test_data[test_ids[i]][labels][index+d/2]:\n",
    "            Ytest[7500*i + j,label.data[1]] = 1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Raw MLP"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "tf.reset_default_graph()\n",
    "tf.set_random_seed(999)\n",
    "\n",
    "k = 2500\n",
    "\n",
    "x = tf.placeholder(tf.float32, shape=[None,d])\n",
    "y_ = tf.placeholder(tf.float32, shape=[None, m])\n",
    "\n",
    "scale = 10e-7\n",
    "\n",
    "w = tf.Variable(scale*tf.random_normal([d,k],seed=999))\n",
    "beta = tf.Variable(scale*tf.random_normal([k,m],seed=999)) # x100\n",
    "\n",
    "zx = tf.log(tf.constant(1.) + tf.nn.relu(tf.matmul(x,w)))\n",
    "y = tf.matmul(zx,beta)\n",
    "L = tf.reduce_mean(tf.nn.l2_loss(y-y_))\n",
    "\n",
    "init = tf.global_variables_initializer()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "if restore_weights:\n",
    "    weights = np.load(open(folder + 'w.npy','rb'))\n",
    "    coefficients = np.load(open(folder + 'beta.npy','rb'))\n",
    "    average_precision = list(np.load(open(folder + 'ap.npy','rb')))\n",
    "    square_error = list(np.load(open(folder + 'loss.npy','rb')))\n",
    "    weights_top = list(np.load(open(folder + 'wtop.npy','rb')))\n",
    "    weights_bottom = list(np.load(open(folder + 'wbot.npy','rb')))\n",
    "\n",
    "    sess = tf.Session()\n",
    "    sess.run(w.assign(weights))\n",
    "    sess.run(beta.assign(coefficients))\n",
    "    init_lr = .01/16\n",
    "else:\n",
    "    square_error = []\n",
    "    weights_top = []\n",
    "    weights_bottom = []\n",
    "    average_precision = []\n",
    "    sess = tf.Session()\n",
    "    sess.run(init)\n",
    "    init_lr = .01"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false,
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "lr =  init_lr\n",
    "opt = tf.train.GradientDescentOptimizer(lr)\n",
    "train_step = opt.minimize(L)\n",
    "Xmb = np.empty([len(train_data),d])\n",
    "np.random.seed(999)\n",
    "start = time()\n",
    "print 'iter\\tsquare_loss\\tweights_top\\tweights_bottom\\tavg_precision\\ttime'\n",
    "for i in xrange(500000):\n",
    "    if i % 100 == 0 and (i != 0 or len(square_error) == 0):\n",
    "        square_error.append(sess.run(L, feed_dict={x: Xtest, y_: Ytest})/Xtest.shape[0])\n",
    "        weights_top.append(np.mean(np.linalg.norm(beta.eval(session=sess),axis=0)))\n",
    "        weights_bottom.append(np.mean(np.linalg.norm(w.eval(session=sess),axis=0)))\n",
    "        \n",
    "        Yhattestbase = sess.run(y,feed_dict={x: Xtest})\n",
    "        yflat = Ytest.reshape(Ytest.shape[0]*Ytest.shape[1])\n",
    "        yhatflat = Yhattestbase.reshape(Yhattestbase.shape[0]*Yhattestbase.shape[1])\n",
    "        average_precision.append(average_precision_score(yflat, yhatflat))\n",
    "        \n",
    "        end = time()\n",
    "        print i,'\\t', round(square_error[-1],8),\\\n",
    "                '\\t', round(weights_top[-1],8),\\\n",
    "                '\\t', round(weights_bottom[-1],8),\\\n",
    "                '\\t', round(average_precision[-1],8),\\\n",
    "                '\\t', round(end-start,8)\n",
    "        start = time()\n",
    "    \n",
    "    Ymb = np.zeros([len(train_data),m])\n",
    "    for j in range(len(train_ids)):\n",
    "        s = np.random.randint(d/2,len(train_data[train_ids[j]][features])-d/2)\n",
    "        Xmb[j] = train_data[train_ids[j]][features][s-d/2:s+d/2]\n",
    "        Xmb[j] /= np.linalg.norm(Xmb[j]) + 10e-6\n",
    "        for label in train_data[train_ids[j]][labels][s]:\n",
    "            Ymb[j,label.data[1]] = 1\n",
    "    \n",
    "    sess.run(train_step, feed_dict={x: Xmb, y_: Ymb})"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "burnin=1\n",
    "fig, ((ax1, ax2),(ax3,ax4)) = plt.subplots(2, 2)\n",
    "fig.set_figwidth(12)\n",
    "fig.set_figheight(10)\n",
    "ax1.set_title('average precision')\n",
    "ax1.plot(average_precision[burnin:],color='g')\n",
    "ax2.set_title('square loss')\n",
    "ax2.plot(square_error[burnin:],color='g')\n",
    "ax3.set_title('weights top')\n",
    "ax3.plot(weights_top[burnin:],color='g')\n",
    "ax4.set_title('weights bottom')\n",
    "ax4.plot(weights_bottom[burnin:],color='g')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "offset = 20\n",
    "window = 2048\n",
    "f, ax = plt.subplots(20,3, sharey=False)\n",
    "f.set_figheight(20)\n",
    "f.set_figwidth(20)\n",
    "weights = w.eval(session=sess)\n",
    "for i in range(20):\n",
    "    ax[i,0].plot(weights[10:-10,i+offset], color=(41/255.,104/255.,168/255.))\n",
    "    ax[i,0].set_xlim([-20,d+20])\n",
    "    ax[i,0].set_ylim([np.min(weights[10:-10,i+offset]),np.max(weights[10:-10,i+offset])])\n",
    "    ax[i,0].set_xticklabels([])\n",
    "    ax[i,0].set_yticklabels([])\n",
    "    ax[i,1].plot(weights[d/2-1024:d/2+1024,i+offset], color=(41/255.,104/255.,168/255.))\n",
    "    ax[i,1].set_xticklabels([])\n",
    "    ax[i,1].set_yticklabels([])\n",
    "    ax[i,1].set_xlim([0,d])\n",
    "    ax[i,2].plot(np.abs(fft(weights[:,i+offset]))[0:200], color=(41/255.,104/255.,168/255.))\n",
    "    ax[i,2].set_xticklabels([])\n",
    "    ax[i,2].set_yticklabels([])\n",
    "    \n",
    "for i in range(ax.shape[0]):\n",
    "    for j in range(ax.shape[1]):\n",
    "        ax[i,j].set_xticks([])\n",
    "        ax[i,j].set_yticks([])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "folder = 'mlp_small_2500/'\n",
    "\n",
    "weights = w.eval(session=sess)\n",
    "coefficients = beta.eval(session=sess)\n",
    "\n",
    "np.save(open(folder + 'w.npy','wb'),weights)\n",
    "np.save(open(folder + 'beta.npy','wb'),coefficients)\n",
    "np.save(open(folder + 'ap.npy','wb'),average_precision)\n",
    "np.save(open(folder + 'loss.npy','wb'),square_error)\n",
    "np.save(open(folder + 'wtop.npy','wb'),weights_top)\n",
    "np.save(open(folder + 'wbot.npy','wb'),weights_bottom)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# MIREX evaluation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import mir_eval"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "Xvalidation = np.zeros([50*len(train_data),d])\n",
    "Yvalidation = np.zeros([50*len(train_data),m])\n",
    "for i in range(len(train_data)):\n",
    "    # 50 random samples from each recording\n",
    "    s = np.random.randint(d/2,len(train_data[train_ids[i]][features])-d/2,size=100)\n",
    "    for j in range(50):\n",
    "        Xvalidation[50*i+j] = train_data[train_ids[i]][features][s[j]-d/2:s[j]+d/2]\n",
    "        Xvalidation[50*i+j] /= np.linalg.norm(Xvalidation[50*i+j]) + 10e-6\n",
    "        # label stuff that's on in the center of the window\n",
    "        for label in train_data[train_ids[i]][labels][s[j]]:\n",
    "            Yvalidation[50*i+j,label.data[1]] = 1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "Yhattestbase = sess.run(y,feed_dict={x: Xtest})\n",
    "yflat = Ytest.reshape(Ytest.shape[0]*Ytest.shape[1])\n",
    "yhatflat = Yhattestbase.reshape(Yhattestbase.shape[0]*Yhattestbase.shape[1])\n",
    "print average_precision_score(yflat, yhatflat)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "Yhatbase = estimate(Xvalidation)\n",
    "\n",
    "# single threshold\n",
    "density = 500\n",
    "P = np.empty(density)\n",
    "R = np.empty(density)\n",
    "F = np.empty(density)\n",
    "for i in np.arange(density):\n",
    "    if i % 100 == 0: print '.',\n",
    "    c = i/float(density)\n",
    "    Yhat = Yhatbase>c\n",
    "    true_positives = np.sum(Yhat*Yvalidation)\n",
    "    P[i] = true_positives/np.sum(Yhat)\n",
    "    R[i] = true_positives/np.sum(Yvalidation)\n",
    "    F[i] = 2*(P[i]*R[i])/(P[i]+R[i])\n",
    "\n",
    "plt.plot(F)\n",
    "i = np.argmax(F)\n",
    "c = i/float(density)\n",
    "print c"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "Yhatbase = estimate(Xtest)\n",
    "\n",
    "Yhat = Yhatbase>c\n",
    "Yhatlist = []\n",
    "Ytestlist = []\n",
    "for i in range(len(Yhat)):\n",
    "    fhat = []\n",
    "    ftest = []\n",
    "    for note in range(128):\n",
    "        if Yhat[i][note] == 1:\n",
    "            fhat.append(440.*2**((note - 69.)/12.))\n",
    "        if Ytest[i][note] == 1:\n",
    "            ftest.append(440.*2**((note - 69.)/12.))\n",
    "    Yhatlist.append(np.array(fhat))\n",
    "    Ytestlist.append(np.array(ftest))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "P,R,Acc,Esub,Emiss,Efa,Etot,cP,cR,cAcc,cEsub,cEmiss,cEfa,cEtot = \\\n",
    "mir_eval.multipitch.metrics(np.arange(len(Ytestlist))/100.,Ytestlist,np.arange(len(Yhatlist))/100.,Yhatlist)\n",
    "\n",
    "print P\n",
    "print R\n",
    "print Acc\n",
    "print Etot\n",
    "print Esub\n",
    "print Emiss\n",
    "print Efa\n",
    "\n",
    "print '-----'\n",
    "\n",
    "print cP\n",
    "print cR\n",
    "print cAcc\n",
    "print cEtot\n",
    "print cEsub\n",
    "print cEmiss\n",
    "print cEfa"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Precision/Recall"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "fs = 44100\n",
    "m = 128\n",
    "\n",
    "composition = 2\n",
    "            \n",
    "# create the test set\n",
    "Xtest = np.empty([7500,d])\n",
    "Ytest = np.zeros([7500,m])\n",
    "for i in range(1):\n",
    "    for j in range(7500):\n",
    "        index = fs+j*512 # start from one second to give us some wiggle room for larger segments\n",
    "        Xtest[j] = test_data[test_ids[composition]][features][index:index+d]\n",
    "        Xtest[j] /= np.linalg.norm(Xtest[7500*i + j]) + 10e-6\n",
    "        \n",
    "        # label stuff that's on in the center of the window\n",
    "        for label in test_data[test_ids[composition]][labels][index+d/2]:\n",
    "            Ytest[7500*i + j,label.data[1]] = 1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "print test_ids[composition]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "Yhatbase = estimate(Xvalidation)\n",
    "\n",
    "# single threshold\n",
    "density = 500\n",
    "P = np.empty(density)\n",
    "R = np.empty(density)\n",
    "F = np.empty(density)\n",
    "for i in np.arange(density):\n",
    "    if i % 100 == 0: print '.',\n",
    "    c = i/float(density)\n",
    "    Yhat = Yhatbase>c\n",
    "    true_positives = np.sum(Yhat*Yvalidation)\n",
    "    P[i] = true_positives/np.sum(Yhat)\n",
    "    R[i] = true_positives/np.sum(Yvalidation)\n",
    "    F[i] = 2*(P[i]*R[i])/(P[i]+R[i])\n",
    "\n",
    "plt.plot(F)\n",
    "i = np.argmax(F)\n",
    "c = i/float(density)\n",
    "print c"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "plt.plot(R,P)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "Yhattestbase = estimate(Xtest)\n",
    "Yhat = Yhattestbase>c\n",
    "true_positives = np.sum(Yhat*Ytest)\n",
    "P = true_positives/(np.sum(Yhat))\n",
    "R = true_positives/(np.sum(Ytest))\n",
    "F = 2*(P*R)/(P+R)\n",
    "print P\n",
    "print R\n",
    "print F"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "yflat = Ytest.reshape(Ytest.shape[0]*Ytest.shape[1])\n",
    "yhatflat = Yhattestbase.reshape(Yhattestbase.shape[0]*Yhattestbase.shape[1])\n",
    "\n",
    "precision, recall, _ = precision_recall_curve(yflat, yhatflat)\n",
    "ap = average_precision_score(yflat, yhatflat)\n",
    "plt.plot(recall,precision)\n",
    "print ap"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "plt.rcParams.update({'font.size': 12})\n",
    "\n",
    "fig = plt.figure()\n",
    "plt.plot(recall,precision)\n",
    "fig.axes[0].set_xlabel('recall')\n",
    "fig.axes[0].set_ylabel('precision')\n",
    "\n",
    "plt.tight_layout()\n",
    "plt.savefig('mlp2500_pr.eps',format='eps', dpi=1000)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}