08_pairing_with_hg38.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "from pyfaidx import Fasta\n",
    "from Bio.Seq import Seq"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "outputs": [
    {
     "data": {
      "text/plain": "         chr     start       end strand  class\n0       chr1     12227     12612      +      1\n1       chr1     12721     13220      +      1\n2       chr1     12057     12178      +      1\n3       chr1     12697     12974      +      1\n4       chr1     13052     13220      +      1\n...      ...       ...       ...    ...    ...\n519029  chrY  25464577  25465486      +      0\n519030  chrY  25513173  25513588      -      0\n519031  chrY  25513745  25516715      -      0\n519032  chrY  25525288  25527646      -      0\n519033  chrY  57213357  57213855      -      0\n\n[519034 rows x 5 columns]",
      "text/html": "<div>\n<style scoped>\n    .dataframe tbody tr th:only-of-type {\n        vertical-align: middle;\n    }\n\n    .dataframe tbody tr th {\n        vertical-align: top;\n    }\n\n    .dataframe thead th {\n        text-align: right;\n    }\n</style>\n<table border=\"1\" class=\"dataframe\">\n  <thead>\n    <tr style=\"text-align: right;\">\n      <th></th>\n      <th>chr</th>\n      <th>start</th>\n      <th>end</th>\n      <th>strand</th>\n      <th>class</th>\n    </tr>\n  </thead>\n  <tbody>\n    <tr>\n      <th>0</th>\n      <td>chr1</td>\n      <td>12227</td>\n      <td>12612</td>\n      <td>+</td>\n      <td>1</td>\n    </tr>\n    <tr>\n      <th>1</th>\n      <td>chr1</td>\n      <td>12721</td>\n      <td>13220</td>\n      <td>+</td>\n      <td>1</td>\n    </tr>\n    <tr>\n      <th>2</th>\n      <td>chr1</td>\n      <td>12057</td>\n      <td>12178</td>\n      <td>+</td>\n      <td>1</td>\n    </tr>\n    <tr>\n      <th>3</th>\n      <td>chr1</td>\n      <td>12697</td>\n      <td>12974</td>\n      <td>+</td>\n      <td>1</td>\n    </tr>\n    <tr>\n      <th>4</th>\n      <td>chr1</td>\n      <td>13052</td>\n      <td>13220</td>\n      <td>+</td>\n      <td>1</td>\n    </tr>\n    <tr>\n      <th>...</th>\n      <td>...</td>\n      <td>...</td>\n      <td>...</td>\n      <td>...</td>\n      <td>...</td>\n    </tr>\n    <tr>\n      <th>519029</th>\n      <td>chrY</td>\n      <td>25464577</td>\n      <td>25465486</td>\n      <td>+</td>\n      <td>0</td>\n    </tr>\n    <tr>\n      <th>519030</th>\n      <td>chrY</td>\n      <td>25513173</td>\n      <td>25513588</td>\n      <td>-</td>\n      <td>0</td>\n    </tr>\n    <tr>\n      <th>519031</th>\n      <td>chrY</td>\n      <td>25513745</td>\n      <td>25516715</td>\n      <td>-</td>\n      <td>0</td>\n    </tr>\n    <tr>\n      <th>519032</th>\n      <td>chrY</td>\n      <td>25525288</td>\n      <td>25527646</td>\n      <td>-</td>\n      <td>0</td>\n    </tr>\n    <tr>\n      <th>519033</th>\n      <td>chrY</td>\n      <td>57213357</td>\n      <td>57213855</td>\n      <td>-</td>\n      <td>0</td>\n    </tr>\n  </tbody>\n</table>\n<p>519034 rows × 5 columns</p>\n</div>"
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "sha512sum = !sha512sum data/04_train_set.tsv\n",
    "assert(\"5006dfaea420edf91d4e86a72bc2428a57ad3b184a8e3b9cc55ec048d17e001247a36322b18b48fa19b24a784b8bd7fe33eb2cfa074d4625138f4728b56d3324\" in sha512sum[0])\n",
    "df = pd.read_csv(\"data/04_train_set.tsv\", sep=\"\\t\")\n",
    "df"
   ],
   "metadata": {
    "collapsed": false
   }
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "outputs": [],
   "source": [
    "# This is downloaded in https://hgdownload.soe.ucsc.edu/goldenPath/hg38/bigZips/\n",
    "genome = Fasta('data/human_ref_hg38_109/GRCh38.primary_assembly.genome.fa', sequence_always_upper=True)\n",
    "\n",
    "# A simple lambda function for matching the chromosome, start and end\n",
    "coords_to_dna = lambda start_c, end_c, chr_c: genome[chr_c][start_c-1:end_c]"
   ],
   "metadata": {
    "collapsed": false
   }
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "outputs": [],
   "source": [
    "# def find_sequence(row):\n",
    "#     chromosome = row[\"chr\"]\n",
    "#     start = row[\"start\"]\n",
    "#     stop = row[\"end\"]\n",
    "#     strand = row[\"strand\"]\n",
    "#\n",
    "#     seq_find_sequence = str(coords_to_dna(int(start)+1, int(stop), chromosome))\n",
    "#\n",
    "#     if strand == \"-\": #strand\n",
    "#         seq_find_sequence = Seq(seq_find_sequence)  # Encode the seqeunce into\n",
    "#         seq_find_sequence = seq_find_sequence.reverse_complement()\n",
    "#\n",
    "#\n",
    "#     return str(seq_find_sequence)\n",
    "\n",
    "def find_sequence_maxentscan(row):\n",
    "    chromosome = row[\"chr\"]\n",
    "    start = row[\"start\"]\n",
    "    stop = row[\"end\"]\n",
    "    strand = row[\"strand\"]\n",
    "\n",
    "    seq_find_sequence = str(coords_to_dna(int(start)+1-3, int(stop)+3, chromosome))\n",
    "\n",
    "    if strand == \"-\": #strand\n",
    "        seq_find_sequence = Seq(seq_find_sequence)  # Encode the sequence into\n",
    "        seq_find_sequence = seq_find_sequence.reverse_complement()\n",
    "\n",
    "    seq_find_sequence = str(seq_find_sequence)\n",
    "    first_three = seq_find_sequence[:3].lower()\n",
    "    last_three = seq_find_sequence[-3:].lower()\n",
    "    middle_part = seq_find_sequence[3:-3].upper()\n",
    "\n",
    "    return first_three + middle_part + last_three"
   ],
   "metadata": {
    "collapsed": false
   }
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "outputs": [],
   "source": [
    "# df_temp = df.copy().head(5000)\n",
    "# df_temp[\"sequence\"] = df_temp.apply(find_sequence, axis=1)\n",
    "# df_temp[\"sequence_maxentscan\"] = df_temp.apply(find_sequence_maxentscan, axis=1)\n",
    "# df_temp"
   ],
   "metadata": {
    "collapsed": false
   }
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "outputs": [],
   "source": [
    "df[\"maxentscan_sequence\"] = df.apply(find_sequence_maxentscan, axis=1)"
   ],
   "metadata": {
    "collapsed": false
   }
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "outputs": [
    {
     "data": {
      "text/plain": "         chr     start       end strand  class  \\\n0       chr1     12227     12612      +      1   \n1       chr1     12721     13220      +      1   \n2       chr1     12057     12178      +      1   \n3       chr1     12697     12974      +      1   \n4       chr1     13052     13220      +      1   \n...      ...       ...       ...    ...    ...   \n519029  chrY  25464577  25465486      +      0   \n519030  chrY  25513173  25513588      -      0   \n519031  chrY  25513745  25516715      -      0   \n519032  chrY  25525288  25527646      -      0   \n519033  chrY  57213357  57213855      -      0   \n\n                                      maxentscan_sequence start_ss end_ss  \n0       ccaGTAAGTAGTGCTTGTGCTCATCTCCTTGGCTGTGATACGTGGC...       GT     AG  \n1       gagGTGAGAGGAGAGTAGACAGTGAGTGGGAGTGGCGTCGCCCCTA...       GT     AG  \n2       gagCACTGGAGTGGAGTTTTCCTGTGGAGAGGAGCCATGCCTAGAG...       CA     GA  \n3       cttGTGAGTGTCCCCAGTGTTGCAGAGGTGAGAGGAGAGTAGACAG...       GT     AG  \n4       tagGCAAGCCTGGCTGCCTCCAGCTGGGTCGACAGACAGGGGCTGG...       GC     AG  \n...                                                   ...      ...    ...  \n519029  cacGTTGAGGCGCCCAGTGGCGGCCTCACGGGGCAGGGCGAGGGCG...       GT     AG  \n519030  aggGTAAGATCAGTGCTATTGTCAGAGGAAAAACTCCTGGCCATCA...       GT     AG  \n519031  gctGTAAGTTCCACATTGATTATCATAGGCTAACCATGGGCCAGGC...       GT     AG  \n519032  tgtGTGTATAAATATCTGGACTTTTTGGTTAAGTAATTATAGTTAA...       GT     AG  \n519033  gagGTGAGAGGAGAGTAGACAGTGAGTGGGAGTGGCGTCGCCCCTA...       GT     AG  \n\n[519034 rows x 8 columns]",
      "text/html": "<div>\n<style scoped>\n    .dataframe tbody tr th:only-of-type {\n        vertical-align: middle;\n    }\n\n    .dataframe tbody tr th {\n        vertical-align: top;\n    }\n\n    .dataframe thead th {\n        text-align: right;\n    }\n</style>\n<table border=\"1\" class=\"dataframe\">\n  <thead>\n    <tr style=\"text-align: right;\">\n      <th></th>\n      <th>chr</th>\n      <th>start</th>\n      <th>end</th>\n      <th>strand</th>\n      <th>class</th>\n      <th>maxentscan_sequence</th>\n      <th>start_ss</th>\n      <th>end_ss</th>\n    </tr>\n  </thead>\n  <tbody>\n    <tr>\n      <th>0</th>\n      <td>chr1</td>\n      <td>12227</td>\n      <td>12612</td>\n      <td>+</td>\n      <td>1</td>\n      <td>ccaGTAAGTAGTGCTTGTGCTCATCTCCTTGGCTGTGATACGTGGC...</td>\n      <td>GT</td>\n      <td>AG</td>\n    </tr>\n    <tr>\n      <th>1</th>\n      <td>chr1</td>\n      <td>12721</td>\n      <td>13220</td>\n      <td>+</td>\n      <td>1</td>\n      <td>gagGTGAGAGGAGAGTAGACAGTGAGTGGGAGTGGCGTCGCCCCTA...</td>\n      <td>GT</td>\n      <td>AG</td>\n    </tr>\n    <tr>\n      <th>2</th>\n      <td>chr1</td>\n      <td>12057</td>\n      <td>12178</td>\n      <td>+</td>\n      <td>1</td>\n      <td>gagCACTGGAGTGGAGTTTTCCTGTGGAGAGGAGCCATGCCTAGAG...</td>\n      <td>CA</td>\n      <td>GA</td>\n    </tr>\n    <tr>\n      <th>3</th>\n      <td>chr1</td>\n      <td>12697</td>\n      <td>12974</td>\n      <td>+</td>\n      <td>1</td>\n      <td>cttGTGAGTGTCCCCAGTGTTGCAGAGGTGAGAGGAGAGTAGACAG...</td>\n      <td>GT</td>\n      <td>AG</td>\n    </tr>\n    <tr>\n      <th>4</th>\n      <td>chr1</td>\n      <td>13052</td>\n      <td>13220</td>\n      <td>+</td>\n      <td>1</td>\n      <td>tagGCAAGCCTGGCTGCCTCCAGCTGGGTCGACAGACAGGGGCTGG...</td>\n      <td>GC</td>\n      <td>AG</td>\n    </tr>\n    <tr>\n      <th>...</th>\n      <td>...</td>\n      <td>...</td>\n      <td>...</td>\n      <td>...</td>\n      <td>...</td>\n      <td>...</td>\n      <td>...</td>\n      <td>...</td>\n    </tr>\n    <tr>\n      <th>519029</th>\n      <td>chrY</td>\n      <td>25464577</td>\n      <td>25465486</td>\n      <td>+</td>\n      <td>0</td>\n      <td>cacGTTGAGGCGCCCAGTGGCGGCCTCACGGGGCAGGGCGAGGGCG...</td>\n      <td>GT</td>\n      <td>AG</td>\n    </tr>\n    <tr>\n      <th>519030</th>\n      <td>chrY</td>\n      <td>25513173</td>\n      <td>25513588</td>\n      <td>-</td>\n      <td>0</td>\n      <td>aggGTAAGATCAGTGCTATTGTCAGAGGAAAAACTCCTGGCCATCA...</td>\n      <td>GT</td>\n      <td>AG</td>\n    </tr>\n    <tr>\n      <th>519031</th>\n      <td>chrY</td>\n      <td>25513745</td>\n      <td>25516715</td>\n      <td>-</td>\n      <td>0</td>\n      <td>gctGTAAGTTCCACATTGATTATCATAGGCTAACCATGGGCCAGGC...</td>\n      <td>GT</td>\n      <td>AG</td>\n    </tr>\n    <tr>\n      <th>519032</th>\n      <td>chrY</td>\n      <td>25525288</td>\n      <td>25527646</td>\n      <td>-</td>\n      <td>0</td>\n      <td>tgtGTGTATAAATATCTGGACTTTTTGGTTAAGTAATTATAGTTAA...</td>\n      <td>GT</td>\n      <td>AG</td>\n    </tr>\n    <tr>\n      <th>519033</th>\n      <td>chrY</td>\n      <td>57213357</td>\n      <td>57213855</td>\n      <td>-</td>\n      <td>0</td>\n      <td>gagGTGAGAGGAGAGTAGACAGTGAGTGGGAGTGGCGTCGCCCCTA...</td>\n      <td>GT</td>\n      <td>AG</td>\n    </tr>\n  </tbody>\n</table>\n<p>519034 rows × 8 columns</p>\n</div>"
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "df[\"start_ss\"] = df[\"maxentscan_sequence\"].str[3:5]\n",
    "df[\"end_ss\"] = df[\"maxentscan_sequence\"].str[-5:-3:]\n",
    "df"
   ],
   "metadata": {
    "collapsed": false
   }
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "outputs": [],
   "source": [
    "df.to_csv(\"data/08_trainset_with_seq.tsv\", sep=\"\\t\", index=False)"
   ],
   "metadata": {
    "collapsed": false
   }
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "159da72ef546de9f15bfc662d51b70987134035499c64674f4240279ed9d120d34b0312bf98f1cda8c4d43222df7175f818c6d53c9dde2487ea7cf747b51ab05  data/08_trainset_with_seq.tsv\r\n"
     ]
    }
   ],
   "source": [
    "!sha512sum data/08_trainset_with_seq.tsv"
   ],
   "metadata": {
    "collapsed": false
   }
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "outputs": [],
   "source": [],
   "metadata": {
    "collapsed": false
   }
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "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.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}