5.2 Selecting Control sgRNAs.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# SgRNA Controls\n",
    "\n",
    "We are also including 1,000 control sgRNAs which do not target any sequences in the human genome. \n",
    "\n",
    "## Extracting from Literature\n",
    "\n",
    "If we need to test the efficacy of the knockout, we may want to use EGFP or puromycin. Therefore, we included 4 control sgRNAs targeting each. The EGFP targeting control sgRNAs were taken from the supliment of <a href=\"https:\\\\doi.org\\10.1126\\science.1247005\">Shalem et al., 2014</a> (1, 2, 5 and 6): \n",
    "\n",
    "<img src=\"Published Libraries/Shalem et al EGFP sgRNA.png\">\n",
    "\n",
    "\n",
    "The puromycin control sgRNAs were CGGCGTCTCGCCCGACCACC, GTCGGGCGAGACGCCGACGG, CGTGGTCCAGACCGCCACCG and ACGCGCGTCGGGCTCGACAT.\n",
    "\n",
    "We wanted the rest of the controls to not target any sequence in the human genome. These other control sequences were taken from <a href=\"https:\\\\doi.org\\10.1126\\science.1246981\">Wang et al., 2014</a> (pLX-sgRNA), <a href=\"https:\\\\doi.org\\10.1126\\science.aac7041\">Wang et al., 2015</a> (Essential-gene) and <a href=\"https:\\\\doi.org\\10.1038\\nmeth.3047\">Sanjana et al., 2014</a> (Human-GeCKOv2; library sequences downloaded from <a href=\"http://genome-engineering.org/gecko/?page_id=15\">here</a>). The control sgRNA sequences were extracted from the list of sgRNAs in the pool and compiled into a single file with out replicating any sgRNAs:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "\n",
    "with open(\"Control sgRNAs/All Literature Control sgRNAs.csv\", \"w\") as fout:\n",
    "    fout.write(\"sgID,sgRNA Sequence,Library,Lab\\n\")\n",
    "    fout.write(\"EGFP_1,GGGCGAGGAGCTGTTCACCG,EGFP,Lander\\n\")\n",
    "    fout.write(\"EGFP_2,GAGCTGGACGGCGACGTAAA,EGFP,Lander\\n\")\n",
    "    fout.write(\"EGFP_5,GAAGTTCGAGGGCGACACCC,EGFP,Lander\\n\")\n",
    "    fout.write(\"EGFP_6,GGTGAACCGCATCGAGCTGA,EGFP,Lander\\n\")\n",
    "    fout.write(\"Puro_1,CGGCGTCTCGCCCGACCACC,Puro,\\n\")\n",
    "    fout.write(\"Puro_2,GTCGGGCGAGACGCCGACGG,Puro,\\n\")\n",
    "    fout.write(\"Puro_3,CGTGGTCCAGACCGCCACCG,Puro,\\n\")\n",
    "    fout.write(\"Puro_4,ACGCGCGTCGGGCTCGACAT,Puro,\\n\")\n",
    "    \n",
    "    sg_list = [\"GGGCGAGGAGCTGTTCACCG\", \"GAGCTGGACGGCGACGTAAA\", \"GAAGTTCGAGGGCGACACCC\", \"GGTGAACCGCATCGAGCTGA\", \n",
    "               \"CGGCGTCTCGCCCGACCACC\", \"GTCGGGCGAGACGCCGACGG\", \"CGTGGTCCAGACCGCCACCG\", \"ACGCGCGTCGGGCTCGACAT\"]\n",
    "    \n",
    "    # Look at GeCKOv2 A (controls are identical in library A and B)\n",
    "    # due to use of \\r as line break, the entire file is interpreted as a single line by python, so use pd\n",
    "    df = pd.read_csv(\"Published Libraries/Human_GeCKOv2_Library_A_09Mar2015.csv\", header=0)\n",
    "    df_con = df[df[\"gene_id\"].str.contains(\"Control\")]\n",
    "    for i in df_con.index.tolist():\n",
    "        sgID = df.loc[i, \"UID\"]\n",
    "        seq = df.loc[i, \"seq\"]\n",
    "        if seq not in sg_list:\n",
    "            sg_list += [seq]\n",
    "            out_str = \"{},{},Human_GeCKOv2_Library,Zhang\\n\".format(sgID, seq)\n",
    "            fout.write(out_str)\n",
    "        else:\n",
    "            print \"SgRNA ID {} with sequence {} already in library\".format(sgID, seq)\n",
    "            \n",
    "    # Look at pLX-sgRNA\n",
    "    with open(\"Published Libraries/Wang et al 2014 Supplementary Table 1 sgRNA Sequences.csv\", \"r\") as pLX_f:\n",
    "        for line in pLX_f:\n",
    "            if \"control\" not in line:\n",
    "                continue\n",
    "            ele = line.split(\",\")\n",
    "            sgID = ele[0]\n",
    "            seq = ele[6]\n",
    "            if seq not in sg_list:\n",
    "                sg_list += [seq]\n",
    "                out_str = \"{},{},pLX-sgRNA_Library,Lander\\n\".format(sgID, seq)\n",
    "                fout.write(out_str)\n",
    "            else:\n",
    "                print \"SgRNA ID {} with sequence {} already in library\".format(sgID, seq)\n",
    "            \n",
    "    #Look at Essential gene screen\n",
    "    with open(\"Published Libraries/Wang et al 2015 Supplementary Table 1 sgRNA Sequences.csv\", \"r\") as ess_f:\n",
    "        for line in ess_f:\n",
    "            if \"CTRL\" not in line:\n",
    "                continue\n",
    "            ele = line.split(\",\")\n",
    "            sgID = ele[0]\n",
    "            seq = ele[5]\n",
    "            if seq not in sg_list:\n",
    "                sg_list += [seq]\n",
    "                out_str = \"{},{},Essential-gene_Library,Sabatini\\n\".format(sgID, seq)\n",
    "                fout.write(out_str)\n",
    "            else:\n",
    "                print \"SgRNA ID {} with sequence {} already in library\".format(sgID, seq)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Making fastq File\n",
    "\n",
    "We then wanted to align these control sequences to the human genome to ensure no matches were found. A fastq file containing the sgRNA sequences was created:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "control_file = \"Control sgRNAs/All Literature Control sgRNAs.csv\"\n",
    "fastq_file = \"Control sgRNAs/lit_controls.fastq\"\n",
    "with open(fastq_file, \"w\") as fout, open(control_file, \"r\") as fin:\n",
    "    fin.next() # skip header line\n",
    "    for line in fin:\n",
    "        ele = line.split(\",\")\n",
    "        sg = ele[1]\n",
    "        sglen = len(sg)\n",
    "        out = \"@{}\\n{}\\n+\\n{}\\n\".format(ele[0], sg, \"I\"*sglen)\n",
    "        fout.write(out)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Align to hg19\n",
    "\n",
    "The control sgRNAs were then aligned to the hg19 version of the human genome using <a href=\"http://bowtie-bio.sourceforge.net/index.shtml\">bowtie</a>."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Server user name: jkurata\n",
      "Server password: ········\n"
     ]
    }
   ],
   "source": [
    "import data_processing.trim_align as ta\n",
    "\n",
    "# Create TrimAndAlign object\n",
    "taObj = ta.TrimAndAlign(\"Control sgRNAs/Literature Control Alignment Log.log\")\n",
    "# move fastq file to server\n",
    "taObj.fileToServer(\"Control sgRNAs/lit_controls.fastq\", \"lit_controls\", ext=\".fastq\")\n",
    "# align using bowtie\n",
    "taObj.align_bowtie(\"lit_controls\", \"../UCSC/hg19/bowtie-indexes/hg19\", options=\"-v 0 -a -p 4\")\n",
    "# move aligned sam file back from server\n",
    "taObj.fileFromServer(\"Control sgRNAs/\", \"lit_controls_bowtie-aligned\", ext=\".sam\")\n",
    "# delete files on server\n",
    "taObj.cleanUp(\"lit_controls_aligned\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Filter by alignment\n",
    "\n",
    "The sgRNAs were then filtered to remove those which aligned to the genome. The sgRNAs removed may not have cut the genome, because our alignment did not require the PAM sequence to be present, but are removed to be safe."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "with open(\"Control sgRNAs/filtered_lit_controls.csv\", \"w\") as fout:\n",
    "    with open(\"Control sgRNAs/lit_controls_bowtie-aligned.sam\", \"r\") as alCon:\n",
    "        for line in alCon:\n",
    "            # Pass on comment lines\n",
    "            if line[0] == \"@\":\n",
    "                continue\n",
    "            ele = line.split(\"\\t\")\n",
    "            # Skip if it does align to the genome\n",
    "            if ele[2] != \"*\" or ele[3] != \"0\":\n",
    "                continue\n",
    "            outStr = \"{},{}\\n\".format(ele[0], ele[9])\n",
    "            fout.write(outStr)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Sort\n",
    "\n",
    "To ensure everytime this code is run, the same control sgRNAs are selected, the sgRNAs were sorted before selection."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "\n",
    "def sort_df(df):\n",
    "    \"\"\"\n",
    "        Sorts the dataframe based on the value of the numbers which are part of the sgName string\n",
    "        from \n",
    "        http://stackoverflow.com/questions/37693600/how-to-sort-dataframe-based-on-particular-stringcolumns-using-python-pandas\n",
    "    \"\"\"\n",
    "    name_ser = df.loc[:, \"LiteratureSgRNAID\"].str.extract(\"(\\d+)\", expand=False)\n",
    "    df = df.assign(sort=pd.to_numeric(name_ser)) # Add new column with name 'sort'\n",
    "    df.sort_values(\"sort\", inplace=True)\n",
    "    df = df.drop(\"sort\", axis=1) # remove 'sort' column\n",
    "    return df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "df = pd.read_csv(\"Control sgRNAs/filtered_lit_controls.csv\", header=None, \n",
    "                 names=[\"LiteratureSgRNAID\", \"SgRNA\"])\n",
    "# divide the dataframe into parts by origin\n",
    "egfp_df = df[df[\"LiteratureSgRNAID\"].str.contains(\"EGFP\")]\n",
    "puro_df = df[df[\"LiteratureSgRNAID\"].str.contains(\"Puro\")]\n",
    "geckoA_df = df[df[\"LiteratureSgRNAID\"].str.contains(\"HGLibA\")]\n",
    "geckoB_df = df[df[\"LiteratureSgRNAID\"].str.contains(\"HGLibB\")]\n",
    "pLX_df = df[(df[\"LiteratureSgRNAID\"].str.contains(\"CTRL\")) & (df[\"LiteratureSgRNAID\"].str.len()==8)]\n",
    "ess_df = df[(df[\"LiteratureSgRNAID\"].str.contains(\"CTRL\")) & (df[\"LiteratureSgRNAID\"].str.len()==9)]\n",
    "\n",
    "part_dfs = [(egfp_df, \"EGFP\"), (puro_df, \"Puro\"), (geckoA_df, \"GeCKOv2\"), (geckoB_df, \"GeCKOv2\"), (pLX_df, \"pLX-sgRNA\"), \n",
    "            (ess_df, \"Essential\")]\n",
    "part_dfs_sorted = []\n",
    "n = 1\n",
    "for part_df, part_name in part_dfs:\n",
    "    sorted_df = sort_df(part_df)\n",
    "    sorted_df[\"SgRNAName\"] = \"\"    # create new column for name\n",
    "    for ind in sorted_df.index.tolist():\n",
    "        sorted_df.set_value(ind, \"SgRNAName\", \"control_{}_{}\".format(n, part_name))\n",
    "        n += 1\n",
    "    part_dfs_sorted += [sorted_df]\n",
    "    \n",
    "# put all of the partial dataframes back together\n",
    "full_df = pd.concat(part_dfs_sorted, ignore_index=True)\n",
    "# output to file\n",
    "full_df.to_csv(\"Control sgRNAs/Filtered and Sorted Literature Control sgRNAs.csv\", index=False)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Making ControlSgRNA Table\n",
    "\n",
    "To hold the control RNA sequences, along with the literature sgRNA ID, the ControlSgRNA table was created. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "import data_processing as dp\n",
    "\n",
    "def create_control_table(db_name, sql_version=\"MySQL\", firewall=False):\n",
    "    \"\"\"\n",
    "        Creates the ControlSgRNA table\n",
    "    \"\"\"\n",
    "    db_con = dp.DatabaseConnection(sql_version, db_name=db_name, firewall=firewall)\n",
    "    db_con.make_table(\"ControlSgRNA\", {\"SgRNAName\": [\"VARCHAR(50)\", \"NOT NULL\"], \"SgRNA\": [\"CHAR(20)\"], \n",
    "                                       \"LiteratureSgRNA\": [\"VARCHAR(200)\"]}, [\"PRIMARY KEY (SgRNAName)\"])\n",
    "    db_con.close_cursor()\n",
    "    db_con.close_connection()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Creating table ControlSgRNA\n",
      "CREATE TABLE ControlSgRNA (LiteratureSgRNA VARCHAR(200),\n",
      "    SgRNA CHAR(20),\n",
      "    SgRNAName VARCHAR(50) NOT NULL,\n",
      "    PRIMARY KEY (SgRNAName)) ENGINE=InnoDB;\n",
      "Sucessfully created table ControlSgRNA\n"
     ]
    }
   ],
   "source": [
    "create_control_table(\"miR-test\", firewall=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Fill ControlSgRNA Table\n",
    "\n",
    "The first 1000 control sgRNAs were then imported into the new table."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "import data_processing as dp\n",
    "\n",
    "def import_controls(db_name, sql_version=\"MySQL\", firewall=False):\n",
    "    df = pd.read_csv(\"Control sgRNAs/Filtered and Sorted Literature Control sgRNAs.csv\", \n",
    "                     header=0)\n",
    "    first_df = df.head(1000) # take first 1,000 rows\n",
    "    \n",
    "    insert_dict = {\"SgRNAName\": [], \"SgRNA\": [], \"LiteratureSgRNA\": []}\n",
    "    for i in first_df.index.tolist():\n",
    "        insert_dict[\"SgRNAName\"] += [first_df.get_value(i, \"SgRNAName\")]\n",
    "        insert_dict[\"SgRNA\"] += [first_df.get_value(i, \"SgRNA\")]\n",
    "        insert_dict[\"LiteratureSgRNA\"] += [first_df.get_value(i, \"LiteratureSgRNAID\")]\n",
    "    \n",
    "    db_con = dp.DatabaseConnection(sql_version, db_name=db_name, firewall=firewall)\n",
    "    db_con.make_many_rows(insert_dict, \"ControlSgRNA\")\n",
    "    db_con.close_cursor()\n",
    "    db_con.close_connection()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "import_controls(\"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
}