16_11_22 perplexity test.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Populating the interactive namespace from numpy and matplotlib\n"
     ]
    }
   ],
   "source": [
    "%pylab inline"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Perplexity experiment\n",
    "\n",
    "Trying to get my head around perplexity as a measure. This notebook does the following:\n",
    "\n",
    "* Generate data from a ground truth unigram distribution (p(0)=0.9, p(1)=0.09, p(2)=0.01)\n",
    "* Estimate the ground truth n-gram distributions from the data\n",
    "* Interpret the n-gram distribution as a posterior distribution conditioned on n-1 (for 3-gram, P(x[n] | x[n-1], x[n-2]))\n",
    "* Measure the perplexity of the sequence as described by each posterior"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 78,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0\n",
      " 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 1 0 0 0 0 0 0 0\n",
      " 1 1 0 0 0 0 0 0 0 0 0 2 0 0 0 0 1 0 0 0 0 1 0 0 0 0]\n"
     ]
    }
   ],
   "source": [
    "# Create sparse data\n",
    "true_unigram = [0.9, 0.09, 0.01]\n",
    "seq_len = 10000\n",
    "sequence = np.random.choice(len(true_unigram), seq_len, p=true_unigram)\n",
    "print sequence[:100]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 79,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "1-gram distribution\n",
      "(0,): 0.8972\n",
      "(1,): 0.0922\n",
      "(2,): 0.0106\n",
      "2-gram distribution\n",
      "(0, 0): 0.804680468047\n",
      "(1, 0): 0.0831083108311\n",
      "(0, 1): 0.0828082808281\n",
      "(0, 2): 0.00970097009701\n",
      "(2, 0): 0.00940094009401\n",
      "(1, 1): 0.00820082008201\n",
      "(2, 1): 0.001200120012\n",
      "(1, 2): 0.000900090009001\n",
      "3-gram distribution\n",
      "(0, 0, 0): 0.72274454891\n",
      "(0, 1, 0): 0.0746149229846\n",
      "(0, 0, 1): 0.0739147829566\n",
      "(1, 0, 0): 0.0736147229446\n",
      "(0, 2, 0): 0.00850170034007\n",
      "(2, 0, 0): 0.00830166033207\n",
      "(1, 0, 1): 0.00820164032807\n",
      "(0, 0, 2): 0.00800160032006\n",
      "(0, 1, 1): 0.00740148029606\n",
      "(1, 1, 0): 0.00730146029206\n",
      "(1, 0, 2): 0.00130026005201\n",
      "(2, 1, 0): 0.00120024004801\n",
      "(0, 2, 1): 0.00120024004801\n",
      "(1, 2, 0): 0.000900180036007\n",
      "(0, 1, 2): 0.000800160032006\n",
      "(1, 1, 1): 0.000800160032006\n",
      "(2, 0, 1): 0.000700140028006\n",
      "(2, 0, 2): 0.000400080016003\n",
      "(1, 1, 2): 0.000100020004001\n"
     ]
    }
   ],
   "source": [
    "# Calculate n-gram distrbutions from sample of ground truth unigram distribution\n",
    "max_ngram = 6\n",
    "n_to_distribution = {}\n",
    "for n in xrange(1, max_ngram + 1):\n",
    "    ngrams = []\n",
    "    for i in xrange(seq_len):\n",
    "        ngram = sequence[i:i + n]\n",
    "        if len(ngram) != n:\n",
    "            break\n",
    "        ngrams.append(ngram)\n",
    "    ngram_to_count = {}\n",
    "    for ngram in ngrams:\n",
    "        ngram = tuple(ngram)\n",
    "        if ngram not in ngram_to_count:\n",
    "            ngram_to_count[ngram] = 0\n",
    "        ngram_to_count[ngram] += 1\n",
    "    ngram_distribution = {i:float(c) / len(ngrams) for i, c in ngram_to_count.items()}\n",
    "    n_to_distribution[n] = ngram_distribution\n",
    "\n",
    "for n in xrange(1, min(4, max_ngram + 1)):\n",
    "    print '{}-gram distribution'.format(n)\n",
    "    for ngram, proportion in sorted(n_to_distribution[n].items(), key=lambda x: -x[1]):\n",
    "        print '{}: {}'.format(ngram, n_to_distribution[n][ngram])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 80,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# Calculate perplexity of n-gram distribution against the sample\n",
    "n_to_perplexity = {}\n",
    "for n in xrange(1, max_ngram + 1):\n",
    "    distribution = n_to_distribution[n]\n",
    "    ncontext = n - 1\n",
    "    log_prob_cum = 0.0\n",
    "    for i in xrange(seq_len):\n",
    "        ngram = tuple(sequence[i:i + n])\n",
    "        if len(ngram) != n:\n",
    "            break\n",
    "        posterior = distribution[ngram]\n",
    "        log_prob_cum += np.log(posterior)\n",
    "    negative_log_prob = -log_prob_cum / (len(sequence) - (n - 1))\n",
    "    perplexity = np.exp(negative_log_prob)\n",
    "    n_to_perplexity[n] = perplexity"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 81,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.text.Text at 0x90f76a0>"
      ]
     },
     "execution_count": 81,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlgAAAFHCAYAAAB9DjIPAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3XeUXXW5//H3kwQCCU1FQIUrBEF6R2rCcEFEingRFAtR\nroCgIEXlKj9L1AXXAiIIBgidS28iCoR2BwhCgHQggFIuRYq0hCSUJPP8/tgnpJBMTpKzZ58z836t\nNWvOzNnn7AfOInzy/T772ZGZSJIkqXF6VV2AJElSd2PAkiRJajADliRJUoMZsCRJkhrMgCVJktRg\nBixJkqQGKz1gRcRRETGh9vXdss8nSZJUtVIDVkRsCHwT2ArYDNgrIgaUeU5JkqSqlb2CtT4wMjPf\nycyZwF3AviWfU5IkqVJlB6yHgIER8YGI6AfsAaxR8jklSZIq1afMN8/MRyPi18CtwBRgDDCzzHNK\nkiRVLbryXoQRcQLwbGaeOc/vvSGiJElqGZkZnT1f6goWQER8ODP/FRH/BvwHsO38jvOm061pyJAh\nDBkypOoytJj8/FqXn11r8/NrbRGdZiugCwIWcE1EfBCYDnw7Myd3wTklSZIqU3rAysxBZZ9DkiSp\nmTjJXUukra2t6hK0BPz8WpefXWvz82tNU6bAMcfUd2yXNrkvsIiIbIY6JEmS5ueZZ+Bzn4MttoDz\nz4+FNrm7giVJktSJ++6D7baDAw+Ec8+t7zVd0eQuSZLUki67DI46Cs47D/baq/7XGbAkSZLm0dEB\nP/sZ/M//wO23w8YbL9rrDViSJElzmDoVvv51eOEFGDkSVlll0d/DHixJkqSa55+HQYOgf3+4447F\nC1dgwJIkSQLgwQdhm23gi1+ECy6Avn0X/73cIpQkST3elVfCd74Dw4bB5z+/5O9nwJIkST1WJvzi\nF8VVgrfeCptt1pj3NWBJkqQe6a234KCD4Omni2b21VZr3HvbgyVJknqcF16AnXaC3r2hvb2x4QoM\nWJIkqYcZPbpoZt9nn2LO1TLLNP4cbhFKkqQe45pr4LDD4Mwz4QtfKO88BixJktTtZcKJJxbBavjw\n4qbNZTJgSZKkbu3tt+Hgg+Gxx4pm9o9+tPxz2oMlSZK6rZdegp13hunT4c47uyZcgQFLkiR1U+PG\nFc3su+8Ol18O/fp13bndIpQkSd3O9dcX24Knnw5f+lLXn9+AJUmSuo1M+M1v4A9/gBtvhK23rqYO\nA5YkSeoW3nkHvvUtGD8e7rsPVl+9ulrswZIkSS3v5Zdhl13gzTfh7rurDVdgwJIkSS1uwoSimb2t\nDa66Cvr3r7oitwglSVIL+8tfihs2n3oqfOUrVVczmwFLkiS1nEw45RQ46SS44QbYdtuqK5qbAUuS\nJLWUd9+Fww+HUaOKZvZ/+7eqK3q/0nuwIuKYiHgoIsZHxCURsXTZ55QkSd3TK6/Apz8Nr74KI0Y0\nZ7iCkgNWRHwUOBLYIjM3oVgxO6DMc0qSpO7pkUeKZvbtt4drr4Xllqu6ogXrii3C3kD/iOgA+gH/\n7IJzSpKkbuTmm2Hw4KLnavDgqqtZuFJXsDLzn8DJwDPA88AbmXlbmeeUJEndRyacdlpxpeB117VG\nuILytwhXAvYBPg58FFguIproIkpJktSspk+Hww6DYcPg3nthhx2qrqh+ZW8R7go8mZmvAUTEtcD2\nwKXzHjhkyJD3Hre1tdHW1lZyaZIkqVm99hrst18xNPRvf4Pll6+ulvb2dtrb2xfpNZGZ5VQDRMSn\ngHOBrYF3gPOBBzLzjHmOyzLrkCRJreOxx2CvveDzn4df/Qp69666orlFBJkZnR1Tdg/W/cDVwBhg\nHBDA2WWeU5Ikta5bb4VBg+BHP4Lf/rb5wlW9Sl3BqrsIV7AkSerxzjgDfvlLuPLKImQ1q3pWsJzk\nLkmSKjVjBhx1FLS3F/1WAwZUXdGSM2BJkqTKvP46fPGL0KdPEa5WXLHqihqj9FvlSJIkzc/f/w7b\nbQcbbVTcsLm7hCswYEmSpArccQfsuCMccwycckqxgtWddLN/HEmS1OzOOgt++lO4/HLYeeeqqymH\nAUuSJHWJGTPg+98v7is4YgSss07VFZXHgCVJkko3aRIccADMnAn33QcrrVR1ReWyB0uSJJXqiSeK\nZva114Ybb+z+4QoMWJIkqUR33lncpPmII+D007tfM/uC9JB/TEmS1NXOPReOPx4uuQR23bXqarqW\nAUuSJDXUzJnwX/8Ff/4z3HUXfPKTVVfU9QxYkiSpYSZPhq98BaZNK5rZP/jBqiuqhj1YkiSpIZ5+\nuui3Wn11GD6854YrMGBJkqQGGDGiuFLw0ENh6FBYaqmqK6qWW4SSJGmJXHgh/OAHcNFFsPvuVVfT\nHAxYkiRpsXR0FFcJXnUVtLfDBhtUXVHzMGBJkqRFNmUKfO1r8MYbMHIkrLxy1RU1F3uwJEnSInnm\nmaKZfeWV4ZZbDFfzY8CSJEl1u/de2HZb+MY3YNgwWHrpqitqTm4RSpKkulxyCRxzDJx/Puy5Z9XV\nNDcDliRJ6lRHB/z0p0XAuuMO2GijqitqfgYsSZK0QFOnwuDB8NJLRTP7KqtUXVFrsAdLkiTN13PP\nwcCBsPzycPvthqtFYcCSJEnvc//9RTP7AQcUPVd9+1ZdUWtxi1CSJM3liivgiCPgnHNgn32qrqY1\nGbAkSRJQNLP/4hdwwQVw222w6aZVV9S6Sg1YEbEucAWQQAADgJ9k5mllnleSJC2aadPgoIPg2WeL\nZvZVV626otZWasDKzMeBzQEiohfwHHBdmeeUJEmL5p//LLYCP/nJYgzDMstUXVHr68om912BJzLz\n2S48pyRJ6sSoUbDNNvAf/wEXX2y4apSu7MH6EnBZF55PkiR14uqr4fDD4ayzYN99q66me4nMLP8k\nEUsB/wQ2yMx/zef57Io6JEkSZMIJJ8DZZ8P118Pmm1ddUWuJCDIzOjumq1awPguMml+4mmXIkCHv\nPW5ra6Otra38qiRJ6mHeegu++U144omimf0jH6m6oubX3t5Oe3v7Ir2mq1awLgNuzswLF/C8K1iS\nJJXsxRfh85+HNdcshocuu2zVFbWmelawSm9yj4h+FA3u15Z9LkmSNH9jxxbN7HvsAZddZrgqW5es\nYC20CFewJEkqzZ/+BIccAn/8I+y/f9XVtL5m6sGSJEldLBN+/Ws4/XS46SbYaquqK+o5DFiSJHVD\n77xTrFo9/HDRzP6xj1VdUc/SlYNGJUlSF3j5Zfj3fy+uGLz7bsNVFQxYkiR1I+PHF83su+wCV1wB\n/fpVXVHP5BahJEndxA03FDOuTj0Vvvzlqqvp2QxYkiS1uEw4+WQ45ZQiZG2zTdUVyYAlSVILe/dd\nOOwwGD0a7rsP1lij6ooE9mBJktSy/vUv2HVXeP11GDHCcNVMDFiSJLWghx8utgIHDoRrroHllqu6\nIs3JLUJJklrMTTfB179e9F0deGDV1Wh+DFiSJLWIzOIKwd/8prj9zfbbV12RFsSAJUlSC3j3XTji\niKKR/d574eMfr7oidcaAJUlSk3v1VdhvP1h+ebjnnuK7mptN7pIkNbFHHy2a2bfeGq67znDVKgxY\nkiQ1qVtugUGD4P/9v6LvqnfvqitSvdwilCSpCZ1+OpxwQjGCYeDAqqvRojJgSZLURKZPh6OOgrvu\ngr/9DdZaq+qKtDgMWJIkNYnXX4f994e+fYtwtcIKVVekxWUPliRJTeDxx2HbbWGTTeDPfzZctToD\nliRJFbv99qLP6gc/gN/9zmb27sCAJUlShc48E776VbjiCjj44KqrUaPYgyVJUgVmzIBjj4Vbb4UR\nI+ATn6i6IjWSAUuSpC72xhvwpS8Vj++9F1Zaqdp61HhuEUqS1IUefxy22w4++Un4618NV92VAUuS\npC4wZUoxkX377eGYY+C006CP+0jdlgFLkqQSZcJll8H668Ozz8L48XDooVVXpbKVnp0jYkXgHGAj\noAP4z8wcWfZ5JUmq2tixcOSRMG0aXH457LBD1RWpq3TFCtapwI2ZuT6wKTCxC84pSVJlXnkFDj8c\ndt8dBg+G++83XPU0pQasiFgBGJiZ5wNk5ozMnFzmOSVJqsqMGXDGGbDBBrD00jBxIhxyiINDe6Ky\ntwjXAl6JiPMpVq8eBI7KzLdKPq8kSV2qvR2++11YeWW44w7YaKOqK1KVyg5YfYAtgO9k5oMR8Xvg\nh8DP5j1wyJAh7z1ua2ujra2t5NIkSVpyzz4L3/8+jBwJJ58M++4LEVVXpUZqb2+nvb19kV4TmVlO\nNUBErArcm5kDaj/vCPxXZu49z3FZZh2SJDXa22/DSSfB738PRxwBxx0H/fpVXZW6QkSQmZ3G6FJX\nsDLzpYh4NiLWzczHgV2AR8o8pyRJZcqE668vbnOzxRbw4IOw5ppVV6VmU+oKFkBEbEoxpmEp4Eng\noMycNM8xrmBJkprexIlw1FHw/PPFoNBddqm6IlWhnhWs0gNWPQxYkqRmNmkS/PzncPHF8JOfFCMY\nllqq6qpUlXoClpPcJUlagI4OOO88WG89ePNNeOSR4kpBw5UWxrsgSZI0HyNHFlPY+/SBG26Arbaq\nuiK1ElewJEmaw4svwkEHFeMWjjwSRowwXGnRGbAkSQLefbeYY7XxxrDKKvDoo3DggdDL/1NqMbhF\nKEnq8YYPL64OHDAA7rkH1l236orU6gxYkqQe68kni3lWDz9cDAzdc8+qK1J34cKnJKnHmToVfvxj\n+NSnYNtt4aGHDFdqLAOWJKnHyITLLy/GLjz9NIwbBz/8IfTtW3Vl6m7q2iKMiJOB8zLz4ZLrkSSp\nFOPGFVcFTpkCl10GO+5YdUXqzupdwZoInB0RIyPisIhYscyiJElqlFdfhW9/G3bbDb72NXjgAcOV\nyldXwMrMczJzB2AwsCYwPiIujYidyyxOkqTFNXMmDB0KG2wAvXsX9xE89NDisVS2uq8ijIjewHq1\nr1eAccCxEfGtzDygpPokSVpkd91V3NLmAx+A224rZltJXamumz1HxCnAXsAdwLmZef8czz2WmZ9c\noiK82bMkqQGefRaOOw7+9jc46STYbz+ITm/JKy26Rt7seTywWWZ+a85wVfOpxapOkqQGefttOOEE\n2HzzYkjoxImw//6GK1Wn3oD1tcycOucvIuJ2gMyc1PCqJEmqQyZcfz1suCGMHl00sP/859CvX9WV\nqafrtAcrIpYB+gErR8QHgFl/F1gB+FjJtUmStECPPgpHHw3PPANnnQW77lp1RdJsC1vB+hYwiqKx\nfXTt8SjgeuD0ckuTJOn9Jk+G738fBg6E3Xcv5lsZrtRs6m1yPzIz/1BaETa5S5IWoqMDLroIjj8e\n9tgDTjwRVlml6qrUE9XT5L6wLcJ/z8w7gOcjYt95n8/Ma5ewRkmSFur++4sp7BFFz9XWW1ddkdS5\nhc3B2oliNMPe83kuAQOWJKk0L70EP/oR3Hwz/OpXxST2Xt5FVy2gri3C0otwi1CSNIfp0+EPf4D/\n/m846CD48Y9hhRWqrkoqNGwOVkRcPOf9ByPi47PGNEiS1Ei33gqbbAK33AIjRsBvfmO4Uuup91Y5\nI4CREXEsxXiGHwDfK60qSVKP8+ST8L3vwYQJcMopsNdeDgpV66p7izAidgT+l+I+hJtn5osNK8It\nQknqsaZOLfqrhg4tAtYxx8Ayy1RdlbRgjdwiPBA4DxgMXADcGBGbLnGFkqQeKxOuuALWX79YvRo7\ntmhoN1ypO6h3DtafgEMz8+Xaz58Czs7MzRpShCtYktSjjBsH3/0uTJpUNLMPHFh1RVL96lnBWuyr\nCCNi6cx8t47jngYmAR3A9Mx8382hDViS1DO89hr85Cdw9dXFPQMPOQR69666KmnRNHKLcN2IuD0i\nHqr9vAlwXJ11dABtmbn5/MKVJKn7mzkTzjyz2A6MgIkT4bDDDFfqvuod1zYM+BEwHSAzxwMH1Pna\nWITzSJK6mbvvhi23hMsvL0YwnH46fPCDVVcllaveMQ39MvP+mPt62Rl1vjaBWyNiJkXf1rBFKVCS\n1Jqeew6OO66YZXXSSbD//o5dUM9Rb8B6JSLWpghLRMR+wAt1vnaHzHwhIj5MEbQmZuaIeQ8aMmTI\ne4/b2tpoa2ur8+0lSc3k7bfhd78rvg4/HIYNg/79q65KWnzt7e20t7cv0mvqvYpwAHA2sD3wOvAU\n8LXMfHqRThbxM+DNzPzdPL+3yV2SWlwm/OUvxRyrjTeGk0+GAQOqrkpqvIZfRRgR/YFemflmncf3\nqx0/pfbaW4CfZ+Yt8xxnwJKkFvbYY3D00fD003DqqbDbblVXJJWnnoDV6RZh7dY4831jgHlXouZj\nVeC6iMjauS6ZN1xJklrX5Mnwy1/CBRfA8cfDEUfAUktVXZVUvYX1YC2/JG+emU8BDRlGKklqHh0d\ncPHFxeT13XeHhx6CVVetuiqpeSz2oNGGFuEWoSS1jAcegCOPLHqu/vAH+JQTDtXDNHLQ6ICIuCEi\n/hURL0fE9bXGd0lSD/Hyy3DwwbDPPsWQ0HvvNVxJC1LvANBLgSuBjwAfBa4CLiurKElS85g+HX7/\ne9hwQ1hppWIK+ze+Ab0cIS0tUL1jGsZn5ibz/G5cZm7akCLcIpSkpnTbbcVNmddYo7g6cL31qq5I\nql7DxjRExK8p5l9dTjFs9EvAB4DfAmTma0tYqAFLkprIU0/B974H48bBKafA3ns7hV2apZEB66lO\nns7MXKJ+LAOWJDWHadPgV7+CM86AY48tQtYyy1RdldRclngOVu1NelFMbb+nYZVJkppKJlx9NXz/\n+7D99jB2bLEtKGnx1LuCNSYzNy+tCFewJKkyEyYUfVavvVaMXRg0qOqKpObWsDENwO0R8YUId+Al\nqbt47bVintUuu8AXvwijRhmupEapN2B9i2I0w7sRMTki3oyIySXWJUkqycyZcNZZsP76xUT2iRPh\n8MOhz0KbRiTVq67/nDJziW6ZI0lqDvfcU6xaLbccDB8Om3kzM6kU9U5yj4j4WkT8pPbzGhHh/F5J\nahHPPw9f/SoccAAcdxzceafhSipTvVuEfwS2A75S+3kKcEYpFUmSGuadd4qxC5tuCmutBY8+WoQs\nO2qlctW7475NZm4REWMAMvP1iFi6xLokSUsgE/76Vzj6aNhoIxg5EtZeu+qqpJ6j3oA1PSJ6U0xx\nJyI+DHSUVpUkabE9/ngRrJ58shgY+pnPVF2R1PPUu0V4GnAdsEpEnACMAE4srSpJ0iKbPLnor9p+\n+2L0wvjxhiupKvVeRXhJRIwCdgEC+HxmTiy1MklSXTo64JJL4Ic/hN12g4cegtVWq7oqqWfrNGBF\nxDLAYcAngAnAWZk5oysKkyQt3IMPFlPYZ8yAa6+FbbapuiJJsPAtwguBrSjC1WeBk0qvSJK0UC+/\nDIccAnvvXXy/7z7DldRMFrZFuEFmbgwQEecC95dfkiRpQf7+dzjzTLjoIhg8uBi7sOKKVVclaV4L\nW8GaPuuBW4OSVI1Z23+f/jTssENxS5uRI+Hkkw1XUrOKzFzwkxEzgamzfgSWBabVHmdmrtCQIiKy\nszokqSd6/nkYNgzOOQfWXLO4X+B++0HfvlVXJvVsEUFmdjqut9Mtwszs3diSJEmd6eiA22+HoUOh\nvb2Yun7TTbDxxlVXJmlReO90SWoCr74KF1wAZ50Fyy5brFZdeCEsv3zVlUlaHAYsSapIZtFLNXQo\nXH99cUXgBRfAdtt5r0Cp1XXag9VlRdiDJakHmToVLr20CFaTJsFhh8FBB8HKK1ddmaR61NOD1SUB\nKyJ6AQ8Cz2Xm5+bzvAFLUrf38MPFiIVLL4WBA4ttwE9/GnrVe9MySU1hiZvcG+go4BGgIVcdSlKr\nePfdYsTC0KHFDKuDD4axY2GNNaquTFKZSg9YEbE6sAdwAnBs2eeTpGbw9NNw9tlw3nmw4YZw5JGw\nzz6w1FJVVyapK3TFwvQpwA8A9wAldWszZ8Jf/wp77QVbbQVvvQV33lmMXdhvP8OV1JOUuoIVEXsC\nL2Xm2IhooxhQOl9Dhgx573FbWxttbW1lliZJDfPSS3DuucWK1SqrFL1VV14J/fpVXZmkRmhvb6e9\nvX2RXlNqk3tEnAh8DZhBMQV+eeDazBw8z3E2uUtqKZlw111Fb9Xw4fCFLxTBasstq65MUtma5irC\nWjE7Ad/zKkJJrWzSJLj44iJYdXQUoWrwYFhppaork9RVmukqQklqaWPGFKHqqqtgt93gjDNgp50c\nCCpp/rosYGXmncCdXXU+SVpSb71V9FINHQovvACHHgoTJ8Jqq1VdmaRm5yR3SZrH3/9eDAS96CLY\neutiG3CPPaB376ork9QM3CKUpDrNmAF//nOxWjV+fHHrmpEjYcCAqiuT1IoMWJJ6tOefh2HDiq8B\nA4rVqi98Afr2rboySa3MgCWpx+noKIZ/Dh0K7e3w5S/DzTfDxhtXXZmk7sKAJanHePVVuOCCor+q\nf/9iterCC2H55auuTFJ3Y8CS1K1lFr1UQ4cWPVZ77100r2+7rSMWJJXHqwgldUtTpsCllxbB6s03\n4bDD4BvfgJVXrroySa2uqSa5d1qEAUtSgzz8cBGqLr20GAR6+OGw667QqytubS+pR3BMg6Qe4Z13\n4Npri2D1xBNw8MHFqIXVV6+6Mkk9lQFLUst6+mk46yw477ziCsCjjoLPfQ6WWqrqyiT1dC6aS2op\nM2fCX/4Ce+4JW21VrF7ddRfcdlsxv8pwJakZuIIlqSW89BKcey6cfTasumrRW3X11bDsslVXJknv\nZ8CS1LQyi9WpoUNh+HDYbz+45hrYcsuqK5OkznkVoaSmM2lSMavqzDOLkHX44XDggbDSSlVXJkle\nRSipxYweXaxWXX01fOYz8Mc/wqBBDgSV1HoMWJIq9dZbcMUVRbB68UX41rfg0UeLPitJalVuEUqq\nxOOPF1uAF10E22xTbAN+9rPQu3fVlUlS59wilNRUpk8v7gc4dChMmAD/+Z/wwAOw1lpVVyZJjWXA\nklS6556DYcPgnHNg7bWL1ap994W+fauuTJLKYcCSVIqOjmL459ChcOed8OUvF6MWNtqo6sokqXwG\nLEkN9eqrcP75xS1s+vcvVqsuvhiWW67qyiSp6xiwJC2xTLjvvmK16oYbivsBXnxx0bzuiAVJPZFX\nEUpabFOmwCWXFMFq6lQ47DD4xjfgQx+qujJJKk89VxEasCQtsoceKkLVZZdBW1uxDbjLLtDL28dL\n6gEc0yCpYd55p7gP4NCh8OSTcMghMH48rL561ZVJUvMpNWBFRF/gLmDp2rmuzsyfl3lOSY311FNF\nw/r558Mmm8Axx8Dee8NSS1VdmSQ1r1IDVma+ExE7Z+a0iOgN3BMRN2Xm/WWeV9KSmTkTbryxWK26\n/34YPBjuvhvWXbfqyiSpNZS+RZiZ02oP+9bOZ7OV1KRefBHOPRfOPhtWW63orbrmGlh22aork6TW\nUnpLakT0iogxwIvArZn5QNnnlFS/TGhvhy99CdZfH55+Gq67DkaOLK4INFxJ0qLrihWsDmDziFgB\n+FNEbJCZj5R9Xkmde+ON4kbLZ55Z/Hz44cXK1YorVluXJHUHXXYVYWZOjoj/BXYH3hewhgwZ8t7j\ntrY22trauqo0qUfILMYr3HRT8fXgg7DnnkWf1aBBDgSVpAVpb2+nvb19kV5T6hysiFgZmJ6ZkyJi\nWWA48KvMvHGe45yDJZVg8uTifoA33QQ331xc+ffZzxZfO+9c3MpGkrRoKh80GhEbAxdS9Hr1Aq7I\nzBPmc5wBS2qAWatUN95YhKpRo2D77WeHqnXXdaVKkpZU5QGrXgYsafHNuUp1002w9NKuUklSmQxY\nUjeUCRMmzA5UrlJJUtcyYEndxKRJc/dS9e07O1C1tblKJUldyYAltaj5rVLtsMPsULXOOq5SSVJV\nDFhSC+lslWrnnaFfv6orlCSBAUtqapkwfvzsVarRo12lkqRWYMCSmsykSXDrrbNXqZZddu5eKlep\nJKn5GbCkis1apZo1l2rMGNhxx7lXqSRJrcWAJVXgjTfm7qVylUqSuhcDltQFMmHcuNm9VK5SSVL3\nZsCSSvLGG3P3UvXvP/cq1bLLVl2hJKksBiypQeZcpbrxxuLxnKtUn/hE1RVKkrqKAUtaAgtapdpj\nD9hpJ1epJKmnMmBJiyATxo6d3UvlKpUkaX4MWNJCvP763KtUyy8/O1C5SiVJmh8DljSPWatUs+ZS\njRsHgwbNDlVrr111hZKkZmfAknCVSpLUWAYs9UgdHXP3Uo0fDwMHukolSWoMA5Z6jNdfh1tumb1K\nteKKswPVoEGuUkmSGseApW6ro6OYmD5rlWrChLl7qQYMqLpCSVJ3ZcBSt/Laa0Uv1Y03wvDhs1ep\n9tijCFfLLFN1hZKknsCApZbmKpUkqRkZsNRyXnttdi/V8OGw0kpz91K5SiVJqpoBS01v1irVrLlU\nDz1UjE6YFarWWqvqCiVJmpsBS01pzlWqm2+GD35wdqAaONBVKklSczNgqSl0dMDo0bN7qVylkiS1\nMgOWKvPqq3P3UrlKJUnqLioPWBGxOnARsCrQAQzLzNPmc5wBq8V1dMCoUbNXqR5+GNraikC1++6u\nUkmSuo9mCFirAatl5tiIWA4YBeyTmY/Oc5wBq8VkwnPPwYgR71+l2mOPYpWqb9+qq5QkqfEqD1jv\nO1nEn4A/ZObt8/zegNXEMuGZZ4o+qlGjZn9FwLbbFitUn/0srLlm1ZVKklS+pgpYEbEm0A5slJlT\n5nnOgNUkMuH//m/uIDV6NPTpA1tuCVtsUXzfckv42MeKkCVJUk/SNAGrtj3YDvwyM6+fz/MGrApk\nwlNPvT9MLbPM3EFqyy3hox+tulpJkppDPQGrTxcU0Qe4Grh4fuFqliFDhrz3uK2tjba2trJL61Ey\n4Ykn3h+mlltudpg6+uji+2qrVV2tJEnNo729nfb29kV6TekrWBFxEfBKZh7byTGuYDVQRwf84x9z\nh6kxY2CFFeZeldpiC1h11aqrlSSptVS+RRgROwB3AROArH0dn5k3z3OcAWsxdXTA44/PvSo1Zgx8\n4APvD1NnDiBqAAAIl0lEQVQf/nDV1UqS1PoqD1j1MmDVZ+ZMeOyxucPU2LGw8srvD1Mf+lDV1UqS\n1D0ZsFrYjBnw6KNzh6lx44otvTnD1OabF/OnJElS1zBgtYgZM+CRR2YHqVGjYPz44sq9ecPUSitV\nXa0kST2bAasJTZ9e3EZmzqGdEybAGmvMvcW3+eaw4opVVytJkuZlwKrYu+/CQw/NHaYefhg+/vG5\nh3Zuvjksv3zV1UqSpHoYsLrQO+8UYWrO0QiPPAIDBswdpjbbrJg9JUmSWpMBqyRvv11s680Zph59\nFD7xibknoG+6KfTvX3W1kiSpkQxYDfDWW0XD+Zxh6vHHYd115w5Tm2wC/fpVXa0kSSqbAWsRTZtW\njEKYM0z94x+w3nrvD1PLLFN1tZIkqQoGrE5MnVoM6ZxzztQTT8D66889GmHjjaFv3y4tTZIkNTED\nVs2UKcXtY+YMU089BRtuOHeY2nBDw5QkSepcjwxYkyfPDlOzxiM88wxstNHcYWqDDWDppRtySkmS\n1IN0+4A1adLsEDXr+3PPFT1Scw7t3GADWGqpEgqXJEk9TrcKWK+/XoSoOYd2vvBCMQphzjC1/vrQ\np08XFS5Jknqclg1Yr70296rUqFHw8svFkM45h3autx707l1h4ZIkqcdpqYB14on5Xph69dXi9jFz\njkZYd13DlCRJql5LBaxjj833wtQ660CvXlVXJUmS9H4tFbCaoQ5JkqSFqSdguU4kSZLUYAYsSZKk\nBjNgSZIkNZgBS5IkqcEMWJIkSQ1mwJIkSWowA5YkSVKDGbAkSZIazIAlSZLUYKUGrIg4NyJeiojx\nZZ5HkiSpmZS9gnU+8JmSz6EKtbe3V12CloCfX+vys2ttfn7dX6kBKzNHAK+XeQ5Vyz8kWpufX+vy\ns2ttfn7dnz1YkiRJDWbAkiRJarDIzHJPEPFx4IbM3KSTY8otQpIkqYEyMzp7vk8X1BC1rwVaWJGS\nJEmtpOwxDZcCfwPWjYhnIuKgMs8nSZLUDErfIpQkSeppKm1ydxBp64qI1SPijoh4OCImRMR3q65J\n9YmIvhExMiLG1D67n1VdkxZdRPSKiNER8eeqa9GiiYinI2Jc7b/B+6uuR/WLiBUj4qqImFj7/982\nCzy2yhWsiNgRmAJc1FkTvJpPRKwGrJaZYyNiOWAUsE9mPlpxaapDRPTLzGkR0Ru4B/huZvoHfQuJ\niGOALYEVMvNzVdej+kXEk8CWmemcyBYTERcAd2bm+RHRB+iXmZPnd2ylK1gOIm1dmfliZo6tPZ4C\nTAQ+Vm1VqldmTqs97EtxsYu9Ai0kIlYH9gDOqboWLZbAMUktJyJWAAZm5vkAmTljQeEK/IDVABGx\nJrAZMLLaSlSv2vbSGOBF4NbMfKDqmrRITgF+gMG4VSVwa0Q8EBGHVF2M6rYW8EpEnF/bnj87IpZd\n0MEGLC2R2vbg1cBRtZUstYDM7MjMzYHVgW0iYoOqa1J9ImJP4KXaCvJCx+CoKe2QmVtQrEJ+p9Yu\no+bXB9gCOKP2+U0Dfriggw1YWmy1/eergYsz8/qq69Giqy1v/y+we9W1qG47AJ+r9fFcBuwcERdV\nXJMWQWa+UPv+L+A64FPVVqQ6PQc8m5kP1n6+miJwzVczBCz/Bta6zgMeycxTqy5E9YuIlSNixdrj\nZYFPA16c0CIy8/jM/LfMHAAcANyRmYOrrkv1iYh+tZV/IqI/sBvwULVVqR6Z+RLwbESsW/vVLsAj\nCzq+Kya5L1BtEGkb8KGIeAb42azmMTW3iNgB+CowodbLk8DxmXlztZWpDh8BLoyIXhR/yboiM2+s\nuCapp1gVuK52i7g+wCWZeUvFNal+3wUuiYilgCeBBQ5Qd9CoJElSgzXDFqEkSVK3YsCSJElqMAOW\nJElSgxmwJEmSGsyAJUmS1GAGLEmSpAYzYEnqEhHRERG/nePn70XET6usaWEiYsuI+H3VdUhqPQYs\nSV3lHWDfiPhgI94sIno34n06k5mjMvPoss8jqfsxYEnqKjOAs4FjF3ZgRHwzIh6LiPtqd6w/rfb7\n8yNiaETcB/w6IraOiL9FxKiIGBER69SO+3pEXBcRt0TEkxHxnYg4JiJG145faT7n3D8iJkTEmIho\nr/1up4i4ofb4r7XXj4mINyLiwIjoFRG/iYiRETE2Ig5p3L8uSa2s0lvlSOpREjiD4vZKv17QQRHx\nEeDHwGbAFIqbUY+d45CPZea2tWOXA3bMzI6I2AX4b2C/2nEb1t6jH/AP4AeZuUVE/A4YDJw2z6l/\nAuyWmS9ExArz1E1m7lk75xYU9+H8E/BN4I3M3CYilgbuiYhbMvP/FuVfjKTux4Alqctk5pSIuBA4\nCnhrAYd9CmjPzEkAEXEVsM4cz181x+OVgItqK1ez7u02y/9m5jRgWkS8Afyl9vsJwMbzOe8Iins0\nXglcO7/CImJl4GJgv8x8MyJ2AzaOiP1rh6xQq9WAJfVwBixJXe1UYDRwPkDtptOjKALSn4ExQHTy\n+qlzPP4lcEdm7hsRH6dY7ZrlnTke5xw/dzCfP/sy89sRsTWwFzCqtlL1nlqdlwFDMnPirF8DR2bm\nrZ3UK6kHsgdLUlcJgMx8HbiSYnuNzOzIzM0zc4vMHAI8AAyKiBUjog/whU7ecwXg+drjBd7Vvq7i\nIgZk5gOZ+TPgZWCNeQ75NTAuM+dcQRsOfLtWJxGxTkQsuyR1SOoeDFiSukrO8fhk4EPz/K44KPOf\nwInA/cDdwFPApPm8B8BvgV9FxCg6//PsfeeZj99GxPiIGA/ck5nj53n+e8ButSb30RGxV2YOAx4B\nRkfEBOBM3BmQBERmPX/uSFLXiYj+mTm1NorhOuDczLy+6rokqV6uYElqRkMiYgxFQ/qThitJrcYV\nLEmSpAZzBUuSJKnBDFiSJEkNZsCSJElqMAOWJElSgxmwJEmSGsyAJUmS1GD/H++b38QGNzYzAAAA\nAElFTkSuQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x898b7b8>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Plot perplexity values\n",
    "plt.figure(figsize=(10, 5))\n",
    "ns = range(1, max_ngram + 1)\n",
    "plt.plot(ns, [n_to_perplexity[n] for n in ns])\n",
    "plt.xlabel('N-gram size')\n",
    "plt.ylabel('Perplexity')"
   ]
  }
 ],
 "metadata": {
  "anaconda-cloud": {},
  "kernelspec": {
   "display_name": "Python [Root]",
   "language": "python",
   "name": "Python [Root]"
  },
  "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
}