hw-week7-fake-genome.ipynb

{
 "metadata": {
  "name": "hw-week6-fake-genome"
 },
 "nbformat": 3,
 "nbformat_minor": 0,
 "worksheets": [
  {
   "cells": [
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "## Making a fake genome\n",
      "\n",
      "For class, we've been working with a bunch of reads from the E. coli genome for a week or so.  But the E. coli genome is annoyingly large to play with, and you need lots of reads to get good coverage for assembly, so that's even more annoying.\n",
      "\n",
      "Let's make a simple simulated genome, at random, to play with."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "import random\n",
      "import string\n",
      "\n",
      "SEED=0                 # set SEED to a number between 1 and 9.\n",
      "random.seed(SEED)\n",
      "\n",
      "genome = [\"A\"] * 2500 + [\"C\"] * 2500 + [\"G\"] * 2500 + [\"T\"] * 2500\n",
      "\n",
      "# 'genome' is now a list with 10,000 As, Cs, Gs, and Ts in it.  Randomize, by shuffling the list:\n",
      "\n",
      "random.shuffle(genome)\n",
      "\n",
      "# Now that we've shuffled the order of nucleotides in 'genome', make it into a string:\n",
      "\n",
      "genome = string.join(genome, \"\")\n"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 1
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# genome is now a 10k long string.  let's save it to a file.\n",
      "outfp = open('/mnt/fake-genome.fa', 'w')\n",
      "outfp.write(genome)\n",
      "outfp.close()"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 2
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# look at that file using the shell command 'ls': exactly 10k characters!\n",
      "!ls -l /mnt/fake-genome.fa"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "-rw-r--r-- 1 root root 10000 2012-10-16 01:34 /mnt/fake-genome.fa\r\n"
       ]
      }
     ],
     "prompt_number": 3
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "## Sequencing a fake genome\n",
      "\n",
      "OK, so we've got a fake genome.  Let's make some fake reads, of length 100, to a coverage of 100x."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# first, make a function that returns a random read from the genome (we'll skip the reverse complement stuff for now)\n",
      "def make_read(genome):\n",
      "    READLEN=100\n",
      "    \n",
      "    # pick a random starting position between 0 and the length of the genome, minus the read length\n",
      "    pos = random.randint(0, len(genome) - READLEN) \n",
      "    \n",
      "    # extract the read sequence\n",
      "    read = genome[pos:pos + READLEN]\n",
      "    \n",
      "    # return!\n",
      "    return read\n"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 4
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# ok, next, we want to put errors in the reads.  Let's just hardcode the error rate to be 1 in 100.\n",
      "# note -- an error rate of 1 in 100 means we *expect* approximately 1 error every 100 bp, but it doesn't guarantee 1 error/read!\n",
      "\n",
      "# the error rate is set by the random.randint function call, below -- random.randint returns a number between\n",
      "# the first and second parameter, inclusive.  For a 10% error rate you could use 1, 10.\n",
      "\n",
      "def make_substitution_error(read):\n",
      "    for i in range(len(read)):\n",
      "        if random.randint(1, 100) == 1:\n",
      "            new_nt = random.choice(\"acgt\")    # make it lowercase for easy identification\n",
      "            read = read[:i] + new_nt + read[i + 1:]\n",
      "    return read\n",
      "\n"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 5
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# now that we've got these functions defined, let's try 'em out:\n",
      "r = make_read(genome)\n",
      "r2 = make_substitution_error(r)\n",
      "print r2"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "TAATTACCAGAAGAATCTATATATTCAGTGCTGGGGCTGCTTAGTCATCAACGTTGCAATGGATCACAGGATgGCTGAAAACGCGGGTTTTCACCAACTT\n"
       ]
      }
     ],
     "prompt_number": 6
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "# some quick calculations: genome size is 10k, 100x coverage would be 100 x 10k bases, at 100 bp per read it would be 10k reads total.\n",
      "# so, make a fake read sequencing data set containing 10k reads.\n",
      "\n",
      "# here, coverage is set by the number of reads (10,000)\n",
      "\n",
      "outfp2 = open('/mnt/fake-reads.fa', 'w')\n",
      "for i in range(10000):\n",
      "    read = make_read(genome)\n",
      "    read = make_substitution_error(read)\n",
      "    outfp2.write('>read%d\\n%s\\n' % (i, read))\n",
      "outfp2.close()"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 7
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "## OK, now let's look at the k-mer abundance distribution in the read data set.\n",
      "\n",
      "Note, for this, you'll need to have installed khmer and screed as per the homework instructions.\n",
      "\n",
      "First, go to the /mnt working directory."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "cd /mnt"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "/mnt\n"
       ]
      }
     ],
     "prompt_number": 8
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Next, ask khmer to count the k-mers contained by the reads in 'fake-reads.fa' and save the resulting data into a file 'fake-reads.kh'.\n",
      "\n",
      "Note that the '-k 26' parameter specifies the k-mer size to use -- 26.  You can change this to any number between 1 and 32."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "!/usr/local/share/khmer/scripts/load-into-counting.py -x 1e8 -k 26 fake-reads.kh fake-reads.fa"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "\r\n",
        "PARAMETERS:\r\n",
        " - kmer size =    26 \t\t(-k)\r\n",
        " - n hashes =     4 \t\t(-N)\r\n",
        " - min hashsize = 1e+08 \t(-x)\r\n",
        "\r\n",
        "Estimated memory usage is 4e+08 bytes (n_hashes x min_hashsize)\r\n",
        "--------\r\n",
        "Saving hashtable to fake-reads.kh\r\n",
        "Loading kmers from sequences in ['fake-reads.fa']\r\n",
        "making hashtable\r\n"
       ]
      },
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "consuming input fake-reads.fa\r\n"
       ]
      },
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "saving fake-reads.kh\r\n"
       ]
      },
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "fp rate estimated to be 0.000\r\n",
        "DONE.\r\n"
       ]
      }
     ],
     "prompt_number": 9
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "ls"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "fake-genome.fa  fake-reads.hist  fake-reads.kh.info\r\n",
        "fake-reads.fa   fake-reads.kh    \u001b[0m\u001b[01;34mlost+found\u001b[0m/\r\n"
       ]
      }
     ],
     "prompt_number": 10
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Now, use khmer to get the abundance distribution of the k-mers in fake-reads.fa using the counts from fake-reads.kh, and put the results in file called fake-reads.hist."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "!/usr/local/share/khmer/scripts/abundance-dist.py fake-reads.kh fake-reads.fa fake-reads.hist\n"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "hashtable from fake-reads.kh\r\n"
       ]
      },
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "K: 26\r\n",
        "HT sizes: [100000007, 100000037, 100000039, 100000049]\r\n",
        "outputting to fake-reads.hist\r\n",
        "ERROR: fake-reads.hist exists; not squashing.\r\n"
       ]
      }
     ],
     "prompt_number": 11
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "fake-reads.hist is a file that contains four columns.  The columns are:\n",
      "\n",
      "1. a k-mer count\n",
      "2. the number of k-mers in the data set with that count.  For example, below, 7,793 k-mers appear exactly twice in the read data set, while 309 appear three times.\n",
      "3. the total number of k-mers at or below this k-mer count.  The third column is the sum of the second column.\n",
      "4. the fraction of the total k-mers seen so far -- this is the third column divided by the sum of the entire second column."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "!cat fake-reads.hist"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "0 0 0 0.0\r\n",
        "1 130315 130315 0.878\r\n",
        "2 7793 138108 0.931\r\n",
        "3 309 138417 0.933\r\n",
        "4 21 138438 0.933\r\n",
        "5 2 138440 0.933\r\n",
        "6 3 138443 0.933\r\n",
        "7 1 138444 0.933\r\n",
        "8 4 138448 0.933\r\n",
        "9 3 138451 0.933\r\n",
        "10 1 138452 0.933\r\n",
        "11 3 138455 0.933\r\n",
        "12 1 138456 0.933\r\n",
        "13 2 138458 0.933\r\n",
        "14 1 138459 0.933\r\n",
        "15 4 138463 0.933\r\n",
        "16 3 138466 0.933\r\n",
        "17 4 138470 0.933\r\n",
        "18 2 138472 0.933\r\n",
        "19 2 138474 0.933\r\n",
        "20 5 138479 0.933\r\n",
        "21 1 138480 0.933\r\n",
        "22 4 138484 0.933\r\n",
        "23 6 138490 0.933\r\n",
        "24 2 138492 0.933\r\n",
        "25 1 138493 0.933\r\n",
        "26 3 138496 0.933\r\n",
        "27 2 138498 0.933\r\n",
        "28 3 138501 0.933\r\n",
        "29 3 138504 0.933\r\n",
        "30 4 138508 0.933\r\n",
        "31 1 138509 0.933\r\n",
        "32 7 138516 0.933\r\n",
        "33 5 138521 0.933\r\n",
        "34 7 138528 0.934\r\n",
        "35 8 138536 0.934\r\n",
        "36 20 138556 0.934\r\n",
        "37 13 138569 0.934\r\n",
        "38 25 138594 0.934\r\n",
        "39 61 138655 0.934\r\n",
        "40 89 138744 0.935\r\n",
        "41 109 138853 0.936\r\n",
        "42 109 138962 0.936\r\n",
        "43 122 139084 0.937\r\n",
        "44 219 139303 0.939\r\n",
        "45 246 139549 0.94\r\n",
        "46 314 139863 0.942\r\n",
        "47 356 140219 0.945\r\n",
        "48 357 140576 0.947\r\n",
        "49 413 140989 0.95\r\n",
        "50 452 141441 0.953\r\n",
        "51 443 141884 0.956\r\n",
        "52 530 142414 0.96\r\n",
        "53 520 142934 0.963\r\n",
        "54 495 143429 0.967\r\n",
        "55 489 143918 0.97\r\n",
        "56 428 144346 0.973\r\n",
        "57 470 144816 0.976\r\n",
        "58 451 145267 0.979\r\n",
        "59 451 145718 0.982\r\n",
        "60 365 146083 0.984\r\n",
        "61 375 146458 0.987\r\n",
        "62 335 146793 0.989\r\n",
        "63 271 147064 0.991\r\n",
        "64 288 147352 0.993\r\n",
        "65 233 147585 0.995\r\n",
        "66 209 147794 0.996\r\n",
        "67 155 147949 0.997\r\n",
        "68 121 148070 0.998\r\n",
        "69 86 148156 0.998\r\n",
        "70 74 148230 0.999\r\n",
        "71 50 148280 0.999\r\n",
        "72 40 148320 0.999\r\n",
        "73 33 148353 1.0\r\n",
        "74 21 148374 1.0\r\n",
        "75 14 148388 1.0\r\n",
        "76 8 148396 1.0\r\n"
       ]
      }
     ],
     "prompt_number": 12
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Now, read in the data from fake-reads.hist and plot it."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "import csv\n",
      "r = csv.reader(open('fake-reads.hist'), delimiter=' ')\n",
      "\n",
      "xx = []\n",
      "yy = []\n",
      "for count, abund, sum_abund, sum_frac in r:\n",
      "    xx.append(int(count))\n",
      "    yy.append(int(abund))"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 13
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "You can see that the plot, below, has an awful lot of low-abundance k-mers..."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "plot(xx, yy)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "pyout",
       "prompt_number": 14,
       "text": [
        "[<matplotlib.lines.Line2D at 0x30f3b90>]"
       ]
      },
      {
       "output_type": "display_data",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAYsAAAD9CAYAAABN7FvjAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3X9UVPed//Hn4A+MG9aqKSTfg4Mh0oHxJ8RhON1aJp7U\ncnZrUFOPuoubFXLMunuO69bTo4meRrse2/zoMbEnYswJ3W1tStUctvklVN0O0bbOwDZbk2E00DVH\nk40Yyi4CghL5fP+YzBWiYsBLHG9fj3PmjPOZ+Xx4X8T78n3vncFljDGIiIgMIOlmFyAiIolPYSEi\nItelsBARketSWIiIyHUpLERE5LoUFiIicl0DhkVpaSlpaWlMnz79iud+8IMfkJSURGtrqzW2fft2\nsrKy8Hq9HDlyxBqPRqPk5eWRmZnJhg0brPGenh7KysrIyMggEAhw5swZ67m9e/fi8XjweDzs27fv\nhjZSRERuzIBhsWLFCqqrq68YP336NAcOHCAjI8MaO3v2LDt27ODQoUOUl5ezevVq67m1a9eybt06\n6urqqK2tpb6+HoCqqira2tqIRqMUFRWxZcsWAHp7e1m/fj0vv/wye/fuZf369bZsrIiIDM2AYTFn\nzhzGjx9/xfi3vvUtnnzyyX5joVCIoqIi3G43hYWFGGPo6OgA4MSJEyxZsoSJEyeyaNEiQqGQNaek\npISxY8eycuVKazwSiTBt2jSmTZvGjBkz8Hq9RCIRWzZYREQGb9DnLH7xi1+Qnp7OjBkz+o2Hw2Fy\ncnKsxx6Ph1AoRFNTE6mpqda41+vl6NGj1hyv1wvAhAkTaG5upru7m1AoZI1/eo6IiHz+Rg7mxefP\nn2fr1q0cOHDAGot/WsjVPjXE5XJdMWaMscaNMf3mDfTJI1dba6BxERG5tsF+0tOgOos//OEPvPfe\ne8ycOZO7776b999/n3vvvZfm5mb8fj8NDQ3Wa48fP47P52PKlCk0Nzdb4w0NDfj9foB+c1pbW0lL\nS2PMmDFXrNV3ztXEQydRb48//vhNr0F1qk7VqTrjt6EYVFhMnz6d5uZmTp48ycmTJ0lPT+d3v/sd\naWlp5OfnU1NTw6lTpwgGgyQlJZGSkgJAdnY2lZWVtLS0UFVV1S8sdu/eTWdnJ7t27aKgoACIHXZ6\n5513ePvttzl27BiRSISpU6cOaQNFROTGDXgYatmyZdTW1vLHP/6RSZMm8d3vfpcVK1ZYz/c9BJSW\nlsaqVauYO3cuo0eP5vnnn7eee/rppykpKeHRRx9l6dKlzJ49G4CFCxdSXV1NTk4OmZmZVFZWAjBi\nxAi+973v8eCDDwLw/e9/374tFhGRQXOZofYkCcLlcg25rfq8BINBAoHAzS7julSnvVSnvVSnfYay\n31RYiIj8iRnKflMf9yEiItelsBARketSWIiIyHUpLERE5LoUFiIicl0KCxERuS6FhYiIXJfCQkRE\nrstRYfHmm9DefrOrEBFxHkeFxXe+A7/97c2uQkTEeRwVFj09sZuIiNjLUWFx8aLCQkRkODgqLHp6\nYoEhIiL2clxYqLMQEbGf48JCnYWIiP0cFxbqLERE7Oe4sFBnISJiP0eFha6GEhEZHo4KC3UWIiLD\nw3Fhoc5CRMR+jgsLdRYiIvZzXFiosxARsd+AYVFaWkpaWhrTp0+3xr797W+Tk5NDXl4ea9asoaur\ny3pu+/btZGVl4fV6OXLkiDUejUbJy8sjMzOTDRs2WOM9PT2UlZWRkZFBIBDgzJkz1nN79+7F4/Hg\n8XjYt2/fdTfk0iUwRp2FiMhwGDAsVqxYQXV1db+xefPmEYlEqK+vp7Ozk5deegmAs2fPsmPHDg4d\nOkR5eTmrV6+25qxdu5Z169ZRV1dHbW0t9fX1AFRVVdHW1kY0GqWoqIgtW7YA0Nvby/r163n55ZfZ\nu3cv69evv+6GxENCnYWIiP0GDIs5c+Ywfvz4fmNf+9rXSEpKIikpia9//evU1tYCEAqFKCoqwu12\nU1hYiDGGjo4OAE6cOMGSJUuYOHEiixYtIhQKWXNKSkoYO3YsK1eutMYjkQjTpk1j2rRpzJgxA6/X\nSyQSGXBD4iGhzkJExH4jb2TyCy+8wMMPPwxAOBwmJyfHes7j8RAKhcjIyCA1NdUa93q9/PSnP+Uf\n//EfCYfDPPLIIwBMmDCB5uZmuru7CYVCeL3efnOOHj3K1KlTr1rHpk2biB8NO3UqAARuZLNERBwl\nGAwSDAZvaI0hh8V3v/tdUlJSWLx4MQDGmCte43K5rhgzxljjxph+8662xkBrxW3atIkzZ+DJJ+GO\nOz7zJoiI/EkIBAIEAgHr8ebNmwe9xpCuhvrXf/1Xampq2L17tzXm9/tpaGiwHh8/fhyfz8eUKVNo\nbm62xhsaGvD7/VfMaW1tJS0tjTFjxlyxVt8516LDUCIiw2fQYVFdXc1TTz3FK6+8wpgxY6zx/Px8\nampqOHXqFMFgkKSkJFJSUgDIzs6msrKSlpYWqqqq+oXF7t276ezsZNeuXRQUFACxw07vvPMOb7/9\nNseOHSMSiVzzEFRcPCx0gltEZBiYASxdutTcddddZtSoUSY9Pd28+OKLZsqUKcbtdptZs2aZWbNm\nmVWrVlmvf+aZZ8w999xjcnJyzJtvvmmNRyIRk5ubayZPnmzWr19vjV+8eNGsWLHCTJo0yRQWFpoP\nP/zQeu7nP/+5ycrKMllZWWbPnj3XrDG+CdGoMWDMN74x0BaJiMh1dv1X5fpk4i3L5XJhjOHtt2HG\nDPj61+FTV/uKiEgf8f3mYDjmHdw6ZyEiMnwcFxY6ZyEiYj9HhcWIEeosRESGg6PCYuxYdRYiIsPB\nUWHxZ3+mzkJEZDg4JiwuXoyFhToLERH7OSYs1FmIiAwfR4WFzlmIiAwPR4WFOgsRkeHhqLBQZyEi\nMjwcFRbqLEREhodjwkJXQ4mIDB/HhEXfw1C39kcjiogkHkeFxejRMHIkfPzxza5GRMRZHBUWo0bF\nbjpvISJiL8eFxejROm8hImI3x4WFOgsREfs5JiwuXox1FeosRETs55iwUGchIjJ8HBcW6ixEROzn\nuLAYNUphISJiN8eFxejROgwlImK3AcOitLSUtLQ0pk+fbo21t7dTXFyM2+1mwYIFdHR0WM9t376d\nrKwsvF4vR44cscaj0Sh5eXlkZmayYcMGa7ynp4eysjIyMjIIBAKcOXPGem7v3r14PB48Hg/79u27\n7oaosxARGT4DhsWKFSuorq7uN1ZeXo7b7aaxsZH09HR27twJwNmzZ9mxYweHDh2ivLyc1atXW3PW\nrl3LunXrqKuro7a2lvr6egCqqqpoa2sjGo1SVFTEli1bAOjt7WX9+vW8/PLL7N27l/Xr1193Q/pe\nDaXOQkTEXgOGxZw5cxg/fny/sXA4TFlZGcnJyZSWlhIKhQAIhUIUFRXhdrspLCzEGGN1HSdOnGDJ\nkiVMnDiRRYsW9ZtTUlLC2LFjWblypTUeiUSYNm0a06ZNY8aMGXi9XiKRyIAbos5CRGT4DPqcRV1d\nHdnZ2QBkZ2cTDoeB2I4/JyfHep3H4yEUCtHU1ERqaqo17vV6OXr0KBALHq/XC8CECRNobm6mu7ub\nUChkjX96zrXonIWIyPAZOdgJZhAf6epyua46Pz5ujOm33kBrX22tuE2bNhGNwp49cO5cgJ6ewGeu\nUUTE6YLBIMFg8IbWGHRY+Hw+otEoubm5RKNRfD4fAH6/n4MHD1qvO378OD6fj5SUFJqbm63xhoYG\n/H6/NaehoQGPx0NraytpaWmMGTMGv9/Pq6++2m/O8uXLr1nTpk2bOHwY/u7v4Lnn1FmIiPQVCAQI\nBALW482bNw96jUEfhvL7/VRUVNDV1UVFRQUFBQUA5OfnU1NTw6lTpwgGgyQlJZGSkgLEDldVVlbS\n0tJCVVVVv7DYvXs3nZ2d7Nq1y1rL6/Xyzjvv8Pbbb3Ps2DEikQhTp04dsC69KU9EZBiZASxdutTc\nddddZvTo0SY9Pd1UVFSYc+fOmQceeMBMmjTJFBcXm/b2duv1zzzzjLnnnntMTk6OefPNN63xSCRi\ncnNzzeTJk8369eut8YsXL5oVK1aYSZMmmcLCQvPhhx9az/385z83WVlZJisry+zZs+eaNcY3oaDA\nmF//2pi/+RtjfvzjgbZKRORP23V2/Vfl+mTiLcvlcmGM4d574fnnYccO+MpXoLT0ZlcmIpKY4vvN\nwXDcO7j1QYIiIvZzXFjonIWIiP0cFxbqLERE7Oe4sFBnISJiP8eFhToLERH7OSYs9GtVRUSGj2PC\nQp2FiMjwcVxYqLMQEbGf48JCnYWIiP0cERbGwMcfq7MQERkujgiLjz+GESPA5VJnISIyHBwRFvEr\noUCdhYjIcHBEWMTPV4A6CxGR4eC4sFBnISJiP8eFxahRCgsREbs5LixGj9ZhKBERuzkuLNRZiIjY\nzzFh0fdqKHUWIiL2ckRYXLyozkJEZDg5Iix0zkJEZHg5LizUWYiI2M9xYaHOQkTEfo4LC3UWIiL2\nG3JYvPDCC3z5y1/m3nvvZc2aNQC0t7dTXFyM2+1mwYIFdHR0WK/fvn07WVlZeL1ejhw5Yo1Ho1Hy\n8vLIzMxkw4YN1nhPTw9lZWVkZGQQCAQ4c+bMNWvRx32IiAyvIYVFa2srW7du5cCBA9TV1fHuu+9S\nU1NDeXk5brebxsZG0tPT2blzJwBnz55lx44dHDp0iPLyclavXm2ttXbtWtatW0ddXR21tbXU19cD\nUFVVRVtbG9FolKKiIrZs2XLNevRBgiIiw2tIYXHbbbdhjKGtrY2uri7Onz/PF77wBcLhMGVlZSQn\nJ1NaWkooFAIgFApRVFSE2+2msLAQY4zVdZw4cYIlS5YwceJEFi1a1G9OSUkJY8eOZeXKldb41aiz\nEBEZXkMOi/LyciZPnsydd97JX/zFX+D3+6mrqyM7OxuA7OxswuEwENvx5+TkWPM9Hg+hUIimpiZS\nU1Otca/Xy9GjRwEIh8N4vV4AJkyYQHNzMxcuXLhqPfogQRGR4TVyKJM++ugjVq1aRUNDA+PHj2fx\n4sW89tprGGM+8xoul+uKMWOMNW6M6bfeQGvv2bOJpibYtAm+8pUAFy8GPnMdIiJOFwwGCQaDN7TG\nkMIiHA5TUFDAlClTAFi8eDGHDx/G5/MRjUbJzc0lGo3i8/kA8Pv9HDx40Jp//PhxfD4fKSkpNDc3\nW+MNDQ34/X5rTkNDAx6Ph9bWVtLS0khOTr5qPd/4xiYOH46FhTFw6RL09kKSI671EhG5MYFAgEAg\nYD3evHnzoNcY0u50zpw51NfX09rayoULF9i/fz/z5s3D7/dTUVFBV1cXFRUVFBQUAJCfn09NTQ2n\nTp0iGAySlJRESkoKEDtcVVlZSUtLC1VVVf3CYvfu3XR2drJr1y5rravpexgq/qtVdShKRMQ+Q+os\n/vzP/5yNGzeycOFCzp8/T1FREffddx/5+fmUlJTg8XjIy8vjiSeeACAtLY1Vq1Yxd+5cRo8ezfPP\nP2+t9fTTT1NSUsKjjz7K0qVLmT17NgALFy6kurqanJwcMjMzqaysvGY9fa+GgsvnLa7RiIiIyCC5\nzGBONCQgl8vFD35gOH0atm2LjY0fD3/4A0yYcHNrExFJRC6Xa1DnmMGB7+AGXRElImI3R4aF3msh\nImIvR4aFOgsREXs5MizUWYiI2MsRYXGtq6FERMQejggLdRYiIsPLkWGhzkJExF6ODAu9g1tExF6O\nDAv9alUREXs5MizUWYiI2MsxYfHpq6HUWYiI2McRYXHxojoLEZHh5Iiw0DkLEZHh5ciwUGchImIv\nR4aFOgsREXs5MizUWYiI2MuRYaHOQkTEXo4Ii09/kKA6CxERezkiLNRZiIgML0eGhToLERF7OTYs\n1FmIiNjHkWGhjygXEbGXI8NCnYWIiL2GHBadnZ089NBDfOlLX8Lr9RIKhWhvb6e4uBi3282CBQvo\n6OiwXr99+3aysrLwer0cOXLEGo9Go+Tl5ZGZmcmGDRus8Z6eHsrKysjIyCAQCHDmzJlr1qJfqyoi\nMryGHBaPP/44brebY8eOcezYMbKzsykvL8ftdtPY2Eh6ejo7d+4E4OzZs+zYsYNDhw5RXl7O6tWr\nrXXWrl3LunXrqKuro7a2lvr6egCqqqpoa2sjGo1SVFTEli1brlmLOgsRkeE15LA4ePAgjz32GGPG\njGHkyJGMGzeOcDhMWVkZycnJlJaWEgqFAAiFQhQVFeF2uyksLMQYY3UdJ06cYMmSJUycOJFFixb1\nm1NSUsLYsWNZuXKlNX41OmchIjK8hhQW77//Pt3d3axatQq/388TTzxBV1cXdXV1ZGdnA5CdnU04\nHAZiO/6cnBxrvsfjIRQK0dTURGpqqjXu9Xo5evQoAOFwGK/XC8CECRNobm7mwoULV61HnYWIyPAa\nOZRJ3d3dvPvuuzz11FPcf//9PPLII+zZswdjzGdew+VyXTFmjLHGjTH91hto7d7eTcSPUgUCAUaP\nDqizEBH5RDAYJBgM3tAaQwqLKVOm4PF4mD9/PgDLli3jxz/+MT6fj2g0Sm5uLtFoFJ/PB4Df7+fg\nwYPW/OPHj+Pz+UhJSaG5udkab2howO/3W3MaGhrweDy0traSlpZGcnLyVesZNWoTmzZdfvzv/67O\nQkQkLhAIEAgErMebN28e9BpDPmeRlZVFKBSit7eX119/nfvvvx+/309FRQVdXV1UVFRQUFAAQH5+\nPjU1NZw6dYpgMEhSUhIpKSlA7HBVZWUlLS0tVFVV9QuL3bt309nZya5du6y1rqbvlVDxx+osRETs\nM6TOAuDpp5/mb//2b+nu7ub+++9n6dKl9Pb2UlJSgsfjIS8vjyeeeAKAtLQ0Vq1axdy5cxk9ejTP\nP/98v3VKSkp49NFHWbp0KbNnzwZg4cKFVFdXk5OTQ2ZmJpWVldespe/5ivhjdRYiIvZxmcGcaEhA\nLpeLL37RcPbs5bHaWvjOd2L3IiLSn8vlGtQ5ZnDIO7iv1lnoMJSIiH0cGRb6iHIREXs5MizUWYiI\n2MsRYXG1q6HUWYiI2McRYaHOQkRkeDkyLNRZiIjYy5Fhoc5CRMRejgwLdRYiIvZyZFiosxARsZcj\nwkJXQ4mIDC9HhIU6CxGR4eXIsBgxInZ/6dLnX4uIiBM5MiziY+ouRETs4diw0HkLERH7ODYs1FmI\niNjH0WGhzkJExB6OCItPXzobH1NnISJiD0eEhToLEZHh5diwUGchImIfx4aFOgsREfs4NizUWYiI\n2MexYaHOQkTEPo4IC10NJSIyvIYcFpcuXSI3N5f58+cD0N7eTnFxMW63mwULFtDR0WG9dvv27WRl\nZeH1ejly5Ig1Ho1GycvLIzMzkw0bNljjPT09lJWVkZGRQSAQ4MyZMwPWos5CRGR4DTksnn32Wbxe\nLy6XC4Dy8nLcbjeNjY2kp6ezc+dOAM6ePcuOHTs4dOgQ5eXlrF692lpj7dq1rFu3jrq6Ompra6mv\nrwegqqqKtrY2otEoRUVFbNmyZcBadM5CRGR4DSks3n//fd544w0efvhhjDEAhMNhysrKSE5OprS0\nlFAoBEAoFKKoqAi3201hYSHGGKvrOHHiBEuWLGHixIksWrSo35ySkhLGjh3LypUrrfFr0cd9iIgM\nryGFxT//8z/z1FNPkZR0eXpdXR3Z2dkAZGdnEw6HgdiOPycnx3qdx+MhFArR1NREamqqNe71ejl6\n9CgQCx6v1wvAhAkTaG5u5sKFC9esRx8kKCIyvEYOdsJrr71Gamoqubm5BINBazzeYXwW8UNXfRlj\nrHFjTL/1rrf2K69s4uTJ2J8DgQCBQECdhYjIJ4LBYL/99VAMOix+85vf8Morr/DGG2/Q3d3NuXPn\nWL58OT6fj2g0Sm5uLtFoFJ/PB4Df7+fgwYPW/OPHj+Pz+UhJSaG5udkab2howO/3W3MaGhrweDy0\ntraSlpZGcnLyNWtatmwT3/xm/zF1FiIiMfH/RMdt3rx50GsM+jDU1q1bOX36NCdPnqSyspK5c+fy\nk5/8BL/fT0VFBV1dXVRUVFBQUABAfn4+NTU1nDp1imAwSFJSEikpKUDscFVlZSUtLS1UVVX1C4vd\nu3fT2dnJrl27rLWuRecsRESG1w2/zyJ+6GjVqlWcOnUKj8fDBx98wN///d8DkJaWxqpVq5g7dy7/\n8A//wLPPPmvNffrpp3nyySfx+XzMmTOH2bNnA7Bw4ULGjRtHTk4O1dXVbNy4ccAadM5CRGR4ucxg\nTjYkIJfLxS9/afja1/qPr1kDkyfH7kVE5DKXyzWo88zgkHdwq7MQERlejg0LnbMQEbGPI8LiWp8N\npc5CRMQejggLdRYiIsPLsWGhzkJExD6ODQt1FiIi9nFsWKizEBGxj2PDQp2FiIh9HBsW6ixEROzj\niLC42qWz6ixEROzjiLBQZyEiMrwcGxbqLERE7OPosFBnISJiD0eERdJVtmL0aHUWIiJ2cURYXI06\nCxER+zg2LNRZiIjYx7Fhoc5CRMQ+jg0LdRYiIvZxbFjo0lkREfs4Niz0pjwREfs4NizUWYiI2Mex\nYaHOQkTEPo4NC3UWIiL2GVJYnD59mvvuu4+pU6cSCAR46aWXAGhvb6e4uBi3282CBQvo6Oiw5mzf\nvp2srCy8Xi9HjhyxxqPRKHl5eWRmZrJhwwZrvKenh7KyMjIyMggEApw5c2ZQNaqzEBGxz5DCYtSo\nUWzbto1IJMK+ffvYuHEj7e3tlJeX43a7aWxsJD09nZ07dwJw9uxZduzYwaFDhygvL2f16tXWWmvX\nrmXdunXU1dVRW1tLfX09AFVVVbS1tRGNRikqKmLLli2DrFGdhYiIXYYUFnfeeSezZs0C4I477mDq\n1KnU1dURDocpKysjOTmZ0tJSQqEQAKFQiKKiItxuN4WFhRhjrK7jxIkTLFmyhIkTJ7Jo0aJ+c0pK\nShg7diwrV660xj+r+JvyjBnKFoqISF83fM6iqamJSCRCfn4+dXV1ZGdnA5CdnU04HAZiO/6cnBxr\njsfjIRQK0dTURGpqqjXu9Xo5evQoAOFwGK/XC8CECRNobm7mwoULn33DkmDECLh06Ua3UERERt7I\n5Pb2dpYsWcK2bdu4/fbbMYP4b7zL5bpizBhjjRtj+q030NqbNm2y/hwIBAgEAsDl8xYjb2grRURu\nbcFgkGAweENrDHk32tPTw4MPPsjy5cspLi4GwOfzEY1Gyc3NJRqN4vP5APD7/Rw8eNCae/z4cXw+\nHykpKTQ3N1vjDQ0N+P1+a05DQwMej4fW1lbS0tJITk6+ai19w6IvnbcQEen/n2iAzZs3D3qNIR2G\nMsZQVlbGtGnTWLNmjTXu9/upqKigq6uLiooKCgoKAMjPz6empoZTp04RDAZJSkoiJSUFiB2uqqys\npKWlhaqqqn5hsXv3bjo7O9m1a5e11mDoiigREZuYITh8+LBxuVxm5syZZtasWWbWrFlm//795ty5\nc+aBBx4wkyZNMsXFxaa9vd2a88wzz5h77rnH5OTkmDfffNMaj0QiJjc310yePNmsX7/eGr948aJZ\nsWKFmTRpkiksLDQffvjhVWsZaBPuusuYDz4YyhaKiDjXUHb9rk8m3rJcLtc1z2dMngzBYOxeRERi\nBtpvXotj38ENOmchImIXR4eFzlmIiNjD0WGhzkJExB6ODgt1FiIi9nB0WKizEBGxh6PDQp2FiIg9\nHB0W6ixEROzh+LBQZyEicuMcHRajR6uzEBGxg6PDQp2FiIg9HB0W6ixEROzh6LDQCW4REXs4Oix0\n6ayIiD0cHRbqLERE7OHosJgwAd5//2ZXISJy63N0WCxYAHv3wq39GztERG4+R4dFbm7sUFQ4fLMr\nERG5tTk6LFwu+Ou/hpdeutmViIjc2hz9a1UBGhthzpzYuYuRIz/HwkREEpR+repVZGXBpEnwq1/d\n7EpERG5djg8L0KEoEZEb5fjDUAD/8z8wdSp8+CGMGfM5FSYikqB0GOoa/t//g7w8eOONm12JiMit\nKaHD4s033yQnJ4esrCx++MMf3tBaN/NQVDAYvDlfeJBUp71Up71U582V0GHxT//0Tzz//PMcPHiQ\n5557jpaWliGvtWgRHDgAbW02FvgZ3So/PKrTXqrTXqrz5krYi0nbPtmrf/WrXwVg3rx5hEIh/uqv\n/mpI640fD/fdBz/6ETz8MNx+u22litju0iXo6oKPP479ubc3dt/TE/twzAsXYvc9PbHX9L311dsb\ne0183qdff+lS7E2rzz0Xe1+SyxUb7+qC7u7YfU/P5edcLkhKunwbMSJ273LFvl78ftSo2Ad5xm/Q\n/+v29vZfJ752nMsFycmxc4zxW2MjvP765U9kSEqKXQ4/YkTsPv41k5Nj9yNHxr6OMZfv468ZNSp2\nGzHi8jbE1xk5sn8tEpOwYVFXV0d2drb12Ov1cvTo0SGHBcCaNfDII/DYY7EfmPR0SE2N/XBc64f2\n0z7rOaG+a7z7LtTXD7nsQX3dgXx6uz69ZmMj1NUNrYZP/0O3y9W+5rXqvJHvUd8dyqfXiW9P/Ln4\n7dKlyzvc+H3f20cfwU9/Gls3vvanv2b8Fn9NdzecPx/bQY8Zc3mHFt+ZjRp1eWfYd6cX38mNGHHl\n9z/+uk+/Pr6jbGmBSORyfSNGwG23xW4pKbHXw5U73niIxT+sMz7fGGhvj4VTPNhcrss1xnfGfbe7\nt7d/zb29sbnx0Orujv29/+//Xv432tvbP4DigRi/ffxx/3CDK4Mzvg3xMI6vNXLk5cCJf5/iYdL3\nPv797rvv+OgjePXVK/cpfW/xn6vr3fq+ru+fKyrgrrsG9zN+oxL2aqiDBw/y4osv8rOf/QyAnTt3\n8sEHH/Av//Iv/V7n0n8BREQGbbC7/oTtLHw+H9/+9retx5FIhKKioitel6BZJyLiKAl7gnvcuHFA\n7Iqo9957jwMHDuD3+29yVSIif5oStrMAeOaZZ3jkkUfo6elh9erV3HHHHTe7JBGRP0kJ21kAFBYW\nEo1GaWpqYvXq1f2es/M9GHYpLS0lLS2N6dOnW2Pt7e0UFxfjdrtZsGABHR0dN7HCmNOnT3Pfffcx\ndepUAoGH1iEnAAAFl0lEQVQAL33yBpREq7W7uxu/38+sWbMoKChg27ZtCVln3KVLl8jNzWX+/PlA\nYtY5efJkZsyYQW5uLvn5+UBi1tnZ2clDDz3El770JbxeL6FQKKHqPHHiBLm5udZt3LhxbN++nY6O\njoSpMe6FF17gy1/+Mvfeey9r1qwBhvZ3ntBhMRA734NhlxUrVlBdXd1vrLy8HLfbTWNjI+np6ezc\nufMmVXfZqFGj2LZtG5FIhH379rFx40ba29sTrtYxY8bwq1/9iv/6r/+itraWF198kcbGxoSrM+7Z\nZ5/F6/VaF10kYp0ul4tgMMhbb71F+JNf9JKIdT7++OO43W6OHTvGsWPHyM7OTqg6PR4Pb731Fm+9\n9Rb/+Z//ydixY1m4cCE7duxImBoBWltb2bp1KwcOHKCuro53332XmpqaIX0vb8mw6PsejIyMDOs9\nGDfbnDlzGD9+fL+xcDhMWVkZycnJlJaWJkSdd955J7NmzQLgjjvuYOrUqdTV1SVkrWPHjgWgo6OD\njz/+mOTk5ISs8/333+eNN97g4Ycfti66SMQ64cqLQhKxzoMHD/LYY48xZswYRo4cybhx4xKyTojV\nOmXKFCZNmpRwNd52220YY2hra6Orq4vz58/zhS98YWh1mlvQgQMHzNKlS63H5eXlZuPGjTexostO\nnjxppk2bZj12u92mq6vLGGNMZ2encbvdN6u0q2psbDR33323aW9vT8haL126ZGbMmGFGjBhhfvjD\nHxpjEvN7+s1vftP87ne/M8Fg0HzjG98wxiRmnXfffbeZMWOGKS4uNr/4xS+MMYlX5+nTp43H4zEP\nPfSQyc/PN9///vfN+fPnE67OuBUrVpjnnnvOGJN430tjjHnjjTfMqFGjzO23324ee+wxY8zQ6rwl\nO4tbiUngS3vb29tZsmQJ27Zt4/bbb0/IWpOSkvj9739PU1MTO3bs4K233kq4Ol977TVSU1PJzc3t\nV1ui1Qnw61//mt///vd873vf41vf+hZnzpxJuDq7u7t59913efDBBwkGg0QiEfbs2ZNwdQJcvHiR\nV199lcWLFwOJ93f+0UcfsWrVKhoaGnjvvff47W9/y2uvvTakOm/JsPD5fBw/ftx6HIlEKCgouIkV\nXZvP5yMajQIQjUbx+Xw3uaKYnp4eHnzwQZYvX05xcTGQuLVC7MTsX/7lXxIKhRKuzt/85je88sor\n3H333Sxbtoz/+I//YPny5QlXJ8Bdn7ztNycnhwceeIBXX3014eqcMmUKHo+H+fPnc9ttt7Fs2TKq\nq6sTrk6A/fv3c++99/LFL34RSLx/Q+FwmIKCAqZMmcLEiRNZvHgxhw8fHlKdt2RY3ErvwfD7/VRU\nVNDV1UVFRUVChJoxhrKyMqZNm2ZdHQGJV2tLSwv/93//B8Af//hHfvnLX1JcXJxwdW7dupXTp09z\n8uRJKisrmTt3Lj/5yU8Srs7z58/T3t4OxP7HWVNTQ1FRUcLVCZCVlUUoFKK3t5fXX3+d+++/PyHr\n/NnPfsayZcusx4lW45w5c6ivr6e1tZULFy6wf/9+5s2bN7Q6h+EQ2eciGAya7Oxsc88995hnn332\nZpdjjDFm6dKl5q677jKjR4826enppqKiwpw7d8488MADZtKkSaa4uNi0t7ff7DLN4cOHjcvlMjNn\nzjSzZs0ys2bNMvv370+4Wo8dO2Zyc3PNjBkzzLx588y//du/GWNMwtXZVzAYNPPnzzfGJF6d//3f\n/21mzpxpZs6caebOnWtefPHFhKzTGGNOnDhh/H6/mTlzplm7dq3p6OhIuDo7OjrMxIkTzblz56yx\nRKvRGGN+9KMfma9+9atm9uzZZuPGjebSpUtDqjNhPxtKREQSxy15GEpERD5fCgsREbkuhYWIiFyX\nwkJERK5LYSEiItelsBARkev6/xWik0C1O5aFAAAAAElFTkSuQmCC\n",
       "text": [
        "<matplotlib.figure.Figure at 0x2fed990>"
       ]
      }
     ],
     "prompt_number": 14
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "But if you zoom in on the y axis, you can see that there are a bunch of k-mers that are high abundance -- that is, there are approximately 550 different k-mers  that show up 50 different times in the data set."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "plot(xx, yy)\n",
      "axis(ymax=600)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "pyout",
       "prompt_number": 15,
       "text": [
        "(0.0, 80.0, 0.0, 600)"
       ]
      },
      {
       "output_type": "display_data",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAXgAAAD9CAYAAAC2l2x5AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3Xt4VNW9//F3IhBQKKJoQMMISMiFSwiaS7WYmIMYbTFU\npEIrUhNbjbboUfu0WqzQi5XaIwVrAnqMemzVn8XWK4IEHRTUSaoICkFiBRMQCReVCYRb2L8/lgn3\nJJPMzN575vN6njwhc9vfhOSTle/ae60Yy7IsREQk4sTaXYCIiISGAl5EJEIp4EVEIpQCXkQkQing\nRUQilAJeRCRCtRrwu3btYsqUKQwePJjU1FR8Ph9+v5+CggI8Hg/jxo2jvr6++fFz5swhMTGR1NRU\nli1bFtLiRUTkxFoN+HvuuQePx8OqVatYtWoVycnJlJaW4vF4qK6uJiEhgblz5wJQV1dHSUkJS5Ys\nobS0lKlTp4b8ExARkeNrNeDLy8u566676Nq1K506daJnz55UVFRQVFREXFwchYWF+Hw+AHw+H/n5\n+Xg8HnJycrAsC7/fH/JPQkREjtViwG/cuJE9e/ZQXFxMVlYWM2fOpKGhgcrKSpKTkwFITk6moqIC\nMAGfkpLS/PykpKTm+0REJLw6tXTnnj17WLduHffffz+jR4/mhhtu4NlnnyWQ1Q1iYmLadJuIiLQu\nkPxtcQQ/aNAgkpKSGDt2LN26dWPSpEksXLiQjIwMqqqqAKiqqiIjIwOArKws1qxZ0/z8tWvXNt93\nvCLb+hYXZ9HQ0PbHB+vtnnvuCfsxVaf9b26o0w01qs7gvwWq1R58YmIiPp+PgwcP8sorrzB69Giy\nsrIoKyujoaGBsrIysrOzAcjMzGTRokXU1NTg9XqJjY2lR48eARd1PO343EREolqLLRqAP//5z1x7\n7bXs2bOH0aNHM3HiRA4ePMg111xDUlISI0eOZObMmQDEx8dTXFxMXl4eXbp0Yd68eUEpUh0dEZHA\ntRrwgwcP5t133z3m9hdeeOG4j7/lllu45ZZbOl6ZA+Tm5tpdQpuozuByQ51uqBFUp91irPY0djp6\n0JiYgPpJ3brB9u1w8skhLEpExOECzU5XLFWgFo2ISOBcEfCgSVYRkUC5IuA1ghcRCZwrAl5ERALn\nmoBXi0ZEJDCuCHi1aEREAueKgBcRkcC5JuDVohERCYwrAl4tGhGRwLki4EEjeBGRQLki4DWCFxEJ\nnCsCXkREAueagFeLRkQkMK4IeLVoREQC54qAFxGRwLkm4NWiEREJjCsCXi0aEZHAuSLgQSN4EZFA\nuSLgNYIXEQmcKwJeRAJzxx3wxBN2VyF2c03Aq0Uj0nb/+pcJ+T/8QT870cwVAa8WjUjbNTbCxo1Q\nWQnz50NxMRw4YHdVYgdXBDxoFCLSVhs3whlnQP/+sHQp/Oc/cOWVsHu33ZVJuLki4DWCF2m79eth\nwADz7299C155Bbp0gbvvtrcuCT9XBLyItN3hAQ8m3B94AB57DL7+2r66JPxcE/Bq0Yi0zdEBD+Dx\nwGWXwSOP2FOT2MMVAa8WjUjbHS/gAW6/HebMgf37w1+T2KPVgO/fvz/Dhw8nPT2dzMxMAPx+PwUF\nBXg8HsaNG0d9fX3z4+fMmUNiYiKpqaksW7YsdJWLyHGdKOBHjoRzz4V//CP8NYk9Wg34mJgYvF4v\nK1asoKKiAoDS0lI8Hg/V1dUkJCQwd+5cAOrq6igpKWHJkiWUlpYyderUoBWqFo1I25wo4MGM4v/n\nf/TzFC3a1KKxjvpuqKiooKioiLi4OAoLC/H5fAD4fD7y8/PxeDzk5ORgWRZ+v7/DRapFI9I2e/bA\ntm1w9tnHv//yy83pkkuXhrcusUebRvB5eXmMGzeOF198EYDKykqSk5MBSE5Obh7Z+3w+UlJSmp+b\nlJTUfF9HacQh0roNG6BfPzjppOPfHxsLt90Gf/5zWMsSm3Rq7QHLly+nb9++VFVVMXbsWDIzM48Z\n0bck5gTD7+nTpzf/Ozc3l9zc3BZeo82HE4lqLbVnmkyebM6Jr6qCw8Zj4kBerxev19vu57ca8H37\n9gUgJSWFK664gpdeeomMjAyqqqpIT0+nqqqKjIwMALKysigvL29+7tq1a5vvO9rhAS8iwbF+PQwc\n2PJjunaFm2+GoiIoLYW0tPDUJoE7evA7Y8aMgJ7fYotm9+7dzT30rVu3smjRIvLz88nKyqKsrIyG\nhgbKysrIzs4GIDMzk0WLFlFTU4PX6yU2NpYePXoE+Ckdn1o0Iq1rywge4Fe/gkmT4NJL4Uc/gk8/\nDX1tEn4tBvyWLVsYNWoUI0aMYOLEidx+++3069eP4uJiampqSEpKYtOmTdx4440AxMfHU1xcTF5e\nHjfddBOzZ88OSpFq0Yi0TVsDvnNn+PnPoboakpIgIwN+//vQ1yfhFWMF0lAP1kFjYgLq48fHw6pV\n5r2InNh555m2yzeXrLTZ+vXmPPkdOzSgcrJAs9MVV7KCWjQibdHWEfzRBgwwvfna2uDXJPZxRcBr\nRCHSuq+/hn37oHfv9j0/LQ1WrgxuTWIvVwQ8aAQv0pqm0Xt7B0QK+MjjioDXCF6kde1tzzRJS4MP\nPghePWI/VwS8iLTu0087HvAawUcW1wS8WjQiLevoCD4pCTZtgsMWhz3Cc8/BV1+1//Ul/FwR8GrR\niLSuowHfqROkpsKHHx573969ZomDX/+6/a8v4eeKgAeN4EVa05ZlClpzojbNO++YXx7z56tP7yat\nrkXjBBrBi7TMssxKkh0ZwcOJJ1oXL4Zx46B/f/jZz+Ctt/Rz6QauGcGLyIlt2QKnnALdu3fsdU40\ngi8vh9GjobDQrDn/97937DgSHq4JeLVoRE6so/33Jmlp8NFHcPDgodu+/NIsLXzBBWad+b/+FX75\nS9i5s+PHk9ByRcDrT0GRlgUr4E89FU4/Hf7zn0O3vf46XHghxMWZj7OzYcwY+N3vOn48CS1X9OBF\n5JC6Ovjxj+GBB+CbjdWCFvBwqE2TmGg+bmrPHO6++2DoUNMS6tLl0O2XXmoWLRNncMUIHtSiEWny\nyitmhH3RRfDNLpohCfgmixfDJZcc+Zj4eHjqKdOP37nTvL39tmnfiHO4YgSvFo3IIQsWwJ13mtH7\nhAmwYoUJ/B/8IDivn5YGjz9u/r1+Pfj9ZrR+tEsuOTL4lyzRmvJOoxG8iIvs329aJvn5phdeWWlG\n2F5v8EbwI0YcGsE3tWdi25AUAwZoZyincUXAawQvYrz9Npx7LvTpYz7u08dMgj7xRMcvcmoyYIBZ\nkmDHDvPL4+j++4n06wdffGGWLBZncEXAi4ixYAFcfvmRt3XpAtdea05hDIbYWBg2zFzw9Prrx/bf\nT6RzZzjrLKipCU4d0nGuCXi1aETg1VePDfhQaOrDn3EGJCS0/XkDB5q+vTiDKwJeLRoRs53e5s1m\ng+xQS0uDp59ue3umifrwzuKKgBcRM3ofMyZ4rZiWjBgBBw60vT3TRCN4Z3FNwKtFI9HueP33UBk6\n1FzRmpMT2PM0gncWVwS8WjQS7fbuhTfeMFeKhsMpp5gzYnr2DOx5GsE7iysCHjSCl+j21ltmM47e\nvcN3zE7tuAxSI3hncUXAawQv0S6c7ZmOOOMM89fG11/bXYmASwJeJNqF6/TIjoqJMaN4tWmcwTUB\nrxaNRKtPPzVrsqen211J26gP7xyuCHi1aCSaLV0KeXltWw/GCdSHd442fcs0NjaSnp7O2LFjAfD7\n/RQUFODxeBg3bhz19fXNj50zZw6JiYmkpqaybNmy0FQtEkXeey88FzcFi0bwztGmgJ89ezapqanE\nfDOULi0txePxUF1dTUJCAnPnzgWgrq6OkpISlixZQmlpKVOnTg1aoWrRSLT697/hvPPsrqLtNIJ3\njlYDfuPGjSxYsIDrr78e65uUraiooKioiLi4OAoLC/H5fAD4fD7y8/PxeDzk5ORgWRZ+v7/DRapF\nI9Fq/3748EP39N9BI3gnaTXg//u//5v777+f2MMagJWVlSR/s1dYcnIyFRUVgAn4lJSU5sclJSU1\n39dRGsFLNFqzBjwe6NHD7krarn9/2LDhyI27xR4tXsrw8ssvc+aZZ5Keno7X622+3QogbWNOMPye\nPn16879zc3PJzc1t4TXafDiRiPLee3D++XZXEZhTTjFXwH7xhVk+WNrP6/Uekb2BajHg3377bV58\n8UUWLFjAnj172LlzJ5MnTyYjI4OqqirS09Opqqoi45sZoKysLMrLy5ufv3bt2ub7jnZ4wIvI8bmt\n/96k6Vx4BXzHHD34nTFjRkDPb7FFc++991JbW8v69et55plnyMvL48knnyQrK4uysjIaGhooKysj\nOzsbgMzMTBYtWkRNTQ1er5fY2Fh6BOlvS7VoJBq99547A37gQE20OkFAq000tVuKi4u55pprSEpK\nYuTIkcycOROA+Ph4iouLycvLo0uXLsybNy8oRapFI9Fo/3746CN3TbA20dWszhBjBdJQD9ZBY2IC\n6uMnJsIrr8DgwSEsSsRhPvgAfvhDM9HqNmVl8OabZlcoCZ5As9MV18ZpBC/RyI0TrE00gncGVwS8\nSDRy6wQrqAfvFK4JeE2ySrRx6wQrmI266+rM0sFiH1cEvFo0Em327YPVq83eqG500knQrx989tmR\nt2/dak890coVAS8SbVavNleEdu9udyXtd3gf/uBBuO028zk1NNhaVlRxTcCrRSPRxM0TrE2a+vB7\n98KkSWZOITnZnF0j4eGKgFeLRqKNmydYmwwYACtWmI3CDx6E116DK6+EhQvtrix6uCLgQSN4iS5u\nnmBtMnAgPPIIpKXBM89A166Qn6+ADydXBLxG8BJN3D7B2uSSS0yw/+UvZtIVzFW527eb1SYl9FwR\n8CLR5KOPzOj3lFPsrqRjevWCq68+coAWG2taNosW2VdXNAloLRo7qUUjkWrXLvj4Y9iyxbyVl7t/\ngrUl+fkwfz7ccIPdlUQ+VwS8WjQSySZOhOpqOOcciI+Hvn3huuvsrip0xoyBm282raguXeyuJrK5\nIuBFItWePeD1Qm0tnHqq3dWExxlnmAUE33kHcnLsriayuaYHrxaNRKLly2HIkOgJ9yY6myY8XBHw\natFIpCovN2ebRBsFfHi4IuBBI3iJTIsXR2fAZ2WZUyU3b7a7ksjmioDXCF4i0fbtsG4dfLPjZVTp\n1AlGjzZXt0rouCLgRSLR66/DqFHReyaJ2jSh55qAV4tGIs3ixWYUG60uvdR8DRob7a4kcrki4NWi\nkUgUrROsTRIS4LTToKrK7koilysCHjSCl8jyn/+YddGHDLG7EnsNHw4ffmh3FZHLFQGvEbxEmvJy\n056J9u/tYcMU8KHkioAXiTTRenrk0RTwoeWagFeLRiJFY6M5g+a//svuSuyngA8tVwR8tP8ZK5Hl\n/ffNgmJnn213JfYbONBsxL1zp92VRCZXBLxIJIn2s2cOd9JJkJpq1sCX4HNNwKtFI5Ei2s9/P5ra\nNKHjioBXi0YihWWZDbUvuMDuSpxDAR86LQb8nj17yMrKYsSIEWRnZzNr1iwA/H4/BQUFeDwexo0b\nR319ffNz5syZQ2JiIqmpqSxbtixohWoEL5Fg40azFd9pp9ldiXMo4EOnxYDv2rUrb7zxBh988AFL\nly7l0Ucfpbq6mtLSUjweD9XV1SQkJDB37lwA6urqKCkpYcmSJZSWljJ16tSgFBkXZ3Z/EXG7NWtM\nz1kOaQp4DeKCr9UWzcknnwxAfX09Bw4cIC4ujoqKCoqKioiLi6OwsBCfzweAz+cjPz8fj8dDTk4O\nlmXh9/s7XOTJJ8Pu3R1+GRHbrV6tgD9afLxZXfLzz+2uJPK0GvAHDx4kLS2N+Ph4fvazn+HxeKis\nrCQ5ORmA5ORkKioqABPwKSkpzc9NSkpqvq8jFPASKdas0fIEx6M2TWi0uidrbGwsK1euZMOGDVx+\n+eVceOGFWAH8LRVzghnS6dOnN/87NzeX3NzcE76GAl4ixZo1cO21dlfhPE0Bn59vdyXO4vV68Xq9\n7X5+mzfd7t+/P5dffjk+n4+MjAyqqqpIT0+nqqqKjIwMALKysigvL29+ztq1a5vvO9rhAd8aBbxE\nAstSi+ZEhg2DpUvtrsJ5jh78zpgxI6Dnt9ii2bZtG1999RUA27dv57XXXqOgoICsrCzKyspoaGig\nrKyM7G+2pMnMzGTRokXU1NTg9XqJjY2lR48eAX5Kx1LASyT4/HPo2hV697a7EudRiyY0WhzBb968\nmSlTptDY2EifPn2444476Nu3L8XFxVxzzTUkJSUxcuRIZs6cCUB8fDzFxcXk5eXRpUsX5s2bF5Qi\nFfASCXQGzYkNGQIffwwHDpgJVwmOGCuQhnqwDhoTE1Af/+67zbZmd98dwqJEQmz2bLMH60MP2V2J\nMw0aBC+9BIedp8Ff/mLW7Jkwwb66nCTQ7HTFlawawUskWL1aZ9C05Og2zXvvwR13wPz59tXkdgp4\nkTBRi6Zlhwf8vn1w3XXwy1+apR2kfRTwImFgWQr41hwe8H/4A3g88NvfmuWEd+ywtza3UsCLhMEX\nX5ilcc880+5KnKsp4D/4AEpKYN488zUbOVKj+PZSwIuEgUbvrRs0CDZvhsmT4U9/OrQhyvnnK+Db\nSwEvEgYK+NZ16gTJyZCQAD/+8aHbMzIU8O2lgBcJA51B0zazZ8Pjjx+5B8T550NlpW0luZoCXiQM\nNIJvm1GjzOqShxs4EHbtMvMYEhgFvEiIaQ2ajomJUR++vRTwIiFWV2feHz0ylbZTwLePAl4kxJra\nM9pbuP0U8O2jgBcJMfXfOy4jw0y0alu/wLgi4Lt1g4YG/eeKO+kMmo5LSDDvN260tw63cUXAx8aa\njbf37LG7EpHAaQTfcZpobR9XBDyoTSPupYAPjqY2jbSdAl4khD76CDp3hr597a7E/TSCD5wCXiSE\nSkrgJz/RGTTB0BTwmotrO9dsjqWAF7fZuROeftpMskrH9ekDp5wCn34K555rdzXuoBG8SIg8+SSM\nHg1nnWV3JZFDbZrAKOBFQsCyzN6rN99sdyWRRROtgVHAi4SA12v67jk5dlcSWcaMMXu07t9vdyXu\noIAXCYGSErjpJk2uBtv558OAAfD//p/dlbiDAl4kyDZtgiVLzM5EEny/+hXcdx8cPGh3Jc6ngBcJ\nsocfhkmT4FvfsruSyDRmjLm2YMECuytxPgW8SBDt3w+PPGLaMxIaMTGHRvHSMgW8SBC9847ZLFqL\ni4XW+PFmh6dly+yuxNkU8CJB9NlnMHiw3VVEvk6d4Be/0Ci+NQp4kSCqqQGPx+4qosOUKfDee7Bq\nld2VOFeLAV9bW8vFF1/MkCFDyM3N5amnngLA7/dTUFCAx+Nh3Lhx1NfXNz9nzpw5JCYmkpqayrIg\n/v2kgBc3qK1VwIdL165w663wpz/ZXYlztRjwnTt3ZtasWaxevZr58+czbdo0/H4/paWleDweqqur\nSUhIYO7cuQDU1dVRUlLCkiVLKC0tZerUqUErVAEvblBTA/362V1F9Cgqgpde0imTJ9JiwPfp04cR\nI0YA0Lt3b4YMGUJlZSUVFRUUFRURFxdHYWEhPp8PAJ/PR35+Ph6Ph5ycHCzLwu/3B6VQBby4gVo0\n4dW7t3n7+GO7K3GmNvfgP/nkE1avXk1mZiaVlZUkJycDkJycTEVFBWACPiUlpfk5SUlJzfd1lAJe\n3EABH35an+bE2rRcsN/v5+qrr2bWrFl0794dK4AFmWNOcK329OnTm/+dm5tLbm5ui6+jgBen+/pr\ns8hYz552VxJdMjNNwF97rd2VBJ/X68Xr9bb7+a0G/P79+xk/fjyTJ0+moKAAgIyMDKqqqkhPT6eq\nqoqMjAwAsrKyKC8vb37u2rVrm+872uEB3xYKeHG6ptG71p8Jr4yMyF2b5ujB74wZMwJ6fostGsuy\nKCoqYujQodx6663Nt2dlZVFWVkZDQwNlZWVkZ2cDkJmZyaJFi6ipqcHr9RIbG0uPHj0CKuhEFPDi\ndGrP2GPkSLM14r59dlfiPC0G/PLly/nb3/7G66+/Tnp6Ounp6SxcuJDi4mJqampISkpi06ZN3Hjj\njQDEx8dTXFxMXl4eN910E7Nnzw5aoQp4cbraWp1BY4dTToGBA+HDD+2uxHlirEAa6sE6aExMQH18\ngC+/NP+JX34ZoqJEOuiuu8xAZNo0uyuJPkVFZinh4mK7KwmtQLNTV7KKBIlaNPbRmTTH55qA79IF\nDhwwbyJOpBaNfRTwx+eagI+JMaP4hga7KxE5Po3g7TNsGHz6KRy2aorgooAHtWnEuRob4fPPISHB\n7kqiU5cuJuTff9/uSpxFAS8SBFu2QK9eEBdndyXRS22aYyngRYJA7Rn7ZWZCkFZGiRgKeJEgUMDb\nTyP4YyngRYJAAW+/wYNh+3bYts3uSpxDAS8SBDpF0n6xseZiJ43iD1HAiwSBRvDOoDbNkRTwIkGg\ngHcGBfyRFPAiQaAWjTM0nUkT/hW2nEkBL9JBDQ2wcyeceabdlUhCgskJrSxpKOBFOqi21gRLrKt+\nmiJTTAyMHw/z59tdiTO46ltSAS9OpPaMs0yYAP/4h9o0oIAX6TBNsDpLZibs2gVr1thdif0U8CId\npIB3lpgYuOoqM4qPdgp4kQ5Si8Z5mto00U4BL9JBGsE7T1aWObMp2ts0CniRDlLAO09srM6mAQW8\nSIdYlgl4tWicR20aBbxIh+zYYTb56NHD7krkaN/+tvn/WbvW7krso4AX6QC1Z5xLbRoFvEiHqD3j\nbNHeplHAi7STZcHTT5vNnsWZLrwQtm6FdevsrsQeCniRdrrnHvjsM/jNb+yuRE4kNha+/3345z/t\nrsQergr4bt1MwGuNCbHbY4/B3/8OL7xgvi/FucaPj96Aj7Gs8MdlTEwM7T1sXJy5gCEuLshFibTR\nkiXwox/B0qWQlGR3NdKaAwegTx94/333T4gHmp2dQlhLSDS1aRTwEg6ffQYbNhz6+Kuv4Cc/MRN3\nCnd36NQJrrgC/vUvuOUWu6sJrxZbNIWFhcTHxzPssFkkv99PQUEBHo+HcePGUV9f33zfnDlzSExM\nJDU1lWXLloWkYPXhJVx27DBbwN19t+mz/+Y38MADUFICOTl2VyeBuPLK6GzTtBjw1113HQsXLjzi\nttLSUjweD9XV1SQkJDB37lwA6urqKCkpYcmSJZSWljJ16tSQFKyAl3D5/e9NMLz5pmnHNL1ddZXd\nlUmgRo+GlSthyxa7KwmvFgN+1KhR9OrV64jbKioqKCoqIi4ujsLCQnw+HwA+n4/8/Hw8Hg85OTlY\nloXf7w96wQp4CYfqavi//4MZM+yuRIKha1e47DIzKR5NAj6LprKykuTkZACSk5OpqKgATMCnpKQ0\nPy4pKan5vmBSwEs4/PKXcMcdEB9vdyUSLNHYpgl4kjWQGdyYmJgT3jd9+vTmf+fm5pKbm9um11TA\nS6gtXWrOuHjqKbsrkWC67DIoKoIvv4SjGhOO5fV68Xq97X5+wAGfkZFBVVUV6enpVFVVkZGRAUBW\nVhbl5eXNj1u7dm3zfcdzeMAHQgEvoXTwINx2G9x3n/mzXiJH9+6QlwcvvwyTJ9tdTdscPfidEWDP\nMOAWTVZWFmVlZTQ0NFBWVkZ2djYAmZmZLFq0iJqaGrxeL7GxsfQIwRJ7CngJpb/9Dbp0gauvtrsS\nCYVoa9O0GPCTJk3iggsuYN26dfTr14/HHnuM4uJiampqSEpKYtOmTdx4440AxMfHU1xcTF5eHjfd\ndBOzZ88OScEKeAmV/fvh1782p0K20F0UF/ve98yFart22V1JeLjuStabb4bUVPNeJJjmz4e//hU6\n0PIUFxgzBn7wA7j+ersrCVzUXMkqEmzz5sENN9hdhYTazJlw6aXQv785Pz6SuWqxMVDAS2h88om5\nEObKK+2uREItPR2eew5++EN45x27qwktBbwI8MgjcO21WuMoWowaZS5kGzcOVq2yu5rQUcBL1Nu3\nDx5/HH76U7srkXDKz4cHHzTvq6vtriY01IOXqPf88zBkCAwebHclEm4/+AFs3mx+ub/xht3VBJ9G\n8BL1NLka3W66Cdavh2+W1YooCniJatXV8NFHZls3iU6dO5t1h2bOtLuS4FPAS1R7+GGYMsVcvSrR\nq7AQli+Hqiq7KwkuBbxErb174YknNLkqJld+/nO4/367KwkuTbJK1Covh+RkGDTI7krECW6+2Xwv\n1NZCv352VxMcGsFL1Hr+efXe5ZBeveC662DWLLsrCR7XrUVTWwsXXGDei7RXYyOcdZa5knHgQLur\nEafYtAmGDTOT76efbnc1xwo0OzWCl6j07rvQp4/CXY509tlmuYoHH7S7kuBQD16i0r/+ZS5TFzna\nnXfChRfCOeeYlo2buS7gu3Y1Zz8cPAixrvv7Q5zAskz//R//sLsScaJzzzXbNl52GWzcCNOmuXd/\nANdFZEwMdOsGDQ12VyJutXq12dxjxAi7KxGnSkqCt982f+ndcAMcOGB3Re3juoAHtWmkY55/3rRn\n3Doqk/Do08eM5D/7DAoKYMsWuysKnAJeok5TwIu0pkcPs0n3sGHm7eGHTXvYLRTwElVqa2HDBrMe\nuEhbdO4M991nLox77DHzvfPhh3ZX1TYKeIkqL7xgNl7u5LrTC8Ruw4eb9WquvRby8uCf/7S7ota5\n8ttcAS/t9fzz8LOf2V2FuFVsrJl0zcw0+7p272428XYqjeAlanz5JVRUOPsHUtwhPd2M4H/0I3O2\njVMp4CVq/PGPJtxPPtnuSiQSfOc78OSTZj2jlSvtrub4XBvw9fV2VyFO4/ebNWaO509/MmdDzJ0b\n3poksuXnw1//ai6K+ve/7a7mWK4M+NxceOQRc0WiCJgt15KT4dvfhg8+OPK+hx+G0lJYvBh697an\nPolcEybA7Nnw3e/CrbeagYZTuDLgCwvNCP7pp+2uRJxg82a45BK46y648UYz+XXHHeZ75JlnYMYM\nE+5nn213pRKpJkwwV0jv3AmpqeYKWCcMQF23XHCTt982O6JXVZmLESQ6bd8OOTlmsuvOO81tW7ea\ngH/9ddiaaW3HAAAKDElEQVS3z5y/PGyYvXVK9Fi61Aw0zjnHzPukpwfvtQPNTtcGPJi9NPv0iczN\ncqV1fj+MHm1advfdd+zSA2++CT17QlqaLeVJFNu3z7SRf/97uPhi+O1vg7NzmCPWg3/zzTdJSUkh\nMTGRB0O4sPLMmfDoo/Dxx6F5fa/XG5oXDrJIqHPTJnMRyddfH3tfYyOsW2cuUpo71/yw3HyzuaJw\nxIjjhzvARRe1L9zd8PV0Q40QvXV26WK+R6urTcsmOxvuuSeoh2iTkAT8Lbfcwrx58ygvL+ehhx5i\n27ZtoTgMffqYvuvUqYH1uyzLfOFraszSwycS7P/0zZvNcXfuPLbehgZzCf26dW1bue7AAbMb0b33\nwpw53hOeVbRrF3i9MGcO/OQnkJUFAwaYr5nPF/o+YWOjWaRpzRooL/cec//27fCLX5irBG+91fTJ\nzznHXG364x9DRgZ861vmbIVHHoH33zcrQaakwPTpUFIS/EXD3BBKbqgRVGf37ma54bVrzdxQuAX9\nStavvxmCXXTRRQCMGTMGn8/Hd7/73WAfCjA7of/v/5qrE2NjTZhs2WICMDXVBMfw4WZ7tuXLzWTb\nkiUQF2fCp67OnHYZHw+Jieaxw4aZ93V18NRTZt2JVatMMA8deugx/fubi2eajrljB5x6qnmtM880\nW36tXGmOWV5uAr5XL/PYgwfNYzp3Nh/v3Wue16kTfPGFOSNk2DDzPi7u0Oe7d6/ZjWjpUnP8iy+G\njz6ChAS4/HK45ho47TRzvMWL4b33TL1paTBypGlr9eoFzz1nLrlubISrrzbLo8bHm7feveHzz83n\nvGqV+fz37TP1Nj2mS5dDn3ddnfk6HP7L4sAB2Lbt0NekVy+zKl9l5aG2yuLF8MADcNVV5hhnnWW+\nLuvXm+Nu326uGhw6VPMs4m69e9tzBlfQA76yspLk5OTmj1NTU3n33XdDFvCdO5szJZ591gTQd75j\nAigmxsxqf/ihub+21pxCd8kl8Ic/mHAEE0pffmlCdd06EyzPPWf+nNqxw0zYDR8OP/2pGUl+9JEZ\nRT7+uPkL4PTTDwX6aafBV1+ZwNuyxQRcSoo55hNPmMmWk04yx921yzxm/37z/J49D41E6+sP1X70\niD42FiZONKPZM880t/Xsaf4cfPZZ87nV15sQvfNO08Y45ZRjv25DhsDdd5tfAP/8J7z22qG66+qg\nb99Dv8i+9z3zS/DwQK+vN9vdffvbpv5evY7cgCU21nxDn3HGoXVffvUr8/jFi80vmmHDzGT54MFH\nPu/cc82biHRM0CdZy8vLefTRR3n6m3MY586dy6ZNm/jd73536KBaiFtEpF0Cieygj+AzMjL4xS9+\n0fzx6tWryc/PP+IxNpy4IyISdYI+ydqzZ0/AnEmzYcMGFi9eTFZWVrAPIyIirQjJcsF/+ctfuOGG\nG9i/fz9Tp06lt64PFxEJu5CcJpmTk0NVVRWffPIJU6dObb49XOfHB6qwsJD4+HiGHXa5o9/vp6Cg\nAI/Hw7hx46h3wOpmtbW1XHzxxQwZMoTc3FyeeuopwHm17tmzh6ysLEaMGEF2djazZs1yZJ0AjY2N\npKenM3bsWMCZNfbv35/hw4eTnp5OZmYm4Mw6d+3axZQpUxg8eDCpqan4fD7H1fnxxx+Tnp7e/Naz\nZ0/mzJlDfX29o+oEeOSRR7jgggs477zzuPXWW4HA/9/DuhZNuM6PD9R1113HwoULj7ittLQUj8dD\ndXU1CQkJzHXAMoSdO3dm1qxZrF69mvnz5zNt2jT8fr/jau3atStvvPEGH3zwAUuXLuXRRx+lurra\ncXUCzJ49m9TU1OaJfyfWGBMTg9frZcWKFVRUVADOrPOee+7B4/GwatUqVq1aRXJysuPqTEpKYsWK\nFaxYsYL33nuPk08+me9///uUlJQ4qs4dO3Zw7733snjxYiorK1m3bh2LFi0K+OsZtoA//Pz4c845\np/n8eCcYNWoUvXr1OuK2iooKioqKiIuLo7Cw0BG19unThxEjRgDQu3dvhgwZQmVlpSNrPfmbRdfr\n6+s5cOAAcXFxjqtz48aNLFiwgOuvv7554t9pNTY5+sQEJ9ZZXl7OXXfdRdeuXenUqRM9e/Z0ZJ1N\nysvLGTRoEP369XNcnd26dcOyLL7++msaGhrYvXs3p556auB1WmGyePFia+LEic0fl5aWWtOmTQvX\n4Vu1fv16a+jQoc0fezweq6GhwbIsy9q1a5fl8XjsKu24qqurrQEDBlh+v9+RtTY2NlrDhw+3Tjrp\nJOvBBx+0LMt5X9OrrrrKev/99y2v12t973vfsyzLeTValmUNGDDAGj58uFVQUGC98MILlmU5r87a\n2lorKSnJmjJlipWZmWndd9991u7dux1X5+Guu+4666GHHrIsy3lfT8uyrAULFlidO3e2unfvbt11\n112WZQVepyuXCw4Hy8Gncvr9fq6++mpmzZpF9+7dHVlrbGwsK1eu5JNPPqGkpIQVK1Y4qs6XX36Z\nM888k/T09CPqclKNTZYvX87KlSv54x//yG233cYXX3zhuDr37NnDunXrGD9+PF6vl9WrV/Pss886\nrs4m+/bt46WXXmLChAmA8/7ft27dSnFxMWvWrGHDhg288847vPzyywHXGbaAz8jIYO3atc0fr169\nmuzs7HAdPmAZGRlUVVUBUFVVRUZGhs0VGfv372f8+PFMnjyZgoICwLm1gpkgvPzyy/H5fI6q8+23\n3+bFF19kwIABTJo0iddff53Jkyc7qsYmffv2BSAlJYUrrriCl156yXF1Dho0iKSkJMaOHUu3bt2Y\nNGkSCxcudFydTV599VXOO+88zjjjDMB5P0MVFRVkZ2czaNAgTj/9dCZMmMBbb70VcJ1hC3i3nR+f\nlZVFWVkZDQ0NlJWVOeKXkWVZFBUVMXTo0OZZdXBerdu2beOrr74CYPv27bz22msUFBQ4qs57772X\n2tpa1q9fzzPPPENeXh5PPvmko2oE2L17N/5vtgjaunUrixYtIj8/33F1AiQmJuLz+Th48CCvvPIK\no0ePdmSdAE8//TSTJk1q/thpdY4aNYp///vf7Nixg7179/Lqq68yZsyYwOsMUfvouLxer5WcnGyd\ne+651uzZs8N56BZNnDjR6tu3r9WlSxcrISHBKisrs3bu3GldccUVVr9+/ayCggLL7/fbXab11ltv\nWTExMVZaWpo1YsQIa8SIEdarr77quFpXrVplpaenW8OHD7fGjBljPfHEE5ZlWY6rs4nX67XGjh1r\nWZbzavz000+ttLQ0Ky0tzcrLy7MeffRRR9ZpWZb18ccfW1lZWVZaWpp1++23W/X19Y6ss76+3jr9\n9NOtnTt3Nt/mxDofe+wx66KLLrLOP/98a9q0aVZjY2PAddqy4YeIiISeJllFRCKUAl5EJEIp4EVE\nIpQCXkQkQingRUQilAJeRCRC/X+q9Eh5R8iURAAAAABJRU5ErkJggg==\n",
       "text": [
        "<matplotlib.figure.Figure at 0x30f8110>"
       ]
      }
     ],
     "prompt_number": 15
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "## Questions --\n",
      "\n",
      "1. qualitatively, how does the high-abundance curve shift as you count k-mers using larger or smaller k?  Why?\n",
      "\n",
      "2. qualitatively, how does the mean of this bump shift with higher or lower read coverage?  Why?\n",
      "\n",
      "3. qualitatively, how does the mean of this bump shift with higher or lower error rate?  Why?\n",
      "\n",
      "Answer these questions by copying the code in the examples above, changing the key parameters, and producing three graphs (with associated discussion) for each set of parameters."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 15
    }
   ],
   "metadata": {}
  }
 ]
}