3.1 sgRNA Activity Scoring.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# SgRNA Activity Scoring\n",
    "\n",
    "After the database had been created and all sgRNAs identified, I scored the sgRNAs using published sgRNA scoring methods. \n",
    "\n",
    "## Doench Score\n",
    "\n",
    "The first method I used was the on-target activity scoring method developed by Doench et al. in <a href=\"http://www.nature.com/nbt/journal/v32/n12/full/nbt.3026.html\">this paper</a>. The calc_doench_score function was downloaded as a python script from the Broad Institute's <a href=\"http://www.broadinstitute.org/rnai/public/analysis-tools/sgrna-design-v1\">sgRNA Designer</a>. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import math\n",
    "\n",
    "### Doench Score ###\n",
    "\n",
    "\"\"\"\n",
    "Calculates the on-target score for an sgRNA\n",
    "Input: 30mer\n",
    "Output: On-target score\n",
    "Run as: python on_target_score_calculator.py <30mer>\n",
    "\"\"\"\n",
    "def calc_doench_score(s):\n",
    "    s_list = list(s)\n",
    "    s_20mer = s[4:24]\n",
    "    nuc_hash = {'A':0, 'T':1, 'C':2, 'G':3}\n",
    "    score = 0.597636154\n",
    "    gc = s_20mer.count('G')+s_20mer.count('C')\n",
    "    gc_low = -0.202625894\n",
    "    gc_high = -0.166587752\n",
    "    if gc < 10:\n",
    "        gc_val = abs(gc-10)\n",
    "        score = score+(gc_val*gc_low)\n",
    "    elif gc > 10:\n",
    "        gc_val = gc-10\n",
    "        score = score+(gc_val*gc_high)\n",
    "    #rows[1-30]cols['ATCG']\n",
    "    sing_nuc_hash = {'G2':-0.275377128,'A3':-0.323887456,'C3':0.172128871,'C4':-0.100666209,'C5':-0.20180294, \\\n",
    "                    'G5':0.245956633,'A6':0.036440041,'C6':0.098376835,'C7':-0.741181291,\\\n",
    "                    'G7':-0.393264397,'A12':-0.466099015,'A15':0.085376945,'C15':-0.013813972,\\\n",
    "                    'A16':0.272620512,'C16':-0.119022648,'T16':-0.285944222,'A17':0.097454592,\\\n",
    "                    'G17':-0.17554617,'C18':-0.345795451,'G18':-0.678096426,'A19':0.22508903,\\\n",
    "                    'C19':-0.507794051,'G20':-0.417373597,'T20':-0.054306959,'G21':0.379899366,\\\n",
    "                    'T21':-0.090712644,'C22':0.057823319,'T22':-0.530567296,'T23':-0.877007428,\\\n",
    "                    'C24':-0.876235846,'G24':0.278916259,'T24':-0.403102218,'A25':-0.077300704,\\\n",
    "                    'C25':0.287935617,'T25':-0.221637217,'G28':-0.689016682,'T28':0.117877577,\\\n",
    "                    'C29':-0.160445304,'G30':0.386342585}\n",
    "    #score_mat = np.matrix('0 0 0 0;0 0 0 -0.275377128;-0.323887456 0 0.172128871 0;0 0 -0.100666209 0;\n",
    "    #0 0 -0.20180294 0.245956633;0.036440041 0 0.098376835 0;0 0 -0.741181291 -0.393264397;0 0 0 0;0 0 0 0;\n",
    "    #0 0 0 0;0 0 0 0;-0.466099015 0 0 0;0 0 0 0;0 0 0 0;0.085376945 0 -0.013813972 0;\n",
    "    #0.272620512 -0.285944222 -0.119022648 0;0.097454592 0 0 -0.17554617;0 0 -0.345795451 -0.678096426;\n",
    "    #0.22508903 0 -0.507794051 0;0 -0.054306959 0 -0.417373597;0 -0.090712644 0 0.379899366;\n",
    "    #0 -0.530567296 0.057823319 0;0 -0.877007428 0 0;0 -0.403102218 -0.876235846 0.278916259;\n",
    "    #-0.077300704 -0.221637217 0.287935617 0;0 0 0 0;0 0 0 0;0 0.117877577 0 -0.689016682;\n",
    "    #0 0 -0.160445304 0;0 0 0 0.386342585')\n",
    "    dinuc_hash = {'GT2':-0.625778696,'GC5':0.300043317,'AA6':-0.834836245,'TA6':0.760627772,\n",
    "                  'GG7':-0.490816749,'GG12':-1.516907439,'TA12':0.7092612,'TC12':0.496298609,\n",
    "                  'TT12':-0.586873894,'GG13':-0.334563735,'GA14':0.76384993,'GC14':-0.53702517,\n",
    "                  'TG17':-0.798146133,'GG19':-0.66680873,'TC19':0.353183252,'CC20':0.748072092,\n",
    "                  'TG20':-0.367266772,'AC21':0.568209132,'CG21':0.329072074,'GA21':-0.836456755,\n",
    "                  'GG21':-0.782207584,'TC22':-1.029692957,'CG23':0.856197823,'CT23':-0.463207679,\n",
    "                  'AA24':-0.579492389,'AG24':0.649075537,'AG25':-0.077300704,'CG25':0.287935617,\n",
    "                  'TG25':-0.221637217,'GT27':0.117877577,'GG29':-0.697740024}\n",
    "    for i,nuc in enumerate(s_list):\n",
    "        key = nuc+str(i+1)\n",
    "        if sing_nuc_hash.has_key(key):\n",
    "            nuc_score = sing_nuc_hash[key]\n",
    "        else:\n",
    "            nuc_score = 0\n",
    "        #nuc_score = score_mat[i,nuc_hash[nuc]]\n",
    "        score = score+nuc_score\n",
    "        if i<29:\n",
    "            dinuc = nuc+s[i+1]+str(i+1)\n",
    "            if dinuc in dinuc_hash.keys():\n",
    "                score = score+dinuc_hash[dinuc]\n",
    "    partial_score = math.e**-score\n",
    "    final_score = 1/(1+partial_score)\n",
    "    return final_score"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The scoring function was then run on every long sgRNA."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "import data_processing as dp\n",
    "\n",
    "def run_doench(db_name, sql_version=\"MySQL\", firewall=False):\n",
    "    db_con = dp.DatabaseConnection(sql_version, db_name=db_name, firewall=firewall)\n",
    "    \n",
    "    rows = db_con.fetch_query(\"SELECT SgID, LongSg FROM SgRNATargetInformation\")\n",
    "    \n",
    "    doench_dict = {\"DoenchScore\": []}\n",
    "    sg_dict = {\"SgID\": [], \"LongSg\": []}\n",
    "    for row in rows:\n",
    "        if sql_version == \"MSSQL\":\n",
    "            sgID = row.SgID\n",
    "            longSg = row.LongSg\n",
    "        else:\n",
    "            sgID, longSg = row\n",
    "            longSg = str(longSg)\n",
    "        d_score = calc_doench_score(longSg)\n",
    "        doench_dict[\"DoenchScore\"] += [d_score]\n",
    "        sg_dict[\"SgID\"] += [sgID]\n",
    "        sg_dict[\"LongSg\"] += [longSg]\n",
    "    db_con.update_many_rows(doench_dict, sg_dict, \"SgRNATargetInformation\")\n",
    "    db_con.close_cursor()\n",
    "    db_con.close_connection()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "run_doench(\"miR-test\", firewall=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Add MaxDoenchScore\n",
    "\n",
    "Once the Doench score had been calculated for each extended sgRNA, the maximum Doench score for the sgRNA was determined."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "import data_processing as dp\n",
    "\n",
    "def find_max_doench(db_name, sql_version=\"MySQL\", firewall=False):\n",
    "    db_con = dp.DatabaseConnection(sql_version, db_name=db_name, firewall=firewall)\n",
    "    rows = db_con.fetch_query(\"\"\"SELECT s.SgID, t.DoenchScore \n",
    "FROM SingleGuideRNA AS s \n",
    "JOIN SgRNATargetInformation AS t \n",
    "ON s.SgID = t.SgID\"\"\")\n",
    "    max_dict = {}\n",
    "    for row in rows:\n",
    "        if sql_version == \"MSSQL\":\n",
    "            sg = row.SgID\n",
    "            score = row.DoenchScore\n",
    "        else:\n",
    "            sg, score = row\n",
    "        if sg in max_dict:\n",
    "            if score > max_dict[sg]:\n",
    "                max_dict[sg] = score\n",
    "            else:\n",
    "                pass\n",
    "        else:\n",
    "            max_dict[sg] = score\n",
    "    sg_dict = {\"SgID\": []}\n",
    "    do_dict = {\"MaxDoenchScore\": []}\n",
    "    for key, val in max_dict.iteritems():\n",
    "        sg_dict[\"SgID\"] += [key]\n",
    "        do_dict[\"MaxDoenchScore\"] += [val]\n",
    "        \n",
    "    db_con.update_many_rows(do_dict, sg_dict, \"SingleGuideRNA\")\n",
    "    db_con.close_cursor()\n",
    "    db_con.close_connection()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "find_max_doench(\"miR-test\", firewall=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Azimuth Score\n",
    "\n",
    "This on-target scoring method was updated in 2015 (<a href=\"https://doi.org/10.1101/021568\">Preprint</a>)/2016 (<a href=\"https://doi.org/10.1038/nbt.3437\">Nature Biotechnology Paper</a>) with the Azimuth score. I contacted Microsoft and recieved an API key which enabled me to access the API and submit requests for the scores of my sgRNAs. The documentation for this API is <a href=\"https://studio.azureml.net/apihelp/workspaces/ee5485c1d9814b8d8c647a89db12d4df/webservices/72e5e606de0b4fa0bcde57666f0ddcba/endpoints/c24d128abfaf4832abf1e7ef45db4b54/score\">here</a>. The source code is now up at <a href=\"https://github.com/MicrosoftResearch/Azimuth\">GitHub</a> and there is a <a href=\"http://www.broadinstitute.org/rnai/public/analysis-tools/sgrna-design\">web interface</a> as well. \n",
    "\n",
    "### From flat file\n",
    "\n",
    "It seems they tweaked the scoring over time. The Azimuth scores from the API are now slightly different than the original and the API is no longer being supported. The older score will therefore be loaded from a flat file."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "import data_processing as dp\n",
    "\n",
    "def import_azimuth_score(db_name, sql_version=\"MySQL\", firewall=False):\n",
    "    \"\"\"\n",
    "        Fetches Azimuth score from flat file\n",
    "    \"\"\"\n",
    "    load_dict = {\"SgID\": [], \"PriID\": [], \"SgStart\": []}\n",
    "    az_dict = {\"AzimuthScore\": []}\n",
    "    \n",
    "    df = pd.read_csv(\"sgRNA Scoring/SgRNATargetInformation_table.csv\", header=0, index_col=0)\n",
    "    for sgID, row in df.iterrows():\n",
    "        load_dict[\"SgID\"] += [sgID]\n",
    "        load_dict[\"PriID\"] += [row[\"PriID\"]]\n",
    "        load_dict[\"SgStart\"] += [int(row[\"SgStart\"])]\n",
    "        if pd.isnull(row[\"AzimuthScore\"]):\n",
    "            az_dict[\"AzimuthScore\"] += [None]\n",
    "        else:\n",
    "            az_dict[\"AzimuthScore\"] += [float(row[\"AzimuthScore\"])]\n",
    "    \n",
    "    db_con = dp.DatabaseConnection(sql_version, db_name=db_name, firewall=firewall)\n",
    "    db_con.update_many_rows(az_dict, load_dict, \"SgRNATargetInformation\")\n",
    "    db_con.close_cursor()\n",
    "    db_con.close_connection()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [],
   "source": [
    "import_azimuth_score(\"miR-test\", firewall=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### From Azimuth API\n",
    "\n",
    "The score can be fetched using the API and imported into a new column."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "import urllib2\n",
    "import math\n",
    "import time\n",
    "import json\n",
    "import data_processing as dp\n",
    "\n",
    "### Azimuth (updated Doench) Score ###\n",
    "\n",
    "def get_azimuth_score(out_file, api_key_file, db_name, sql_version=\"MySQL\", firewall=False):\n",
    "        \n",
    "    db_con = dp.DatabaseConnection(sql_version, db_name=db_name, firewall=firewall)\n",
    "    \n",
    "    # create new column\n",
    "    #db_con.add_column(\"AzimuthScorev2\", \"FLOAT\", \"SgRNATargetInformation\")\n",
    "\n",
    "    STEP = 1000\n",
    "    \n",
    "    rows = db_con.fetch_query(\"SELECT SgID, LongSg FROM SgRNATargetInformation\")\n",
    "    \n",
    "    sg_dict = {\"SgID\": [], \"LongSg\": []}\n",
    "    az_dict = {\"AzimuthScorev2\": []}\n",
    "    for row in rows:\n",
    "        if sql_version == \"MSSQL\":\n",
    "            sg_dict[\"SgID\"] += [row.SgID]\n",
    "            sg_dict[\"LongSg\"] += [row.LongSg]\n",
    "        else:\n",
    "            sgID, longSg = row\n",
    "            sg_dict[\"SgID\"] += [sgID]\n",
    "            sg_dict[\"LongSg\"] += [str(longSg)]\n",
    "    with open(out_file, \"a\") as fout:\n",
    "        fout.write(\"SgID,LongSg,AzimuthScorev2\\n\")\n",
    "    \n",
    "    num_sg = len(sg_dict[\"LongSg\"])\n",
    "    num_chunks = int(math.ceil(num_sg/float(STEP)))\n",
    "    for n in range(num_chunks):\n",
    "        longSg_start = n*STEP\n",
    "        longSg_end = min(n*STEP+STEP, num_sg)\n",
    "        longSg_chunk = sg_dict[\"LongSg\"][longSg_start:longSg_end]\n",
    "        \n",
    "        sgRNAs = [[sg.upper(), -1, -1] for sg in longSg_chunk]\n",
    "\n",
    "        data = {\"Inputs\": {\n",
    "                            \"input1\":{\n",
    "                                      \"ColumnNames\": [\"sequence\", \"cutsite\", \"percentpeptide\"],\n",
    "                                      \"Values\": sgRNAs\n",
    "                                     },        \n",
    "                          },\n",
    "                \"GlobalParameters\": {}\n",
    "               }\n",
    "\n",
    "        body = str.encode(json.dumps(data))\n",
    "\n",
    "        url = 'https://ussouthcentral.services.azureml.net/workspaces/ee5485c1d9814b8d8c647a89db12d4df/services/c24d128abfaf4832abf1e7ef45db4b54/execute?api-version=2.0&details=true'\n",
    "        with open(api_key_file, \"r\") as api:\n",
    "            api_key = api.readline()\n",
    "        headers = {'Content-Type':'application/json', 'Authorization':('Bearer '+ api_key)}\n",
    "\n",
    "        req = urllib2.Request(url, body, headers)\n",
    "\n",
    "        try:\n",
    "            response = urllib2.urlopen(req)\n",
    "\n",
    "            result = response.read()\n",
    "            parsed_results = json.loads(result)\n",
    "            listOfResults = parsed_results[\"Results\"][\"output2\"][\"value\"][\"Values\"]\n",
    "            with open(out_file, \"a\") as fout:\n",
    "                for m in range(len(listOfResults)):\n",
    "                    score = float(listOfResults[m][0])\n",
    "                    sgID = sg_dict[\"SgID\"][longSg_start+m]\n",
    "                    longSg = sg_dict[\"LongSg\"][longSg_start+m]\n",
    "                    fout.write(\"{},{},{}\\n\".format(sgID, longSg, score))\n",
    "                    az_dict[\"AzimuthScorev2\"] += [score]\n",
    "                    \n",
    "        except urllib2.HTTPError, error:\n",
    "            print(\"The request failed with status code: \" + str(error.code))\n",
    "\n",
    "            # Print the headers - they include the requert ID and the timestamp, which are useful for debugging the failure\n",
    "            print(error.info())\n",
    "\n",
    "            print(json.loads(error.read()))\n",
    "        time.sleep(30)\n",
    "    db_con.update_many_rows(az_dict, sg_dict, \"SgRNATargetInformation\")\n",
    "    db_con.close_cursor()\n",
    "    db_con.close_connection()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "get_azimuth_score(\"sgRNA Scoring/Azimuth v2 Scores.csv\", \"Azimuth_API_key.txt\", \"miR-test\", firewall=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Max Azimuth Score\n",
    "\n",
    "The maximum azimuth score for all target sites for each sgRNA can then be determined. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "import data_processing as dp\n",
    "\n",
    "def find_max_azimuth(db_name, sql_version=\"MySQL\", firewall=False):\n",
    "    db_con = dp.DatabaseConnection(sql_version, db_name=db_name, firewall=firewall)\n",
    "    \n",
    "    rows = db_con.fetch_query(\"\"\"SELECT s.SgID, t.AzimuthScore \n",
    "FROM SingleGuideRNA AS s \n",
    "JOIN SgRNATargetInformation AS t \n",
    "ON s.SgID = t.SgID\"\"\")\n",
    "    \n",
    "    max_dict = {}\n",
    "    for row in rows:\n",
    "        if sql_version == \"MSSQL\":\n",
    "            sg = row.SgID\n",
    "            score = row.AzimuthScore\n",
    "        else:\n",
    "            sg, score = row\n",
    "        if sg in max_dict:\n",
    "            if score > max_dict[sg]:\n",
    "                max_dict[sg] = score\n",
    "            else:\n",
    "                pass\n",
    "        else:\n",
    "            max_dict[sg] = score\n",
    "    sg_dict = {\"SgID\": []}\n",
    "    do_dict = {\"MaxAzimuthScore\": []}\n",
    "    for key, val in max_dict.iteritems():\n",
    "        sg_dict[\"SgID\"] += [key]\n",
    "        do_dict[\"MaxAzimuthScore\"] += [val]\n",
    "        \n",
    "    db_con.update_many_rows(do_dict, sg_dict, \"SingleGuideRNA\")\n",
    "    db_con.close_cursor()\n",
    "    db_con.close_connection()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "find_max_azimuth(\"miR-test\", firewall=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## sgRNA Scorer\n",
    "\n",
    "The Church lab also came up with an on-target activity scoring algorthm call <a href=\"https://crispr.med.harvard.edu/\">sgRNA Scorer</a>. <a href=\"http://www.nature.com/nmeth/journal/v12/n9/full/nmeth.3473.html\">This paper</a> describes how this method was developed. The input for the scoring is a fasta file with the sgRNA + PAM sequence. This information can be sliced from the longSg sequence using the below code:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import data_processing as dp\n",
    "\n",
    "def get_sgRNAs(out_file, db_name, sql_version=\"MySQL\", firewall=False):\n",
    "    \"\"\"\n",
    "        Fetches the sgRNA and PAM sequences from the database for sgRNA Scorer\n",
    "    \"\"\"\n",
    "    db_con = dp.DatabaseConnection(sql_version, db_name=db_name, firewall=firewall)\n",
    "    rows = db_con.fetch_query(\"SELECT SgID, LongSg FROM SgRNATargetInformation\")\n",
    "    db_con.close_cursor()\n",
    "    db_con.close_connection()\n",
    "    \n",
    "    with open(out_file, \"w\") as fout:\n",
    "        for row in rows:\n",
    "            if sql_version == \"MSSQL\":\n",
    "                sgID = row.SgID\n",
    "                longSg = row.LongSg\n",
    "            else:\n",
    "                sgID, longSg = row\n",
    "                longSg = str(longSg)\n",
    "            header_str = \">{}_{}\\n\".format(sgID, longSg)\n",
    "            fout.write(header_str)\n",
    "            fout.write(\"{}\\n\".format(longSg[4:-3]))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "get_sgRNAs(\"Downloaded Programs/sgRNAScorer1_0/miR_sgRNAs.fa\", \"miR-test\", firewall=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The stand alone code can be downloaded from this <a href=\"https://crispr.med.harvard.edu/\">website</a> and the <a href=\"http://svmlight.joachims.org/\">svm-light binaries</a> are also necessary . The sgRNAs were then scored relative to the human exome sgRNAs included with the package."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "C:\\Users\\jkurata\\Documents\\GitHub\\LX-miR_Design\\Downloaded Programs\\sgRNAScorer1_0\n",
      "Total number of gRNAs: 26376\n",
      "Finished 0 sequences\n",
      "Time: Wed Feb 14 16:41:14 2018\n",
      "Finished 100 sequences\n",
      "Time: Wed Feb 14 16:42:00 2018\n",
      "Finished 200 sequences\n",
      "Time: Wed Feb 14 16:42:45 2018\n",
      "Finished 300 sequences\n",
      "Time: Wed Feb 14 16:43:31 2018\n",
      "Finished 400 sequences\n",
      "Time: Wed Feb 14 16:44:17 2018\n",
      "Finished 500 sequences\n",
      "Time: Wed Feb 14 16:45:01 2018\n",
      "Finished 600 sequences\n",
      "Time: Wed Feb 14 16:45:46 2018\n",
      "Finished 700 sequences\n",
      "Time: Wed Feb 14 16:46:31 2018\n",
      "Finished 800 sequences\n",
      "Time: Wed Feb 14 16:47:16 2018\n",
      "Finished 900 sequences\n",
      "Time: Wed Feb 14 16:48:00 2018\n",
      "Finished 1000 sequences\n",
      "Time: Wed Feb 14 16:48:45 2018\n",
      "Finished 1100 sequences\n",
      "Time: Wed Feb 14 16:49:30 2018\n",
      "Finished 1200 sequences\n",
      "Time: Wed Feb 14 16:50:15 2018\n",
      "Finished 1300 sequences\n",
      "Time: Wed Feb 14 16:51:00 2018\n",
      "Finished 1400 sequences\n",
      "Time: Wed Feb 14 16:51:45 2018\n",
      "Finished 1500 sequences\n",
      "Time: Wed Feb 14 16:52:30 2018\n",
      "Finished 1600 sequences\n",
      "Time: Wed Feb 14 16:53:15 2018\n",
      "Finished 1700 sequences\n",
      "Time: Wed Feb 14 16:54:00 2018\n",
      "Finished 1800 sequences\n",
      "Time: Wed Feb 14 16:54:45 2018\n",
      "Finished 1900 sequences\n",
      "Time: Wed Feb 14 16:55:30 2018\n",
      "Finished 2000 sequences\n",
      "Time: Wed Feb 14 16:56:15 2018\n",
      "Finished 2100 sequences\n",
      "Time: Wed Feb 14 16:57:00 2018\n",
      "Finished 2200 sequences\n",
      "Time: Wed Feb 14 16:57:45 2018\n",
      "Finished 2300 sequences\n",
      "Time: Wed Feb 14 16:58:31 2018\n",
      "Finished 2400 sequences\n",
      "Time: Wed Feb 14 16:59:16 2018\n",
      "Finished 2500 sequences\n",
      "Time: Wed Feb 14 17:00:01 2018\n",
      "Finished 2600 sequences\n",
      "Time: Wed Feb 14 17:00:46 2018\n",
      "Finished 2700 sequences\n",
      "Time: Wed Feb 14 17:01:31 2018\n",
      "Finished 2800 sequences\n",
      "Time: Wed Feb 14 17:02:16 2018\n",
      "Finished 2900 sequences\n",
      "Time: Wed Feb 14 17:03:01 2018\n",
      "Finished 3000 sequences\n",
      "Time: Wed Feb 14 17:03:46 2018\n",
      "Finished 3100 sequences\n",
      "Time: Wed Feb 14 17:04:30 2018\n",
      "Finished 3200 sequences\n",
      "Time: Wed Feb 14 17:05:15 2018\n",
      "Finished 3300 sequences\n",
      "Time: Wed Feb 14 17:06:00 2018\n",
      "Finished 3400 sequences\n",
      "Time: Wed Feb 14 17:06:45 2018\n",
      "Finished 3500 sequences\n",
      "Time: Wed Feb 14 17:07:30 2018\n",
      "Finished 3600 sequences\n",
      "Time: Wed Feb 14 17:08:15 2018\n",
      "Finished 3700 sequences\n",
      "Time: Wed Feb 14 17:08:59 2018\n",
      "Finished 3800 sequences\n",
      "Time: Wed Feb 14 17:09:44 2018\n",
      "Finished 3900 sequences\n",
      "Time: Wed Feb 14 17:10:31 2018\n",
      "Finished 4000 sequences\n",
      "Time: Wed Feb 14 17:11:16 2018\n",
      "Finished 4100 sequences\n",
      "Time: Wed Feb 14 17:12:00 2018\n",
      "Finished 4200 sequences\n",
      "Time: Wed Feb 14 17:12:45 2018\n",
      "Finished 4300 sequences\n",
      "Time: Wed Feb 14 17:13:30 2018\n",
      "Finished 4400 sequences\n",
      "Time: Wed Feb 14 17:14:15 2018\n",
      "Finished 4500 sequences\n",
      "Time: Wed Feb 14 17:15:00 2018\n",
      "Finished 4600 sequences\n",
      "Time: Wed Feb 14 17:15:45 2018\n",
      "Finished 4700 sequences\n",
      "Time: Wed Feb 14 17:16:30 2018\n",
      "Finished 4800 sequences\n",
      "Time: Wed Feb 14 17:17:16 2018\n",
      "Finished 4900 sequences\n",
      "Time: Wed Feb 14 17:18:02 2018\n",
      "Finished 5000 sequences\n",
      "Time: Wed Feb 14 17:18:48 2018\n",
      "Finished 5100 sequences\n",
      "Time: Wed Feb 14 17:19:34 2018\n",
      "Finished 5200 sequences\n",
      "Time: Wed Feb 14 17:20:19 2018\n",
      "Finished 5300 sequences\n",
      "Time: Wed Feb 14 17:21:05 2018\n",
      "Finished 5400 sequences\n",
      "Time: Wed Feb 14 17:21:50 2018\n",
      "Finished 5500 sequences\n",
      "Time: Wed Feb 14 17:22:35 2018\n",
      "Finished 5600 sequences\n",
      "Time: Wed Feb 14 17:23:22 2018\n",
      "Finished 5700 sequences\n",
      "Time: Wed Feb 14 17:24:07 2018\n",
      "Finished 5800 sequences\n",
      "Time: Wed Feb 14 17:24:52 2018\n",
      "Finished 5900 sequences\n",
      "Time: Wed Feb 14 17:25:37 2018\n",
      "Finished 6000 sequences\n",
      "Time: Wed Feb 14 17:26:22 2018\n",
      "Finished 6100 sequences\n",
      "Time: Wed Feb 14 17:27:07 2018\n",
      "Finished 6200 sequences\n",
      "Time: Wed Feb 14 17:27:52 2018\n",
      "Finished 6300 sequences\n",
      "Time: Wed Feb 14 17:28:37 2018\n",
      "Finished 6400 sequences\n",
      "Time: Wed Feb 14 17:29:22 2018\n",
      "Finished 6500 sequences\n",
      "Time: Wed Feb 14 17:30:08 2018\n",
      "Finished 6600 sequences\n",
      "Time: Wed Feb 14 17:30:54 2018\n",
      "Finished 6700 sequences\n",
      "Time: Wed Feb 14 17:31:39 2018\n",
      "Finished 6800 sequences\n",
      "Time: Wed Feb 14 17:32:26 2018\n",
      "Finished 6900 sequences\n",
      "Time: Wed Feb 14 17:33:12 2018\n",
      "Finished 7000 sequences\n",
      "Time: Wed Feb 14 17:33:57 2018\n",
      "Finished 7100 sequences\n",
      "Time: Wed Feb 14 17:34:43 2018\n",
      "Finished 7200 sequences\n",
      "Time: Wed Feb 14 17:35:29 2018\n",
      "Finished 7300 sequences\n",
      "Time: Wed Feb 14 17:36:15 2018\n",
      "Finished 7400 sequences\n",
      "Time: Wed Feb 14 17:37:02 2018\n",
      "Finished 7500 sequences\n",
      "Time: Wed Feb 14 17:37:47 2018\n",
      "Finished 7600 sequences\n",
      "Time: Wed Feb 14 17:38:33 2018\n",
      "Finished 7700 sequences\n",
      "Time: Wed Feb 14 17:39:18 2018\n",
      "Finished 7800 sequences\n",
      "Time: Wed Feb 14 17:40:03 2018\n",
      "Finished 7900 sequences\n",
      "Time: Wed Feb 14 17:40:48 2018\n",
      "Finished 8000 sequences\n",
      "Time: Wed Feb 14 17:41:34 2018\n",
      "Finished 8100 sequences\n",
      "Time: Wed Feb 14 17:42:19 2018\n",
      "Finished 8200 sequences\n",
      "Time: Wed Feb 14 17:43:04 2018\n",
      "Finished 8300 sequences\n",
      "Time: Wed Feb 14 17:43:49 2018\n",
      "Finished 8400 sequences\n",
      "Time: Wed Feb 14 17:44:34 2018\n",
      "Finished 8500 sequences\n",
      "Time: Wed Feb 14 17:45:19 2018\n",
      "Finished 8600 sequences\n",
      "Time: Wed Feb 14 17:46:04 2018\n",
      "Finished 8700 sequences\n",
      "Time: Wed Feb 14 17:46:49 2018\n",
      "Finished 8800 sequences\n",
      "Time: Wed Feb 14 17:47:34 2018\n",
      "Finished 8900 sequences\n",
      "Time: Wed Feb 14 17:48:19 2018\n",
      "Finished 9000 sequences\n",
      "Time: Wed Feb 14 17:49:04 2018\n",
      "Finished 9100 sequences\n",
      "Time: Wed Feb 14 17:49:49 2018\n",
      "Finished 9200 sequences\n",
      "Time: Wed Feb 14 17:50:34 2018\n",
      "Finished 9300 sequences\n",
      "Time: Wed Feb 14 17:51:19 2018\n",
      "Finished 9400 sequences\n",
      "Time: Wed Feb 14 17:52:03 2018\n",
      "Finished 9500 sequences\n",
      "Time: Wed Feb 14 17:52:48 2018\n",
      "Finished 9600 sequences\n",
      "Time: Wed Feb 14 17:53:33 2018\n",
      "Finished 9700 sequences\n",
      "Time: Wed Feb 14 17:54:18 2018\n",
      "Finished 9800 sequences\n",
      "Time: Wed Feb 14 17:55:03 2018\n",
      "Finished 9900 sequences\n",
      "Time: Wed Feb 14 17:55:48 2018\n",
      "Finished 10000 sequences\n",
      "Time: Wed Feb 14 17:56:33 2018\n",
      "Finished 10100 sequences\n",
      "Time: Wed Feb 14 17:57:17 2018\n",
      "Finished 10200 sequences\n",
      "Time: Wed Feb 14 17:58:02 2018\n",
      "Finished 10300 sequences\n",
      "Time: Wed Feb 14 17:58:47 2018\n",
      "Finished 10400 sequences\n",
      "Time: Wed Feb 14 17:59:31 2018\n",
      "Finished 10500 sequences\n",
      "Time: Wed Feb 14 18:00:16 2018\n",
      "Finished 10600 sequences\n",
      "Time: Wed Feb 14 18:01:01 2018\n",
      "Finished 10700 sequences\n",
      "Time: Wed Feb 14 18:01:46 2018\n",
      "Finished 10800 sequences\n",
      "Time: Wed Feb 14 18:02:31 2018\n",
      "Finished 10900 sequences\n",
      "Time: Wed Feb 14 18:03:16 2018\n",
      "Finished 11000 sequences\n",
      "Time: Wed Feb 14 18:04:01 2018\n",
      "Finished 11100 sequences\n",
      "Time: Wed Feb 14 18:04:47 2018\n",
      "Finished 11200 sequences\n",
      "Time: Wed Feb 14 18:05:32 2018\n",
      "Finished 11300 sequences\n",
      "Time: Wed Feb 14 18:06:17 2018\n",
      "Finished 11400 sequences\n",
      "Time: Wed Feb 14 18:07:02 2018\n",
      "Finished 11500 sequences\n",
      "Time: Wed Feb 14 18:07:47 2018\n",
      "Finished 11600 sequences\n",
      "Time: Wed Feb 14 18:08:32 2018\n",
      "Finished 11700 sequences\n",
      "Time: Wed Feb 14 18:09:17 2018\n",
      "Finished 11800 sequences\n",
      "Time: Wed Feb 14 18:10:13 2018\n",
      "Finished 11900 sequences\n",
      "Time: Wed Feb 14 18:11:09 2018\n",
      "Finished 12000 sequences\n",
      "Time: Wed Feb 14 18:12:04 2018\n",
      "Finished 12100 sequences\n",
      "Time: Wed Feb 14 18:12:59 2018\n",
      "Finished 12200 sequences\n",
      "Time: Wed Feb 14 18:13:55 2018\n",
      "Finished 12300 sequences\n",
      "Time: Wed Feb 14 18:14:51 2018\n",
      "Finished 12400 sequences\n",
      "Time: Wed Feb 14 18:15:47 2018\n",
      "Finished 12500 sequences\n",
      "Time: Wed Feb 14 18:16:42 2018\n",
      "Finished 12600 sequences\n",
      "Time: Wed Feb 14 18:17:38 2018\n",
      "Finished 12700 sequences\n",
      "Time: Wed Feb 14 18:18:34 2018\n",
      "Finished 12800 sequences\n",
      "Time: Wed Feb 14 18:19:30 2018\n",
      "Finished 12900 sequences\n",
      "Time: Wed Feb 14 18:20:25 2018\n",
      "Finished 13000 sequences\n",
      "Time: Wed Feb 14 18:21:21 2018\n",
      "Finished 13100 sequences\n",
      "Time: Wed Feb 14 18:22:16 2018\n",
      "Finished 13200 sequences\n",
      "Time: Wed Feb 14 18:23:12 2018\n",
      "Finished 13300 sequences\n",
      "Time: Wed Feb 14 18:24:07 2018\n",
      "Finished 13400 sequences\n",
      "Time: Wed Feb 14 18:25:03 2018\n",
      "Finished 13500 sequences\n",
      "Time: Wed Feb 14 18:25:59 2018\n",
      "Finished 13600 sequences\n",
      "Time: Wed Feb 14 18:26:54 2018\n",
      "Finished 13700 sequences\n",
      "Time: Wed Feb 14 18:27:50 2018\n",
      "Finished 13800 sequences\n",
      "Time: Wed Feb 14 18:28:46 2018\n",
      "Finished 13900 sequences\n",
      "Time: Wed Feb 14 18:29:42 2018\n",
      "Finished 14000 sequences\n",
      "Time: Wed Feb 14 18:30:38 2018\n",
      "Finished 14100 sequences\n",
      "Time: Wed Feb 14 18:31:34 2018\n",
      "Finished 14200 sequences\n",
      "Time: Wed Feb 14 18:32:30 2018\n",
      "Finished 14300 sequences\n",
      "Time: Wed Feb 14 18:33:26 2018\n",
      "Finished 14400 sequences\n",
      "Time: Wed Feb 14 18:34:22 2018\n",
      "Finished 14500 sequences\n",
      "Time: Wed Feb 14 18:35:18 2018\n",
      "Finished 14600 sequences\n",
      "Time: Wed Feb 14 18:36:14 2018\n",
      "Finished 14700 sequences\n",
      "Time: Wed Feb 14 18:37:10 2018\n",
      "Finished 14800 sequences\n",
      "Time: Wed Feb 14 18:38:06 2018\n",
      "Finished 14900 sequences\n",
      "Time: Wed Feb 14 18:39:01 2018\n",
      "Finished 15000 sequences\n",
      "Time: Wed Feb 14 18:39:57 2018\n",
      "Finished 15100 sequences\n",
      "Time: Wed Feb 14 18:40:53 2018\n",
      "Finished 15200 sequences\n",
      "Time: Wed Feb 14 18:41:48 2018\n",
      "Finished 15300 sequences\n",
      "Time: Wed Feb 14 18:42:43 2018\n",
      "Finished 15400 sequences\n",
      "Time: Wed Feb 14 18:43:38 2018\n",
      "Finished 15500 sequences\n",
      "Time: Wed Feb 14 18:44:33 2018\n",
      "Finished 15600 sequences\n",
      "Time: Wed Feb 14 18:45:29 2018\n",
      "Finished 15700 sequences\n",
      "Time: Wed Feb 14 18:46:24 2018\n",
      "Finished 15800 sequences\n",
      "Time: Wed Feb 14 18:47:20 2018\n",
      "Finished 15900 sequences\n",
      "Time: Wed Feb 14 18:48:16 2018\n",
      "Finished 16000 sequences\n",
      "Time: Wed Feb 14 18:49:12 2018\n",
      "Finished 16100 sequences\n",
      "Time: Wed Feb 14 18:50:08 2018\n",
      "Finished 16200 sequences\n",
      "Time: Wed Feb 14 18:51:04 2018\n",
      "Finished 16300 sequences\n",
      "Time: Wed Feb 14 18:52:00 2018\n",
      "Finished 16400 sequences\n",
      "Time: Wed Feb 14 18:52:56 2018\n",
      "Finished 16500 sequences\n",
      "Time: Wed Feb 14 18:53:52 2018\n",
      "Finished 16600 sequences\n",
      "Time: Wed Feb 14 18:54:48 2018\n",
      "Finished 16700 sequences\n",
      "Time: Wed Feb 14 18:55:43 2018\n",
      "Finished 16800 sequences\n",
      "Time: Wed Feb 14 18:56:38 2018\n",
      "Finished 16900 sequences\n",
      "Time: Wed Feb 14 18:57:34 2018\n",
      "Finished 17000 sequences\n",
      "Time: Wed Feb 14 18:58:29 2018\n",
      "Finished 17100 sequences\n",
      "Time: Wed Feb 14 18:59:24 2018\n",
      "Finished 17200 sequences\n",
      "Time: Wed Feb 14 19:00:19 2018\n",
      "Finished 17300 sequences\n",
      "Time: Wed Feb 14 19:01:15 2018\n",
      "Finished 17400 sequences\n",
      "Time: Wed Feb 14 19:02:10 2018\n",
      "Finished 17500 sequences\n",
      "Time: Wed Feb 14 19:03:05 2018\n",
      "Finished 17600 sequences\n",
      "Time: Wed Feb 14 19:04:01 2018\n",
      "Finished 17700 sequences\n",
      "Time: Wed Feb 14 19:04:57 2018\n",
      "Finished 17800 sequences\n",
      "Time: Wed Feb 14 19:05:53 2018\n",
      "Finished 17900 sequences\n",
      "Time: Wed Feb 14 19:06:48 2018\n",
      "Finished 18000 sequences\n",
      "Time: Wed Feb 14 19:07:43 2018\n",
      "Finished 18100 sequences\n",
      "Time: Wed Feb 14 19:08:39 2018\n",
      "Finished 18200 sequences\n",
      "Time: Wed Feb 14 19:09:35 2018\n",
      "Finished 18300 sequences\n",
      "Time: Wed Feb 14 19:10:31 2018\n",
      "Finished 18400 sequences\n",
      "Time: Wed Feb 14 19:11:27 2018\n",
      "Finished 18500 sequences\n",
      "Time: Wed Feb 14 19:12:22 2018\n",
      "Finished 18600 sequences\n",
      "Time: Wed Feb 14 19:13:18 2018\n",
      "Finished 18700 sequences\n",
      "Time: Wed Feb 14 19:14:14 2018\n",
      "Finished 18800 sequences\n",
      "Time: Wed Feb 14 19:15:10 2018\n",
      "Finished 18900 sequences\n",
      "Time: Wed Feb 14 19:16:05 2018\n",
      "Finished 19000 sequences\n",
      "Time: Wed Feb 14 19:17:01 2018\n",
      "Finished 19100 sequences\n",
      "Time: Wed Feb 14 19:17:56 2018\n",
      "Finished 19200 sequences\n",
      "Time: Wed Feb 14 19:18:52 2018\n",
      "Finished 19300 sequences\n",
      "Time: Wed Feb 14 19:19:47 2018\n",
      "Finished 19400 sequences\n",
      "Time: Wed Feb 14 19:20:43 2018\n",
      "Finished 19500 sequences\n",
      "Time: Wed Feb 14 19:21:38 2018\n",
      "Finished 19600 sequences\n",
      "Time: Wed Feb 14 19:22:34 2018\n",
      "Finished 19700 sequences\n",
      "Time: Wed Feb 14 19:23:30 2018\n",
      "Finished 19800 sequences\n",
      "Time: Wed Feb 14 19:24:25 2018\n",
      "Finished 19900 sequences\n",
      "Time: Wed Feb 14 19:25:20 2018\n",
      "Finished 20000 sequences\n",
      "Time: Wed Feb 14 19:26:16 2018\n",
      "Finished 20100 sequences\n",
      "Time: Wed Feb 14 19:27:11 2018\n",
      "Finished 20200 sequences\n",
      "Time: Wed Feb 14 19:28:07 2018\n",
      "Finished 20300 sequences\n",
      "Time: Wed Feb 14 19:29:02 2018\n",
      "Finished 20400 sequences\n",
      "Time: Wed Feb 14 19:29:58 2018\n",
      "Finished 20500 sequences\n",
      "Time: Wed Feb 14 19:30:54 2018\n",
      "Finished 20600 sequences\n",
      "Time: Wed Feb 14 19:31:49 2018\n",
      "Finished 20700 sequences\n",
      "Time: Wed Feb 14 19:32:45 2018\n",
      "Finished 20800 sequences\n",
      "Time: Wed Feb 14 19:33:41 2018\n",
      "Finished 20900 sequences\n",
      "Time: Wed Feb 14 19:34:37 2018\n",
      "Finished 21000 sequences\n",
      "Time: Wed Feb 14 19:35:32 2018\n",
      "Finished 21100 sequences\n",
      "Time: Wed Feb 14 19:36:28 2018\n",
      "Finished 21200 sequences\n",
      "Time: Wed Feb 14 19:37:24 2018\n",
      "Finished 21300 sequences\n",
      "Time: Wed Feb 14 19:38:19 2018\n",
      "Finished 21400 sequences\n",
      "Time: Wed Feb 14 19:39:15 2018\n",
      "Finished 21500 sequences\n",
      "Time: Wed Feb 14 19:40:11 2018\n",
      "Finished 21600 sequences\n",
      "Time: Wed Feb 14 19:41:07 2018\n",
      "Finished 21700 sequences\n",
      "Time: Wed Feb 14 19:42:02 2018\n",
      "Finished 21800 sequences\n",
      "Time: Wed Feb 14 19:42:58 2018\n",
      "Finished 21900 sequences\n",
      "Time: Wed Feb 14 19:43:53 2018\n",
      "Finished 22000 sequences\n",
      "Time: Wed Feb 14 19:44:49 2018\n",
      "Finished 22100 sequences\n",
      "Time: Wed Feb 14 19:45:44 2018\n",
      "Finished 22200 sequences\n",
      "Time: Wed Feb 14 19:46:40 2018\n",
      "Finished 22300 sequences\n",
      "Time: Wed Feb 14 19:47:35 2018\n",
      "Finished 22400 sequences\n",
      "Time: Wed Feb 14 19:48:31 2018\n",
      "Finished 22500 sequences\n",
      "Time: Wed Feb 14 19:49:26 2018\n",
      "Finished 22600 sequences\n",
      "Time: Wed Feb 14 19:50:21 2018\n",
      "Finished 22700 sequences\n",
      "Time: Wed Feb 14 19:51:17 2018\n",
      "Finished 22800 sequences\n",
      "Time: Wed Feb 14 19:52:12 2018\n",
      "Finished 22900 sequences\n",
      "Time: Wed Feb 14 19:53:07 2018\n",
      "Finished 23000 sequences\n",
      "Time: Wed Feb 14 19:54:03 2018\n",
      "Finished 23100 sequences\n",
      "Time: Wed Feb 14 19:54:59 2018\n",
      "Finished 23200 sequences\n",
      "Time: Wed Feb 14 19:55:55 2018\n",
      "Finished 23300 sequences\n",
      "Time: Wed Feb 14 19:56:50 2018\n",
      "Finished 23400 sequences\n",
      "Time: Wed Feb 14 19:57:46 2018\n",
      "Finished 23500 sequences\n",
      "Time: Wed Feb 14 19:58:42 2018\n",
      "Finished 23600 sequences\n",
      "Time: Wed Feb 14 19:59:38 2018\n",
      "Finished 23700 sequences\n",
      "Time: Wed Feb 14 20:00:34 2018\n",
      "Finished 23800 sequences\n",
      "Time: Wed Feb 14 20:01:29 2018\n",
      "Finished 23900 sequences\n",
      "Time: Wed Feb 14 20:02:26 2018\n",
      "Finished 24000 sequences\n",
      "Time: Wed Feb 14 20:03:22 2018\n",
      "Finished 24100 sequences\n",
      "Time: Wed Feb 14 20:04:17 2018\n",
      "Finished 24200 sequences\n",
      "Time: Wed Feb 14 20:05:12 2018\n",
      "Finished 24300 sequences\n",
      "Time: Wed Feb 14 20:06:07 2018\n",
      "Finished 24400 sequences\n",
      "Time: Wed Feb 14 20:07:02 2018\n",
      "Finished 24500 sequences\n",
      "Time: Wed Feb 14 20:07:58 2018\n",
      "Finished 24600 sequences\n",
      "Time: Wed Feb 14 20:08:54 2018\n",
      "Finished 24700 sequences\n",
      "Time: Wed Feb 14 20:09:49 2018\n",
      "Finished 24800 sequences\n",
      "Time: Wed Feb 14 20:10:45 2018\n",
      "Finished 24900 sequences\n",
      "Time: Wed Feb 14 20:11:41 2018\n",
      "Finished 25000 sequences\n",
      "Time: Wed Feb 14 20:12:36 2018\n",
      "Finished 25100 sequences\n",
      "Time: Wed Feb 14 20:13:32 2018\n",
      "Finished 25200 sequences\n",
      "Time: Wed Feb 14 20:14:28 2018\n",
      "Finished 25300 sequences\n",
      "Time: Wed Feb 14 20:15:24 2018\n",
      "Finished 25400 sequences\n",
      "Time: Wed Feb 14 20:16:20 2018\n",
      "Finished 25500 sequences\n",
      "Time: Wed Feb 14 20:17:16 2018\n",
      "Finished 25600 sequences\n",
      "Time: Wed Feb 14 20:18:11 2018\n",
      "Finished 25700 sequences\n",
      "Time: Wed Feb 14 20:19:07 2018\n",
      "Finished 25800 sequences\n",
      "Time: Wed Feb 14 20:20:03 2018\n",
      "Finished 25900 sequences\n",
      "Time: Wed Feb 14 20:20:59 2018\n",
      "Finished 26000 sequences\n",
      "Time: Wed Feb 14 20:21:55 2018\n",
      "Finished 26100 sequences\n",
      "Time: Wed Feb 14 20:22:50 2018\n",
      "Finished 26200 sequences\n",
      "Time: Wed Feb 14 20:23:46 2018\n",
      "Finished 26300 sequences\n",
      "Time: Wed Feb 14 20:24:41 2018\n",
      "Time: Wed Feb 14 16:41:04 2018\n",
      "Generating putative gRNA sites...\n",
      "Time: Wed Feb 14 16:41:05 2018\n",
      "Generating SVM input file from gRNA sequences...\n",
      "Time: Wed Feb 14 16:41:08 2018\n",
      "Running classification using SVM-Light\n",
      "Time: Wed Feb 14 16:41:10 2018\n",
      "Converting scores to ranks based on global SP score distribution\n"
     ]
    }
   ],
   "source": [
    "%cd \"Downloaded Programs/sgRNAScorer1_0\"\n",
    "!python scoreMySites.py miR_sgRNAs.fa Hg SP miR_sgRNAScorer\n",
    "%cd \"../..\""
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The scores for each sgRNA were then imported into the SingleGuideRNA table. For sgRNAs which had multiple possible PAM sites, the highest score is added to the table."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import data_processing as dp\n",
    "\n",
    "def import_sgRNAScorer(in_file, db_name, sql_version=\"MySQL\", firewall=False):\n",
    "    \"\"\"\n",
    "        Imports sgRNA Scorer score\n",
    "    \"\"\"\n",
    "    score_dict = {}\n",
    "    with open(in_file, \"r\") as fin:\n",
    "        # skip header line\n",
    "        fin.next()\n",
    "        for line in fin:\n",
    "            ele = line.strip(\"\\n\").split(\"\\t\")\n",
    "            sgID = int(ele[0].split(\"_\")[0])\n",
    "            if sgID not in score_dict:\n",
    "                score_dict[sgID] = float(ele[2])\n",
    "            else:\n",
    "                if score_dict[sgID] < float(ele[2]):\n",
    "                    score_dict[sgID] = float(ele[2])\n",
    "                else:\n",
    "                    continue\n",
    "    sg_dict = {\"SgID\": []}\n",
    "    sg_scorer = {\"SgRNAScorer\": []}\n",
    "    for key, val in score_dict.iteritems():\n",
    "        sg_dict[\"SgID\"] += [key]\n",
    "        sg_scorer[\"SgRNAScorer\"] += [val]\n",
    "    db_con = dp.DatabaseConnection(sql_version, db_name=db_name, firewall=firewall)\n",
    "    db_con.update_many_rows(sg_scorer, sg_dict, \"SingleGuideRNA\")\n",
    "    db_con.close_cursor()\n",
    "    db_con.close_connection()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "import_sgRNAScorer(\"Downloaded Programs/sgRNAScorer1_0/miR_sgRNAScorer.FinalOutput.txt\",\n",
    "                   \"miR-test\", firewall=True)"
   ]
  }
 ],
 "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.9"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}