examples/CPTAC_Mapping_Workflow.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import clumpsptm\n",
    "import pandas as pd\n",
    "import os\n",
    "import numpy as np\n",
    "import glob\n",
    "from tqdm import tqdm\n",
    "import subprocess\n",
    "import ast\n",
    "from agutil.parallel import parallelize2\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "---\n",
    "\n",
    "__Author__:  Shankara Anand\n",
    "\n",
    "__Email__: sanand@broadinstitute.org\n",
    "\n",
    "This notebook walks through how to take a source dataset of PTM sites and map them to PDB structures to eventually use with PTM-CLUMPS. This is necessary pre-processing for a new dataset to map the feature space onto available structures. This notebook specifically walks through steps used in [CPTAC PanCancer project](https://github.com/getzlab/CPTAC_PanCan_2021).\n",
    "\n",
    "---"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Reference Files\n",
    "UNIPROT_SWISSPROT = \"https://ftp.uniprot.org/pub/databases/uniprot/current_release/knowledgebase/complete/uniprot_sprot.fasta.gz\"\n",
    "SIFTS_DB = \"ftp://ftp.ebi.ac.uk/pub/databases/msd/sifts/flatfiles/tsv/pdb_chain_uniprot.tsv.gz\"\n",
    "\n",
    "# Source Files\n",
    "CPTAC_REFSEQ_FASTA = \"../../CPTAC_PanCan_2021/data/ref/genomic/RefSeq.20180629_Human_ucsc_hg38_cpdbnr_mito_264contams_553smORFs.fasta\"\n",
    "CPTAC_FEATURE_FILE = \"../../CPTAC_PanCan_2021/data/processed/061721/var_map_full.tsv\"\n",
    "\n",
    "# Downloaded Database\n",
    "PDB_DIR = \"../../getzlab-CLUMPS2/clumps/db/ref/pdbs\"\n",
    "\n",
    "# Output Directory\n",
    "REF_DIR = \"example_refseq\"\n",
    "\n",
    "# Other\n",
    "n_proteins_for_demo = 100\n",
    "n_threads = 16"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1. Download Reference Files\n",
    "---\n",
    "\n",
    "* We are mapping our source fasta (from CPTAC) to the Uniprot Swissprot annotated database\n",
    "* We use the Sifts protein database for mapping Uniprot files to annotated PDB protein chains"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "   * Downloading uniprot_sprot.fasta.gz...\n",
      "   * Downloading pdb_chain_uniprot.tsv.gz...\n"
     ]
    }
   ],
   "source": [
    "clumpsptm.mp.dl_ref(REF_DIR, [UNIPROT_SWISSPROT, SIFTS_DB])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2. Split Source Fasta\n",
    "---\n",
    "\n",
    "* Splits source fasta (from CPTAC) into individual protein sequence fastas for easy paralell processing\n",
    "* Will select only 5 sequences to use blasting for example purposes"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "100%|██████████| 42311/42311 [00:03<00:00, 12557.30it/s]\n"
     ]
    }
   ],
   "source": [
    "clumpsptm.mp.split_fastas(\n",
    "    CPTAC_REFSEQ_FASTA,\n",
    "    os.path.join(REF_DIR, \"refseq_fastas\"),\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "['example_refseq/refseq_fastas/NP_061943.2.seq',\n",
       " 'example_refseq/refseq_fastas/NP_001316327.1.seq',\n",
       " 'example_refseq/refseq_fastas/NP_001275683.1.seq',\n",
       " 'example_refseq/refseq_fastas/NP_004205.2.seq',\n",
       " 'example_refseq/refseq_fastas/NP_001159737.1.seq',\n",
       " 'example_refseq/refseq_fastas/NP_066568.3.seq',\n",
       " 'example_refseq/refseq_fastas/NP_776170.2.seq',\n",
       " 'example_refseq/refseq_fastas/NP_001263215.1.seq',\n",
       " 'example_refseq/refseq_fastas/NP_001333622.1.seq',\n",
       " 'example_refseq/refseq_fastas/NP_001191226.1.seq']"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "refseq_individual_fastas = glob.glob(os.path.join(REF_DIR, \"refseq_fastas/*\"))[:n_proteins_for_demo]\n",
    "refseq_individual_fastas[:10]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 3. Blast Protein Sequences\n",
    "---\n",
    "\n",
    "* Use clumps-ptm module to run blastp+ with easy scaling via multi-threading\n",
    "* For full CPTAC source fasta (~22947 sequences to blast), this requires roughly 1 hr of run-time with an `n1-standard-96` on Google Cloud Compute using `384` threads"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "\n",
      "Building a new DB, current time: 07/26/2022 21:07:10\n",
      "New DB name:   /home/sanand/CLUMPS-PTM/examples/example_refseq/uniprot_db/uniprot_db\n",
      "New DB title:  UniprotDB\n",
      "Sequence type: Protein\n",
      "Deleted existing Protein BLAST database named /home/sanand/CLUMPS-PTM/examples/example_refseq/uniprot_db/uniprot_db\n",
      "Keep Linkouts: T\n",
      "Keep MBits: T\n",
      "Maximum file size: 1000000000B\n",
      "Adding sequences from FASTA; added 567483 sequences in 15.1239 seconds.\n",
      "\n",
      "   * Running 100 blast sequences using 16 threads\n",
      "CPU times: user 103 ms, sys: 542 ms, total: 645 ms\n",
      "Wall time: 1min 18s\n"
     ]
    }
   ],
   "source": [
    "%%time\n",
    "\n",
    "clumpsptm.mp.blast_sequences(\n",
    "    os.path.join(REF_DIR, \"uniprot_sprot.fasta\"),\n",
    "    refseq_individual_fastas,\n",
    "    output_dir=REF_DIR,\n",
    "    n_threads=n_threads,\n",
    "    db_title=\"UniprotDB\",\n",
    "    blast_dir_name=\"refseq_to_uniprot_blast\",\n",
    "    seq_db_name=\"uniprot_db\"\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 4. Cross-Reference with Sifts Data-Base\n",
    "---"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "sifts_df = pd.read_csv(os.path.join(REF_DIR, \"pdb_chain_uniprot.tsv\"), comment=\"#\", sep='\\t')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "['example_refseq/refseq_to_uniprot_blast/NP_039236.1.blasted.seq.gz',\n",
       " 'example_refseq/refseq_to_uniprot_blast/NP_001191231.1.blasted.seq.gz',\n",
       " 'example_refseq/refseq_to_uniprot_blast/NP_001289928.1.blasted.seq.gz',\n",
       " 'example_refseq/refseq_to_uniprot_blast/NP_001336478.1.blasted.seq.gz',\n",
       " 'example_refseq/refseq_to_uniprot_blast/NP_001616.1.blasted.seq.gz',\n",
       " 'example_refseq/refseq_to_uniprot_blast/NP_835454.1.blasted.seq.gz',\n",
       " 'example_refseq/refseq_to_uniprot_blast/NP_005083.2.blasted.seq.gz',\n",
       " 'example_refseq/refseq_to_uniprot_blast/NP_001308477.1.blasted.seq.gz',\n",
       " 'example_refseq/refseq_to_uniprot_blast/NP_001307809.1.blasted.seq.gz',\n",
       " 'example_refseq/refseq_to_uniprot_blast/NP_001186782.1.blasted.seq.gz']"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "blasted_files = glob.glob(os.path.join(REF_DIR, \"refseq_to_uniprot_blast/*\"))\n",
    "blasted_files[:10]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "100%|██████████| 100/100 [00:03<00:00, 29.29it/s]\n"
     ]
    },
    {
     "data": {
      "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>Hit_def</th>\n",
       "      <th>Hsp_num</th>\n",
       "      <th>Hsp_bit-score</th>\n",
       "      <th>Hsp_score</th>\n",
       "      <th>Hsp_evalue</th>\n",
       "      <th>Hsp_query-from</th>\n",
       "      <th>Hsp_query-to</th>\n",
       "      <th>Hsp_hit-from</th>\n",
       "      <th>Hsp_hit-to</th>\n",
       "      <th>Hsp_query-frame</th>\n",
       "      <th>...</th>\n",
       "      <th>Hsp_gaps</th>\n",
       "      <th>Hsp_align-len</th>\n",
       "      <th>Hsp_qseq</th>\n",
       "      <th>Hsp_hseq</th>\n",
       "      <th>Hsp_midline</th>\n",
       "      <th>relIdentity</th>\n",
       "      <th>uniprot</th>\n",
       "      <th>species</th>\n",
       "      <th>sifts</th>\n",
       "      <th>blast</th>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>query</th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>NP_039236.1</th>\n",
       "      <td>sp|P23760|PAX3_HUMAN Paired box protein Pax-3 ...</td>\n",
       "      <td>1</td>\n",
       "      <td>751.51</td>\n",
       "      <td>1939</td>\n",
       "      <td>0</td>\n",
       "      <td>1</td>\n",
       "      <td>473</td>\n",
       "      <td>1</td>\n",
       "      <td>478</td>\n",
       "      <td>0</td>\n",
       "      <td>...</td>\n",
       "      <td>19</td>\n",
       "      <td>485</td>\n",
       "      <td>MAALPGTVPRMMRPAPGQNYPRTGFPLEVSTPLGQGRVNQLGGVFI...</td>\n",
       "      <td>MTTLAGAVPRMMRPGPGQNYPRSGFPLEVSTPLGQGRVNQLGGVFI...</td>\n",
       "      <td>M  L G VPRMMRP PGQNYPR+GFPLEVSTPLGQGRVNQLGGVFI...</td>\n",
       "      <td>0.771134</td>\n",
       "      <td>P23760</td>\n",
       "      <td>HUMAN</td>\n",
       "      <td>True</td>\n",
       "      <td>True</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>NP_001191231.1</th>\n",
       "      <td>sp|P05067|A4_HUMAN Amyloid-beta precursor prot...</td>\n",
       "      <td>1</td>\n",
       "      <td>1511.12</td>\n",
       "      <td>3911</td>\n",
       "      <td>0</td>\n",
       "      <td>1</td>\n",
       "      <td>733</td>\n",
       "      <td>1</td>\n",
       "      <td>770</td>\n",
       "      <td>0</td>\n",
       "      <td>...</td>\n",
       "      <td>37</td>\n",
       "      <td>770</td>\n",
       "      <td>MLPGLALLLLAAWTARALEVPTDGNAGLLAEPQIAMFCGRLNMHMN...</td>\n",
       "      <td>MLPGLALLLLAAWTARALEVPTDGNAGLLAEPQIAMFCGRLNMHMN...</td>\n",
       "      <td>MLPGLALLLLAAWTARALEVPTDGNAGLLAEPQIAMFCGRLNMHMN...</td>\n",
       "      <td>0.950649</td>\n",
       "      <td>P05067</td>\n",
       "      <td>HUMAN</td>\n",
       "      <td>True</td>\n",
       "      <td>True</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>NP_001289928.1</th>\n",
       "      <td>sp|P81139|LIPR2_CAVPO Pancreatic lipase-relate...</td>\n",
       "      <td>1</td>\n",
       "      <td>153.295</td>\n",
       "      <td>386</td>\n",
       "      <td>6.28424e-40</td>\n",
       "      <td>17</td>\n",
       "      <td>334</td>\n",
       "      <td>13</td>\n",
       "      <td>349</td>\n",
       "      <td>0</td>\n",
       "      <td>...</td>\n",
       "      <td>45</td>\n",
       "      <td>350</td>\n",
       "      <td>NKRPCLEFSQLSVKDSFRDLFIPRIETILMMYTRNNLNCAEPLFEQ...</td>\n",
       "      <td>DEKPWAGTSQRPIKSLPSDP--KKINTRFLLYTNENQNSYQLITAT...</td>\n",
       "      <td>+++P    SQ  +K    D    +I T  ++YT  N N  + +   ...</td>\n",
       "      <td>0.297143</td>\n",
       "      <td>P81139</td>\n",
       "      <td>CAVPO</td>\n",
       "      <td>True</td>\n",
       "      <td>True</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>NP_001336478.1</th>\n",
       "      <td>sp|Q91YI6|GLD2_MOUSE Poly(A) RNA polymerase GL...</td>\n",
       "      <td>1</td>\n",
       "      <td>926.776</td>\n",
       "      <td>2394</td>\n",
       "      <td>0</td>\n",
       "      <td>1</td>\n",
       "      <td>509</td>\n",
       "      <td>1</td>\n",
       "      <td>484</td>\n",
       "      <td>0</td>\n",
       "      <td>...</td>\n",
       "      <td>25</td>\n",
       "      <td>509</td>\n",
       "      <td>MFPNSILGRPPFTPNHQQHNNFFTLSPTVYSHQQLIDAQFNFQNAD...</td>\n",
       "      <td>MFPNSILGRPPFTPTHQQHNNFFALSPTLYSHQQLIDAQFNFQNVD...</td>\n",
       "      <td>MFPNSILGRPPFTP HQQHNNFF LSPT+YSHQQLIDAQFNFQN D...</td>\n",
       "      <td>0.880157</td>\n",
       "      <td>Q91YI6</td>\n",
       "      <td>MOUSE</td>\n",
       "      <td>True</td>\n",
       "      <td>True</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>NP_001616.1</th>\n",
       "      <td>sp|P54819|KAD2_HUMAN Adenylate kinase 2, mitoc...</td>\n",
       "      <td>1</td>\n",
       "      <td>490.73</td>\n",
       "      <td>1262</td>\n",
       "      <td>5.20366e-176</td>\n",
       "      <td>1</td>\n",
       "      <td>239</td>\n",
       "      <td>1</td>\n",
       "      <td>239</td>\n",
       "      <td>0</td>\n",
       "      <td>...</td>\n",
       "      <td>0</td>\n",
       "      <td>239</td>\n",
       "      <td>MAPSVPAAEPEYPKGIRAVLLGPPGAGKGTQAPRLAENFCVCHLAT...</td>\n",
       "      <td>MAPSVPAAEPEYPKGIRAVLLGPPGAGKGTQAPRLAENFCVCHLAT...</td>\n",
       "      <td>MAPSVPAAEPEYPKGIRAVLLGPPGAGKGTQAPRLAENFCVCHLAT...</td>\n",
       "      <td>1.0</td>\n",
       "      <td>P54819</td>\n",
       "      <td>HUMAN</td>\n",
       "      <td>True</td>\n",
       "      <td>True</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",
       "      <td>...</td>\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",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>NP_937839.1</th>\n",
       "      <td>sp|P40200|TACT_HUMAN T-cell surface protein ta...</td>\n",
       "      <td>1</td>\n",
       "      <td>1210.67</td>\n",
       "      <td>3131</td>\n",
       "      <td>0</td>\n",
       "      <td>1</td>\n",
       "      <td>585</td>\n",
       "      <td>1</td>\n",
       "      <td>585</td>\n",
       "      <td>0</td>\n",
       "      <td>...</td>\n",
       "      <td>0</td>\n",
       "      <td>585</td>\n",
       "      <td>MEKKWKYCAVYYIIQIHFVKGVWEKTVNTEENVYATLGSDVNLTCQ...</td>\n",
       "      <td>MEKKWKYCAVYYIIQIHFVKGVWEKTVNTEENVYATLGSDVNLTCQ...</td>\n",
       "      <td>MEKKWKYCAVYYIIQIHFVKGVWEKTVNTEENVYATLGSDVNLTCQ...</td>\n",
       "      <td>1.0</td>\n",
       "      <td>P40200</td>\n",
       "      <td>HUMAN</td>\n",
       "      <td>True</td>\n",
       "      <td>True</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>NP_001284484.1</th>\n",
       "      <td>sp|O14733|MP2K7_HUMAN Dual specificity mitogen...</td>\n",
       "      <td>1</td>\n",
       "      <td>860.907</td>\n",
       "      <td>2223</td>\n",
       "      <td>0</td>\n",
       "      <td>1</td>\n",
       "      <td>435</td>\n",
       "      <td>1</td>\n",
       "      <td>419</td>\n",
       "      <td>0</td>\n",
       "      <td>...</td>\n",
       "      <td>16</td>\n",
       "      <td>435</td>\n",
       "      <td>MAASSLEQKLSRLEAKLKQENREARRRIDLNLDISPQRPRPIIVIT...</td>\n",
       "      <td>MAASSLEQKLSRLEAKLKQENREARRRIDLNLDISPQRPRP-----...</td>\n",
       "      <td>MAASSLEQKLSRLEAKLKQENREARRRIDLNLDISPQRPRP     ...</td>\n",
       "      <td>0.96092</td>\n",
       "      <td>O14733</td>\n",
       "      <td>HUMAN</td>\n",
       "      <td>True</td>\n",
       "      <td>True</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>NP_696961.2</th>\n",
       "      <td>sp|Q7L0J3|SV2A_HUMAN Synaptic vesicle glycopro...</td>\n",
       "      <td>1</td>\n",
       "      <td>73.1738</td>\n",
       "      <td>178</td>\n",
       "      <td>3.05011e-12</td>\n",
       "      <td>152</td>\n",
       "      <td>366</td>\n",
       "      <td>212</td>\n",
       "      <td>466</td>\n",
       "      <td>0</td>\n",
       "      <td>...</td>\n",
       "      <td>42</td>\n",
       "      <td>256</td>\n",
       "      <td>FAGVLVGAVAFGYLSDRFGRRRLLLVAYVSTLVLGLASAASVSYVM...</td>\n",
       "      <td>YLGMMVGAFLWGGLADRLGRRQCLLISLSVNSVFAFFSSFVQGYGT...</td>\n",
       "      <td>+ G++VGA  +G L+DR GRR+ LL++     V    S+    Y  ...</td>\n",
       "      <td>0.269531</td>\n",
       "      <td>Q7L0J3</td>\n",
       "      <td>HUMAN</td>\n",
       "      <td>True</td>\n",
       "      <td>True</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>NP_963840.2</th>\n",
       "      <td>sp|O43897|TLL1_HUMAN Tolloid-like protein 1 OS...</td>\n",
       "      <td>1</td>\n",
       "      <td>221.476</td>\n",
       "      <td>563</td>\n",
       "      <td>6.14567e-62</td>\n",
       "      <td>30</td>\n",
       "      <td>379</td>\n",
       "      <td>462</td>\n",
       "      <td>895</td>\n",
       "      <td>0</td>\n",
       "      <td>...</td>\n",
       "      <td>94</td>\n",
       "      <td>439</td>\n",
       "      <td>CGGVLSAPSGNFSSPNFPRLYPYNTECSWLIVVAEGSSVLLTFHAF...</td>\n",
       "      <td>CGGEIRKNEGQIQSPNYPDDYRPMKECVWKITVSESYHVGLTFQSF...</td>\n",
       "      <td>CGG +    G   SPN+P  Y    EC W I V+E   V LTF +F...</td>\n",
       "      <td>0.291572</td>\n",
       "      <td>O43897</td>\n",
       "      <td>HUMAN</td>\n",
       "      <td>True</td>\n",
       "      <td>True</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>NP_000057.2</th>\n",
       "      <td>sp|P07358|CO8B_HUMAN Complement component C8 b...</td>\n",
       "      <td>1</td>\n",
       "      <td>1240.71</td>\n",
       "      <td>3209</td>\n",
       "      <td>0</td>\n",
       "      <td>1</td>\n",
       "      <td>591</td>\n",
       "      <td>1</td>\n",
       "      <td>591</td>\n",
       "      <td>0</td>\n",
       "      <td>...</td>\n",
       "      <td>0</td>\n",
       "      <td>591</td>\n",
       "      <td>MKNSRTWAWRAPVELFLLCAALGCLSLPGSRGERPHSFGSNAVNKS...</td>\n",
       "      <td>MKNSRTWAWRAPVELFLLCAALGCLSLPGSRGERPHSFGSNAVNKS...</td>\n",
       "      <td>MKNSRTWAWRAPVELFLLCAALGCLSLPGSRGERPHSFGSNAVNKS...</td>\n",
       "      <td>1.0</td>\n",
       "      <td>P07358</td>\n",
       "      <td>HUMAN</td>\n",
       "      <td>True</td>\n",
       "      <td>True</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>100 rows × 23 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "                                                          Hit_def Hsp_num  \\\n",
       "query                                                                       \n",
       "NP_039236.1     sp|P23760|PAX3_HUMAN Paired box protein Pax-3 ...       1   \n",
       "NP_001191231.1  sp|P05067|A4_HUMAN Amyloid-beta precursor prot...       1   \n",
       "NP_001289928.1  sp|P81139|LIPR2_CAVPO Pancreatic lipase-relate...       1   \n",
       "NP_001336478.1  sp|Q91YI6|GLD2_MOUSE Poly(A) RNA polymerase GL...       1   \n",
       "NP_001616.1     sp|P54819|KAD2_HUMAN Adenylate kinase 2, mitoc...       1   \n",
       "...                                                           ...     ...   \n",
       "NP_937839.1     sp|P40200|TACT_HUMAN T-cell surface protein ta...       1   \n",
       "NP_001284484.1  sp|O14733|MP2K7_HUMAN Dual specificity mitogen...       1   \n",
       "NP_696961.2     sp|Q7L0J3|SV2A_HUMAN Synaptic vesicle glycopro...       1   \n",
       "NP_963840.2     sp|O43897|TLL1_HUMAN Tolloid-like protein 1 OS...       1   \n",
       "NP_000057.2     sp|P07358|CO8B_HUMAN Complement component C8 b...       1   \n",
       "\n",
       "               Hsp_bit-score Hsp_score    Hsp_evalue Hsp_query-from  \\\n",
       "query                                                                 \n",
       "NP_039236.1           751.51      1939             0              1   \n",
       "NP_001191231.1       1511.12      3911             0              1   \n",
       "NP_001289928.1       153.295       386   6.28424e-40             17   \n",
       "NP_001336478.1       926.776      2394             0              1   \n",
       "NP_001616.1           490.73      1262  5.20366e-176              1   \n",
       "...                      ...       ...           ...            ...   \n",
       "NP_937839.1          1210.67      3131             0              1   \n",
       "NP_001284484.1       860.907      2223             0              1   \n",
       "NP_696961.2          73.1738       178   3.05011e-12            152   \n",
       "NP_963840.2          221.476       563   6.14567e-62             30   \n",
       "NP_000057.2          1240.71      3209             0              1   \n",
       "\n",
       "               Hsp_query-to Hsp_hit-from Hsp_hit-to Hsp_query-frame  ...  \\\n",
       "query                                                                ...   \n",
       "NP_039236.1             473            1        478               0  ...   \n",
       "NP_001191231.1          733            1        770               0  ...   \n",
       "NP_001289928.1          334           13        349               0  ...   \n",
       "NP_001336478.1          509            1        484               0  ...   \n",
       "NP_001616.1             239            1        239               0  ...   \n",
       "...                     ...          ...        ...             ...  ...   \n",
       "NP_937839.1             585            1        585               0  ...   \n",
       "NP_001284484.1          435            1        419               0  ...   \n",
       "NP_696961.2             366          212        466               0  ...   \n",
       "NP_963840.2             379          462        895               0  ...   \n",
       "NP_000057.2             591            1        591               0  ...   \n",
       "\n",
       "               Hsp_gaps Hsp_align-len  \\\n",
       "query                                   \n",
       "NP_039236.1          19           485   \n",
       "NP_001191231.1       37           770   \n",
       "NP_001289928.1       45           350   \n",
       "NP_001336478.1       25           509   \n",
       "NP_001616.1           0           239   \n",
       "...                 ...           ...   \n",
       "NP_937839.1           0           585   \n",
       "NP_001284484.1       16           435   \n",
       "NP_696961.2          42           256   \n",
       "NP_963840.2          94           439   \n",
       "NP_000057.2           0           591   \n",
       "\n",
       "                                                         Hsp_qseq  \\\n",
       "query                                                               \n",
       "NP_039236.1     MAALPGTVPRMMRPAPGQNYPRTGFPLEVSTPLGQGRVNQLGGVFI...   \n",
       "NP_001191231.1  MLPGLALLLLAAWTARALEVPTDGNAGLLAEPQIAMFCGRLNMHMN...   \n",
       "NP_001289928.1  NKRPCLEFSQLSVKDSFRDLFIPRIETILMMYTRNNLNCAEPLFEQ...   \n",
       "NP_001336478.1  MFPNSILGRPPFTPNHQQHNNFFTLSPTVYSHQQLIDAQFNFQNAD...   \n",
       "NP_001616.1     MAPSVPAAEPEYPKGIRAVLLGPPGAGKGTQAPRLAENFCVCHLAT...   \n",
       "...                                                           ...   \n",
       "NP_937839.1     MEKKWKYCAVYYIIQIHFVKGVWEKTVNTEENVYATLGSDVNLTCQ...   \n",
       "NP_001284484.1  MAASSLEQKLSRLEAKLKQENREARRRIDLNLDISPQRPRPIIVIT...   \n",
       "NP_696961.2     FAGVLVGAVAFGYLSDRFGRRRLLLVAYVSTLVLGLASAASVSYVM...   \n",
       "NP_963840.2     CGGVLSAPSGNFSSPNFPRLYPYNTECSWLIVVAEGSSVLLTFHAF...   \n",
       "NP_000057.2     MKNSRTWAWRAPVELFLLCAALGCLSLPGSRGERPHSFGSNAVNKS...   \n",
       "\n",
       "                                                         Hsp_hseq  \\\n",
       "query                                                               \n",
       "NP_039236.1     MTTLAGAVPRMMRPGPGQNYPRSGFPLEVSTPLGQGRVNQLGGVFI...   \n",
       "NP_001191231.1  MLPGLALLLLAAWTARALEVPTDGNAGLLAEPQIAMFCGRLNMHMN...   \n",
       "NP_001289928.1  DEKPWAGTSQRPIKSLPSDP--KKINTRFLLYTNENQNSYQLITAT...   \n",
       "NP_001336478.1  MFPNSILGRPPFTPTHQQHNNFFALSPTLYSHQQLIDAQFNFQNVD...   \n",
       "NP_001616.1     MAPSVPAAEPEYPKGIRAVLLGPPGAGKGTQAPRLAENFCVCHLAT...   \n",
       "...                                                           ...   \n",
       "NP_937839.1     MEKKWKYCAVYYIIQIHFVKGVWEKTVNTEENVYATLGSDVNLTCQ...   \n",
       "NP_001284484.1  MAASSLEQKLSRLEAKLKQENREARRRIDLNLDISPQRPRP-----...   \n",
       "NP_696961.2     YLGMMVGAFLWGGLADRLGRRQCLLISLSVNSVFAFFSSFVQGYGT...   \n",
       "NP_963840.2     CGGEIRKNEGQIQSPNYPDDYRPMKECVWKITVSESYHVGLTFQSF...   \n",
       "NP_000057.2     MKNSRTWAWRAPVELFLLCAALGCLSLPGSRGERPHSFGSNAVNKS...   \n",
       "\n",
       "                                                      Hsp_midline relIdentity  \\\n",
       "query                                                                           \n",
       "NP_039236.1     M  L G VPRMMRP PGQNYPR+GFPLEVSTPLGQGRVNQLGGVFI...    0.771134   \n",
       "NP_001191231.1  MLPGLALLLLAAWTARALEVPTDGNAGLLAEPQIAMFCGRLNMHMN...    0.950649   \n",
       "NP_001289928.1  +++P    SQ  +K    D    +I T  ++YT  N N  + +   ...    0.297143   \n",
       "NP_001336478.1  MFPNSILGRPPFTP HQQHNNFF LSPT+YSHQQLIDAQFNFQN D...    0.880157   \n",
       "NP_001616.1     MAPSVPAAEPEYPKGIRAVLLGPPGAGKGTQAPRLAENFCVCHLAT...         1.0   \n",
       "...                                                           ...         ...   \n",
       "NP_937839.1     MEKKWKYCAVYYIIQIHFVKGVWEKTVNTEENVYATLGSDVNLTCQ...         1.0   \n",
       "NP_001284484.1  MAASSLEQKLSRLEAKLKQENREARRRIDLNLDISPQRPRP     ...     0.96092   \n",
       "NP_696961.2     + G++VGA  +G L+DR GRR+ LL++     V    S+    Y  ...    0.269531   \n",
       "NP_963840.2     CGG +    G   SPN+P  Y    EC W I V+E   V LTF +F...    0.291572   \n",
       "NP_000057.2     MKNSRTWAWRAPVELFLLCAALGCLSLPGSRGERPHSFGSNAVNKS...         1.0   \n",
       "\n",
       "               uniprot species sifts blast  \n",
       "query                                       \n",
       "NP_039236.1     P23760   HUMAN  True  True  \n",
       "NP_001191231.1  P05067   HUMAN  True  True  \n",
       "NP_001289928.1  P81139   CAVPO  True  True  \n",
       "NP_001336478.1  Q91YI6   MOUSE  True  True  \n",
       "NP_001616.1     P54819   HUMAN  True  True  \n",
       "...                ...     ...   ...   ...  \n",
       "NP_937839.1     P40200   HUMAN  True  True  \n",
       "NP_001284484.1  O14733   HUMAN  True  True  \n",
       "NP_696961.2     Q7L0J3   HUMAN  True  True  \n",
       "NP_963840.2     O43897   HUMAN  True  True  \n",
       "NP_000057.2     P07358   HUMAN  True  True  \n",
       "\n",
       "[100 rows x 23 columns]"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "mapped_acc_df = clumpsptm.mp.get_blast_hits_with_sifts(blasted_files, sifts_df)\n",
    "#mapped_acc_df = clumpsptm.mp.get_blast_hits_with_sifts(blasted_files, sifts_df, filter_only_human=True)\n",
    "mapped_acc_df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      " 100 / 100 (100.00 %) with blast results.\n",
      " 78 / 100 (78.00 %) with sifts results.\n"
     ]
    }
   ],
   "source": [
    "w_sifts = mapped_acc_df[(mapped_acc_df['blast']) & (mapped_acc_df['sifts'])].shape[0]\n",
    "w_blast = mapped_acc_df[mapped_acc_df['blast']].shape[0]\n",
    "tot = mapped_acc_df.shape[0]\n",
    "\n",
    "print(\" {} / {} ({:.2f} %) with blast results.\".format(w_blast, tot, 100*w_blast/tot))\n",
    "print(\" {} / {} ({:.2f} %) with sifts results.\".format(w_sifts, tot, 100*w_sifts/tot))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 5. Filter SIFTS with Downloaded PDBs\n",
    "---"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "PDBStore\n",
      "   * 112098 PDB files downloaded\n",
      "   * Directory Path: ../../getzlab-CLUMPS2/clumps/db/ref/pdbs\n"
     ]
    }
   ],
   "source": [
    "pdbstore = clumpsptm.PdbStore(PDB_DIR)\n",
    "print(pdbstore)\n",
    "\n",
    "# Sifts that overlap Uniprots\n",
    "sifts_filt_df = sifts_df[sifts_df['SP_PRIMARY'].isin(mapped_acc_df['uniprot'])]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Missing 36 / 560 (6.43 %) PDB files.\n"
     ]
    }
   ],
   "source": [
    "missing_pdbs = set(sifts_filt_df['PDB']) - pdbstore.downloaded_pdbs\n",
    "\n",
    "print(\"Missing {} / {} ({:.2f} %) PDB files.\".format(\n",
    "    len(missing_pdbs), \n",
    "    len(set(sifts_filt_df['PDB'])), \n",
    "    (100*len(missing_pdbs))/(len(set(sifts_filt_df['PDB'])))))\n",
    "\n",
    "# Drop PDBs unable to be downloaded\n",
    "sifts_filt_df = sifts_filt_df[sifts_filt_df[\"PDB\"].isin(pdbstore.downloaded_pdbs)]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 6. Match PTMs from Source Data-Set to Uniprot-Mapped Proteins\n",
    "---\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/home/sanand/anaconda3/envs/py37/lib/python3.7/site-packages/IPython/core/interactiveshell.py:3063: DtypeWarning: Columns (1,4,6,7,8,9,10) have mixed types.Specify dtype option on import or set low_memory=False.\n",
      "  interactivity=interactivity, compiler=compiler, result=result)\n"
     ]
    },
    {
     "data": {
      "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>id.description</th>\n",
       "      <th>geneSymbol</th>\n",
       "      <th>protein_mw</th>\n",
       "      <th>accession_number</th>\n",
       "      <th>feature</th>\n",
       "      <th>variableSites</th>\n",
       "      <th>sequence</th>\n",
       "      <th>sequenceVML</th>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>id</th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>NP_001333374.1_S18s_1_1_18_18</th>\n",
       "      <td>neuroblast differentiation-associated protein...</td>\n",
       "      <td>AHNAK</td>\n",
       "      <td>629603.7</td>\n",
       "      <td>NP_001333374.1</td>\n",
       "      <td>phosphoproteome</td>\n",
       "      <td>['S18s', 'S18s']</td>\n",
       "      <td>['ELLLPNWQGsGSHGLTIAQR', 'ELLLPNWQGsGSHGLTIAQR']</td>\n",
       "      <td>['ELLLPNWQGS(0.99)GS(0.0)HGLT(0.0)IAQR', 'ELLL...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>NP_001333374.1_S41s_1_1_41_41</th>\n",
       "      <td>neuroblast differentiation-associated protein...</td>\n",
       "      <td>AHNAK</td>\n",
       "      <td>629603.7</td>\n",
       "      <td>NP_001333374.1</td>\n",
       "      <td>phosphoproteome</td>\n",
       "      <td>['S41s', 'S41s']</td>\n",
       "      <td>['DDGVFVQEVTQNsPAAR', 'DDGVFVQEVTQNsPAAR']</td>\n",
       "      <td>['DDGVFVQEVT(0.0)QNS(0.99)PAAR', 'DDGVFVQEVT(0...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>NP_001333374.1_S93s_1_1_93_93</th>\n",
       "      <td>neuroblast differentiation-associated protein...</td>\n",
       "      <td>AHNAK</td>\n",
       "      <td>629603.7</td>\n",
       "      <td>NP_001333374.1</td>\n",
       "      <td>phosphoproteome</td>\n",
       "      <td>['S93s', 'S93s', 'S93s', 'S93s']</td>\n",
       "      <td>['GDRsPEPGQTWTR', 'KGDRsPEPGQTWTR', 'KGDRsPEPG...</td>\n",
       "      <td>['GDRS(0.99)PEPGQT(0.0)WT(0.0)R', 'KGDRS(0.99)...</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "                                                                  id.description  \\\n",
       "id                                                                                 \n",
       "NP_001333374.1_S18s_1_1_18_18   neuroblast differentiation-associated protein...   \n",
       "NP_001333374.1_S41s_1_1_41_41   neuroblast differentiation-associated protein...   \n",
       "NP_001333374.1_S93s_1_1_93_93   neuroblast differentiation-associated protein...   \n",
       "\n",
       "                              geneSymbol  protein_mw accession_number  \\\n",
       "id                                                                      \n",
       "NP_001333374.1_S18s_1_1_18_18      AHNAK    629603.7   NP_001333374.1   \n",
       "NP_001333374.1_S41s_1_1_41_41      AHNAK    629603.7   NP_001333374.1   \n",
       "NP_001333374.1_S93s_1_1_93_93      AHNAK    629603.7   NP_001333374.1   \n",
       "\n",
       "                                       feature  \\\n",
       "id                                               \n",
       "NP_001333374.1_S18s_1_1_18_18  phosphoproteome   \n",
       "NP_001333374.1_S41s_1_1_41_41  phosphoproteome   \n",
       "NP_001333374.1_S93s_1_1_93_93  phosphoproteome   \n",
       "\n",
       "                                                  variableSites  \\\n",
       "id                                                                \n",
       "NP_001333374.1_S18s_1_1_18_18                  ['S18s', 'S18s']   \n",
       "NP_001333374.1_S41s_1_1_41_41                  ['S41s', 'S41s']   \n",
       "NP_001333374.1_S93s_1_1_93_93  ['S93s', 'S93s', 'S93s', 'S93s']   \n",
       "\n",
       "                                                                        sequence  \\\n",
       "id                                                                                 \n",
       "NP_001333374.1_S18s_1_1_18_18   ['ELLLPNWQGsGSHGLTIAQR', 'ELLLPNWQGsGSHGLTIAQR']   \n",
       "NP_001333374.1_S41s_1_1_41_41         ['DDGVFVQEVTQNsPAAR', 'DDGVFVQEVTQNsPAAR']   \n",
       "NP_001333374.1_S93s_1_1_93_93  ['GDRsPEPGQTWTR', 'KGDRsPEPGQTWTR', 'KGDRsPEPG...   \n",
       "\n",
       "                                                                     sequenceVML  \n",
       "id                                                                                \n",
       "NP_001333374.1_S18s_1_1_18_18  ['ELLLPNWQGS(0.99)GS(0.0)HGLT(0.0)IAQR', 'ELLL...  \n",
       "NP_001333374.1_S41s_1_1_41_41  ['DDGVFVQEVT(0.0)QNS(0.99)PAAR', 'DDGVFVQEVT(0...  \n",
       "NP_001333374.1_S93s_1_1_93_93  ['GDRS(0.99)PEPGQT(0.0)WT(0.0)R', 'KGDRS(0.99)...  "
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "pmap_df = pd.read_csv(CPTAC_FEATURE_FILE, sep='\\t', index_col=0)\n",
    "ptm_df = pmap_df[pmap_df['feature'].isin(['phosphoproteome','acetylome'])].dropna(1)\n",
    "ptm_df.head(3)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/home/sanand/anaconda3/envs/py37/lib/python3.7/site-packages/pandas/core/indexing.py:1597: SettingWithCopyWarning: \n",
      "A value is trying to be set on a copy of a slice from a DataFrame.\n",
      "Try using .loc[row_indexer,col_indexer] = value instead\n",
      "\n",
      "See the caveats in the documentation: https://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy\n",
      "  self.obj[key] = value\n",
      "/home/sanand/anaconda3/envs/py37/lib/python3.7/site-packages/pandas/core/indexing.py:1676: SettingWithCopyWarning: \n",
      "A value is trying to be set on a copy of a slice from a DataFrame.\n",
      "Try using .loc[row_indexer,col_indexer] = value instead\n",
      "\n",
      "See the caveats in the documentation: https://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy\n",
      "  self._setitem_single_column(ilocs[0], value, pi)\n"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "  * 133619 single PTM sites total in CPTAC dataset\n",
      "  * 41394 multi PTM sites total in CPTAC dataset\n"
     ]
    },
    {
     "data": {
      "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>id.description</th>\n",
       "      <th>geneSymbol</th>\n",
       "      <th>protein_mw</th>\n",
       "      <th>accession_number</th>\n",
       "      <th>feature</th>\n",
       "      <th>variableSites</th>\n",
       "      <th>sequence</th>\n",
       "      <th>sequenceVML</th>\n",
       "      <th>ptmSite</th>\n",
       "      <th>Hit_def</th>\n",
       "      <th>...</th>\n",
       "      <th>Hsp_qseq</th>\n",
       "      <th>Hsp_hseq</th>\n",
       "      <th>Hsp_midline</th>\n",
       "      <th>relIdentity</th>\n",
       "      <th>uniprot</th>\n",
       "      <th>species</th>\n",
       "      <th>sifts</th>\n",
       "      <th>blast</th>\n",
       "      <th>acc_res</th>\n",
       "      <th>acc_res_i</th>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>id</th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>NP_001158136.1_S137s_1_1_137_137</th>\n",
       "      <td>microtubule-associated serine/threonine-prote...</td>\n",
       "      <td>MAST4</td>\n",
       "      <td>286317.4</td>\n",
       "      <td>NP_001158136.1</td>\n",
       "      <td>phosphoproteome</td>\n",
       "      <td>S137s</td>\n",
       "      <td>['KCsNPDVASGPGK', 'KCsNPDVASGPGK']</td>\n",
       "      <td>['KCS(0.99)NPDVAS(0.0)GPGK', 'KCS(0.99)NPDVAS(...</td>\n",
       "      <td>S137s</td>\n",
       "      <td>sp|O15021|MAST4_HUMAN Microtubule-associated s...</td>\n",
       "      <td>...</td>\n",
       "      <td>MGEKVSEAPEPVPRGCSGHGSRTPASALVAASSPGASSAESSSGSE...</td>\n",
       "      <td>MGEKVSEAPEPVPRGCSGHGSRTPASALVAASSPGASSAESSSGSE...</td>\n",
       "      <td>MGEKVSEAPEPVPRGCSGHGSRTPASALVAASSPGASSAESSSGSE...</td>\n",
       "      <td>1.0</td>\n",
       "      <td>O15021</td>\n",
       "      <td>HUMAN</td>\n",
       "      <td>True</td>\n",
       "      <td>True</td>\n",
       "      <td>S</td>\n",
       "      <td>137</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>NP_001158136.1_S1283s_1_1_1283_1283</th>\n",
       "      <td>microtubule-associated serine/threonine-prote...</td>\n",
       "      <td>MAST4</td>\n",
       "      <td>286317.4</td>\n",
       "      <td>NP_001158136.1</td>\n",
       "      <td>phosphoproteome</td>\n",
       "      <td>S1283s</td>\n",
       "      <td>['sFSCLNR', 'sFSCLNR']</td>\n",
       "      <td>['S(0.99)FS(0.0)CLNR', 'S(0.99)FS(0.0)CLNR']</td>\n",
       "      <td>S1283s</td>\n",
       "      <td>sp|O15021|MAST4_HUMAN Microtubule-associated s...</td>\n",
       "      <td>...</td>\n",
       "      <td>MGEKVSEAPEPVPRGCSGHGSRTPASALVAASSPGASSAESSSGSE...</td>\n",
       "      <td>MGEKVSEAPEPVPRGCSGHGSRTPASALVAASSPGASSAESSSGSE...</td>\n",
       "      <td>MGEKVSEAPEPVPRGCSGHGSRTPASALVAASSPGASSAESSSGSE...</td>\n",
       "      <td>1.0</td>\n",
       "      <td>O15021</td>\n",
       "      <td>HUMAN</td>\n",
       "      <td>True</td>\n",
       "      <td>True</td>\n",
       "      <td>S</td>\n",
       "      <td>1283</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>NP_001158136.1_S1309s_1_1_1309_1309</th>\n",
       "      <td>microtubule-associated serine/threonine-prote...</td>\n",
       "      <td>MAST4</td>\n",
       "      <td>286317.4</td>\n",
       "      <td>NP_001158136.1</td>\n",
       "      <td>phosphoproteome</td>\n",
       "      <td>S1309s</td>\n",
       "      <td>['sPTPSYR', 'sPTPSYR']</td>\n",
       "      <td>['S(0.99)PT(0.0)PS(0.0)Y(0.0)R', 'S(0.99)PT(0....</td>\n",
       "      <td>S1309s</td>\n",
       "      <td>sp|O15021|MAST4_HUMAN Microtubule-associated s...</td>\n",
       "      <td>...</td>\n",
       "      <td>MGEKVSEAPEPVPRGCSGHGSRTPASALVAASSPGASSAESSSGSE...</td>\n",
       "      <td>MGEKVSEAPEPVPRGCSGHGSRTPASALVAASSPGASSAESSSGSE...</td>\n",
       "      <td>MGEKVSEAPEPVPRGCSGHGSRTPASALVAASSPGASSAESSSGSE...</td>\n",
       "      <td>1.0</td>\n",
       "      <td>O15021</td>\n",
       "      <td>HUMAN</td>\n",
       "      <td>True</td>\n",
       "      <td>True</td>\n",
       "      <td>S</td>\n",
       "      <td>1309</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>NP_001158136.1_S1370s_1_1_1370_1370</th>\n",
       "      <td>microtubule-associated serine/threonine-prote...</td>\n",
       "      <td>MAST4</td>\n",
       "      <td>286317.4</td>\n",
       "      <td>NP_001158136.1</td>\n",
       "      <td>phosphoproteome</td>\n",
       "      <td>S1370s</td>\n",
       "      <td>['sAGNIPLSPLAR', 'sAGNIPLSPLAR']</td>\n",
       "      <td>['S(0.99)AGNIPLS(0.0)PLAR', 'S(0.99)AGNIPLS(0....</td>\n",
       "      <td>S1370s</td>\n",
       "      <td>sp|O15021|MAST4_HUMAN Microtubule-associated s...</td>\n",
       "      <td>...</td>\n",
       "      <td>MGEKVSEAPEPVPRGCSGHGSRTPASALVAASSPGASSAESSSGSE...</td>\n",
       "      <td>MGEKVSEAPEPVPRGCSGHGSRTPASALVAASSPGASSAESSSGSE...</td>\n",
       "      <td>MGEKVSEAPEPVPRGCSGHGSRTPASALVAASSPGASSAESSSGSE...</td>\n",
       "      <td>1.0</td>\n",
       "      <td>O15021</td>\n",
       "      <td>HUMAN</td>\n",
       "      <td>True</td>\n",
       "      <td>True</td>\n",
       "      <td>S</td>\n",
       "      <td>1370</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>NP_001158136.1_S1384s_1_1_1384_1384</th>\n",
       "      <td>microtubule-associated serine/threonine-prote...</td>\n",
       "      <td>MAST4</td>\n",
       "      <td>286317.4</td>\n",
       "      <td>NP_001158136.1</td>\n",
       "      <td>phosphoproteome</td>\n",
       "      <td>S1384s</td>\n",
       "      <td>['TPsPTPQPTSPQR', 'TPsPTPQPTSPQR']</td>\n",
       "      <td>['T(0.0)PS(0.99)PT(0.0)PQPT(0.0)S(0.0)PQR', 'T...</td>\n",
       "      <td>S1384s</td>\n",
       "      <td>sp|O15021|MAST4_HUMAN Microtubule-associated s...</td>\n",
       "      <td>...</td>\n",
       "      <td>MGEKVSEAPEPVPRGCSGHGSRTPASALVAASSPGASSAESSSGSE...</td>\n",
       "      <td>MGEKVSEAPEPVPRGCSGHGSRTPASALVAASSPGASSAESSSGSE...</td>\n",
       "      <td>MGEKVSEAPEPVPRGCSGHGSRTPASALVAASSPGASSAESSSGSE...</td>\n",
       "      <td>1.0</td>\n",
       "      <td>O15021</td>\n",
       "      <td>HUMAN</td>\n",
       "      <td>True</td>\n",
       "      <td>True</td>\n",
       "      <td>S</td>\n",
       "      <td>1384</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>5 rows × 34 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "                                                                        id.description  \\\n",
       "id                                                                                       \n",
       "NP_001158136.1_S137s_1_1_137_137      microtubule-associated serine/threonine-prote...   \n",
       "NP_001158136.1_S1283s_1_1_1283_1283   microtubule-associated serine/threonine-prote...   \n",
       "NP_001158136.1_S1309s_1_1_1309_1309   microtubule-associated serine/threonine-prote...   \n",
       "NP_001158136.1_S1370s_1_1_1370_1370   microtubule-associated serine/threonine-prote...   \n",
       "NP_001158136.1_S1384s_1_1_1384_1384   microtubule-associated serine/threonine-prote...   \n",
       "\n",
       "                                    geneSymbol  protein_mw accession_number  \\\n",
       "id                                                                            \n",
       "NP_001158136.1_S137s_1_1_137_137         MAST4    286317.4   NP_001158136.1   \n",
       "NP_001158136.1_S1283s_1_1_1283_1283      MAST4    286317.4   NP_001158136.1   \n",
       "NP_001158136.1_S1309s_1_1_1309_1309      MAST4    286317.4   NP_001158136.1   \n",
       "NP_001158136.1_S1370s_1_1_1370_1370      MAST4    286317.4   NP_001158136.1   \n",
       "NP_001158136.1_S1384s_1_1_1384_1384      MAST4    286317.4   NP_001158136.1   \n",
       "\n",
       "                                             feature variableSites  \\\n",
       "id                                                                   \n",
       "NP_001158136.1_S137s_1_1_137_137     phosphoproteome         S137s   \n",
       "NP_001158136.1_S1283s_1_1_1283_1283  phosphoproteome        S1283s   \n",
       "NP_001158136.1_S1309s_1_1_1309_1309  phosphoproteome        S1309s   \n",
       "NP_001158136.1_S1370s_1_1_1370_1370  phosphoproteome        S1370s   \n",
       "NP_001158136.1_S1384s_1_1_1384_1384  phosphoproteome        S1384s   \n",
       "\n",
       "                                                               sequence  \\\n",
       "id                                                                        \n",
       "NP_001158136.1_S137s_1_1_137_137     ['KCsNPDVASGPGK', 'KCsNPDVASGPGK']   \n",
       "NP_001158136.1_S1283s_1_1_1283_1283              ['sFSCLNR', 'sFSCLNR']   \n",
       "NP_001158136.1_S1309s_1_1_1309_1309              ['sPTPSYR', 'sPTPSYR']   \n",
       "NP_001158136.1_S1370s_1_1_1370_1370    ['sAGNIPLSPLAR', 'sAGNIPLSPLAR']   \n",
       "NP_001158136.1_S1384s_1_1_1384_1384  ['TPsPTPQPTSPQR', 'TPsPTPQPTSPQR']   \n",
       "\n",
       "                                                                           sequenceVML  \\\n",
       "id                                                                                       \n",
       "NP_001158136.1_S137s_1_1_137_137     ['KCS(0.99)NPDVAS(0.0)GPGK', 'KCS(0.99)NPDVAS(...   \n",
       "NP_001158136.1_S1283s_1_1_1283_1283       ['S(0.99)FS(0.0)CLNR', 'S(0.99)FS(0.0)CLNR']   \n",
       "NP_001158136.1_S1309s_1_1_1309_1309  ['S(0.99)PT(0.0)PS(0.0)Y(0.0)R', 'S(0.99)PT(0....   \n",
       "NP_001158136.1_S1370s_1_1_1370_1370  ['S(0.99)AGNIPLS(0.0)PLAR', 'S(0.99)AGNIPLS(0....   \n",
       "NP_001158136.1_S1384s_1_1_1384_1384  ['T(0.0)PS(0.99)PT(0.0)PQPT(0.0)S(0.0)PQR', 'T...   \n",
       "\n",
       "                                    ptmSite  \\\n",
       "id                                            \n",
       "NP_001158136.1_S137s_1_1_137_137      S137s   \n",
       "NP_001158136.1_S1283s_1_1_1283_1283  S1283s   \n",
       "NP_001158136.1_S1309s_1_1_1309_1309  S1309s   \n",
       "NP_001158136.1_S1370s_1_1_1370_1370  S1370s   \n",
       "NP_001158136.1_S1384s_1_1_1384_1384  S1384s   \n",
       "\n",
       "                                                                               Hit_def  \\\n",
       "id                                                                                       \n",
       "NP_001158136.1_S137s_1_1_137_137     sp|O15021|MAST4_HUMAN Microtubule-associated s...   \n",
       "NP_001158136.1_S1283s_1_1_1283_1283  sp|O15021|MAST4_HUMAN Microtubule-associated s...   \n",
       "NP_001158136.1_S1309s_1_1_1309_1309  sp|O15021|MAST4_HUMAN Microtubule-associated s...   \n",
       "NP_001158136.1_S1370s_1_1_1370_1370  sp|O15021|MAST4_HUMAN Microtubule-associated s...   \n",
       "NP_001158136.1_S1384s_1_1_1384_1384  sp|O15021|MAST4_HUMAN Microtubule-associated s...   \n",
       "\n",
       "                                     ...  \\\n",
       "id                                   ...   \n",
       "NP_001158136.1_S137s_1_1_137_137     ...   \n",
       "NP_001158136.1_S1283s_1_1_1283_1283  ...   \n",
       "NP_001158136.1_S1309s_1_1_1309_1309  ...   \n",
       "NP_001158136.1_S1370s_1_1_1370_1370  ...   \n",
       "NP_001158136.1_S1384s_1_1_1384_1384  ...   \n",
       "\n",
       "                                                                              Hsp_qseq  \\\n",
       "id                                                                                       \n",
       "NP_001158136.1_S137s_1_1_137_137     MGEKVSEAPEPVPRGCSGHGSRTPASALVAASSPGASSAESSSGSE...   \n",
       "NP_001158136.1_S1283s_1_1_1283_1283  MGEKVSEAPEPVPRGCSGHGSRTPASALVAASSPGASSAESSSGSE...   \n",
       "NP_001158136.1_S1309s_1_1_1309_1309  MGEKVSEAPEPVPRGCSGHGSRTPASALVAASSPGASSAESSSGSE...   \n",
       "NP_001158136.1_S1370s_1_1_1370_1370  MGEKVSEAPEPVPRGCSGHGSRTPASALVAASSPGASSAESSSGSE...   \n",
       "NP_001158136.1_S1384s_1_1_1384_1384  MGEKVSEAPEPVPRGCSGHGSRTPASALVAASSPGASSAESSSGSE...   \n",
       "\n",
       "                                                                              Hsp_hseq  \\\n",
       "id                                                                                       \n",
       "NP_001158136.1_S137s_1_1_137_137     MGEKVSEAPEPVPRGCSGHGSRTPASALVAASSPGASSAESSSGSE...   \n",
       "NP_001158136.1_S1283s_1_1_1283_1283  MGEKVSEAPEPVPRGCSGHGSRTPASALVAASSPGASSAESSSGSE...   \n",
       "NP_001158136.1_S1309s_1_1_1309_1309  MGEKVSEAPEPVPRGCSGHGSRTPASALVAASSPGASSAESSSGSE...   \n",
       "NP_001158136.1_S1370s_1_1_1370_1370  MGEKVSEAPEPVPRGCSGHGSRTPASALVAASSPGASSAESSSGSE...   \n",
       "NP_001158136.1_S1384s_1_1_1384_1384  MGEKVSEAPEPVPRGCSGHGSRTPASALVAASSPGASSAESSSGSE...   \n",
       "\n",
       "                                                                           Hsp_midline  \\\n",
       "id                                                                                       \n",
       "NP_001158136.1_S137s_1_1_137_137     MGEKVSEAPEPVPRGCSGHGSRTPASALVAASSPGASSAESSSGSE...   \n",
       "NP_001158136.1_S1283s_1_1_1283_1283  MGEKVSEAPEPVPRGCSGHGSRTPASALVAASSPGASSAESSSGSE...   \n",
       "NP_001158136.1_S1309s_1_1_1309_1309  MGEKVSEAPEPVPRGCSGHGSRTPASALVAASSPGASSAESSSGSE...   \n",
       "NP_001158136.1_S1370s_1_1_1370_1370  MGEKVSEAPEPVPRGCSGHGSRTPASALVAASSPGASSAESSSGSE...   \n",
       "NP_001158136.1_S1384s_1_1_1384_1384  MGEKVSEAPEPVPRGCSGHGSRTPASALVAASSPGASSAESSSGSE...   \n",
       "\n",
       "                                    relIdentity uniprot species sifts blast  \\\n",
       "id                                                                            \n",
       "NP_001158136.1_S137s_1_1_137_137            1.0  O15021   HUMAN  True  True   \n",
       "NP_001158136.1_S1283s_1_1_1283_1283         1.0  O15021   HUMAN  True  True   \n",
       "NP_001158136.1_S1309s_1_1_1309_1309         1.0  O15021   HUMAN  True  True   \n",
       "NP_001158136.1_S1370s_1_1_1370_1370         1.0  O15021   HUMAN  True  True   \n",
       "NP_001158136.1_S1384s_1_1_1384_1384         1.0  O15021   HUMAN  True  True   \n",
       "\n",
       "                                    acc_res acc_res_i  \n",
       "id                                                     \n",
       "NP_001158136.1_S137s_1_1_137_137          S       137  \n",
       "NP_001158136.1_S1283s_1_1_1283_1283       S      1283  \n",
       "NP_001158136.1_S1309s_1_1_1309_1309       S      1309  \n",
       "NP_001158136.1_S1370s_1_1_1370_1370       S      1370  \n",
       "NP_001158136.1_S1384s_1_1_1384_1384       S      1384  \n",
       "\n",
       "[5 rows x 34 columns]"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "def grab_ptm_site(x):\n",
    "    \"\"\"Grab actual PTM site.\"\"\"\n",
    "    res = list()\n",
    "    \n",
    "    for s in x.split(\" \"):\n",
    "        if s.startswith((\"K\",\"Y\",\"S\",\"T\")):\n",
    "            res.append(s)\n",
    "        else:\n",
    "            pass\n",
    "        \n",
    "    if len(res) == 1:\n",
    "        return res[0]\n",
    "    else:\n",
    "        return res\n",
    "\n",
    "# Expand PTMs in Each Feature\n",
    "ptm_df['variableSites'] = ptm_df['variableSites'].apply(ast.literal_eval)\n",
    "ptm_df = ptm_df.explode('variableSites').drop_duplicates()\n",
    "\n",
    "# Proteins with single PTM sites\n",
    "ptm_sing_df = ptm_df.loc[ptm_df.index.map(lambda x: x.split(\"_\")[-3]==\"1\" and x.split(\"_\")[-4]==\"1\")]\n",
    "ptm_sing_df.loc[:,\"ptmSite\"] = ptm_sing_df['variableSites'].apply(grab_ptm_site)\n",
    "\n",
    "# Proteins with multiple PTM sites\n",
    "ptm_multi_df = pd.concat((\n",
    "    ptm_df.loc[ptm_df.index.map(lambda x: x.split(\"_\")[-3]==\"3\" and x.split(\"_\")[-4]==\"3\")],\n",
    "    ptm_df.loc[ptm_df.index.map(lambda x: x.split(\"_\")[-3]==\"2\" and x.split(\"_\")[-4]==\"2\")]\n",
    "))\n",
    "\n",
    "ptm_multi_df.loc[:,\"ptmSite\"] = ptm_multi_df['variableSites'].apply(grab_ptm_site)\n",
    "ptm_multi_df = ptm_multi_df.explode(\"ptmSite\")\n",
    "\n",
    "# Combine & intersect accession numbers\n",
    "ptm_comb_df = pd.concat((ptm_sing_df, ptm_multi_df))\n",
    "ptm_comb_df = ptm_comb_df[ptm_comb_df[\"accession_number\"].isin(mapped_acc_df.index)]\n",
    "\n",
    "# Join PTM features from source data-set to mapped accession-numbers\n",
    "ptm_comb_df = pd.merge(\n",
    "    ptm_comb_df.reset_index(), \n",
    "    mapped_acc_df.reset_index().rename(columns={\"query\":\"accession_number\"}),\n",
    "    how=\"left\",\n",
    "    on=\"accession_number\"\n",
    ").set_index(\"id\")\n",
    "\n",
    "ptm_comb_df[\"acc_res\"] = ptm_comb_df[\"ptmSite\"].apply(lambda x: x[0])\n",
    "ptm_comb_df[\"acc_res_i\"] = ptm_comb_df[\"ptmSite\"].apply(lambda x: int(x[1:-1]))\n",
    "\n",
    "print(\"  * {} single PTM sites total in CPTAC dataset\".format(ptm_sing_df.shape[0]))\n",
    "print(\"  * {} multi PTM sites total in CPTAC dataset\".format(ptm_multi_df.shape[0]))\n",
    "\n",
    "ptm_comb_df.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Filter\n",
    "ptm_comb_filt_df = ptm_comb_df[ptm_comb_df['Hsp_query-from'].notna()].copy()\n",
    "ptm_comb_filt_df['Hsp_query-from'] = ptm_comb_filt_df['Hsp_query-from'].astype(int)\n",
    "ptm_comb_filt_df['Hsp_hit-from'] = ptm_comb_filt_df['Hsp_hit-from'].astype(int)\n",
    "\n",
    "def get_source_blast_idx(row):\n",
    "    \"\"\"\n",
    "    Get source blast index.\n",
    "    \"\"\"\n",
    "    d = dict()\n",
    "    c = row['Hsp_query-from']\n",
    "\n",
    "    for idx,res in enumerate(row['Hsp_qseq']):\n",
    "        if res is not '-':\n",
    "            d[c] = idx\n",
    "            c = c+1\n",
    "    \n",
    "    if row['acc_res_i'] not in d.keys():\n",
    "        return \"X\"\n",
    "    else:\n",
    "        return d[row['acc_res_i']]\n",
    "\n",
    "def get_uniprot_blast_i(row):\n",
    "    \"\"\"\n",
    "    Get Uniprot Mapped Residue ID.\n",
    "    \"\"\"\n",
    "    d = dict()\n",
    "    c = row['Hsp_hit-from']\n",
    "\n",
    "    for idx,res in enumerate(row['Hsp_hseq']):\n",
    "        if res is not '-':\n",
    "            d[idx] = c\n",
    "            c = c+1\n",
    "\n",
    "    if row['acc_res_idx'] not in d.keys():\n",
    "        return \"X\"\n",
    "    else:\n",
    "        return d[row['acc_res_idx']]\n",
    "\n",
    "def _get_qres(row):\n",
    "    \"\"\"Get query residue\"\"\"\n",
    "    try:\n",
    "        return row[\"Hsp_qseq\"][row[\"acc_res_idx\"]]\n",
    "    except:\n",
    "        return None\n",
    "\n",
    "def _get_hres(row):\n",
    "    \"\"\"Get query residue\"\"\"\n",
    "    try:\n",
    "        return row[\"Hsp_hseq\"][row[\"acc_res_idx\"]]\n",
    "    except:\n",
    "        return None\n",
    "\n",
    "ptm_comb_filt_df[\"acc_res_idx\"] = ptm_comb_filt_df.apply(get_source_blast_idx, 1)\n",
    "ptm_comb_filt_df[\"acc_res\"] = ptm_comb_filt_df.apply(_get_qres, 1)\n",
    "\n",
    "# Get Uniprot Matches\n",
    "ptm_comb_filt_df[\"uniprot_res\"] = ptm_comb_filt_df.apply(_get_hres, 1)\n",
    "ptm_comb_filt_df[\"uniprot_res_i\"] = ptm_comb_filt_df.apply(get_uniprot_blast_i, 1)\n",
    "ptm_comb_filt_df['uniprot_match'] = ptm_comb_filt_df['acc_res']==ptm_comb_filt_df['uniprot_res']"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 7. Align Source PTMs to PDB Residues via SIFTS\n",
    "---"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [],
   "source": [
    "sifts_filt_df['PDB_BEG'] = sifts_filt_df['PDB_BEG'].apply(lambda x: 0 if x==\"None\" else x)\n",
    "sifts_filt_df['len'] = sifts_filt_df['RES_END'] - sifts_filt_df['RES_BEG']\n",
    "sifts_drop_dup = sifts_filt_df.drop_duplicates(subset=['PDB','len'])\n",
    "\n",
    "ptm_comb_filt2_df = ptm_comb_filt_df[ptm_comb_filt_df[\"uniprot_match\"]].copy()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "  2%|▏         | 12/542 [00:00<00:04, 119.27it/s]"
     ]
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "   * Grabbing header info for 542 pdb-chains\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "100%|██████████| 542/542 [00:09<00:00, 59.59it/s]\n"
     ]
    }
   ],
   "source": [
    "def get_pdb_headers(sifts_df, pdbstore):\n",
    "    \"\"\"Get PDB headers.\"\"\" \n",
    "    \n",
    "    pdb_chain_unis = [x for x in zip(sifts_df['PDB'], sifts_df['CHAIN'], sifts_df['SP_PRIMARY'])]\n",
    "    print(\"   * Grabbing header info for {} pdb-chains\".format(len(pdb_chain_unis)))\n",
    "    headers = list()\n",
    "    \n",
    "    def get_header_info(pdb_ch_uni):\n",
    "        \"\"\"Download pdb.\"\"\"\n",
    "        pdb_headers = dict()\n",
    "        \n",
    "        try:\n",
    "            pdb, chain, uniprot = pdb_ch_uni\n",
    "            hd = pdbstore.load_header(pdb)\n",
    "\n",
    "            for polymer in hd['polymers']:\n",
    "                pdb_headers[polymer] = dict()\n",
    "                if polymer.chid == chain:\n",
    "                    for dbref in polymer.dbrefs:\n",
    "                        pdb_headers[polymer][\"PDB\"] = pdb\n",
    "                        pdb_headers[polymer][\"CHAIN\"] = chain\n",
    "                        pdb_headers[polymer][\"db_ref\"] = dbref.database\n",
    "                        pdb_headers[polymer][\"db_accession\"] = dbref.accession\n",
    "                        pdb_headers[polymer][\"db_first_from\"] = dbref.first[0]\n",
    "                        pdb_headers[polymer][\"db_first_to\"] = dbref.first[2]\n",
    "\n",
    "                        if dbref.accession == uniprot:\n",
    "                            continue\n",
    "\n",
    "            headers.append(pd.DataFrame.from_dict(pdb_headers).dropna(1))\n",
    "        except:\n",
    "            print(\"  * Error for {}\".format(pdb_ch_uni))\n",
    "    \n",
    "    for pdb_chain_uni in tqdm(pdb_chain_unis):\n",
    "        get_header_info(pdb_chain_uni)\n",
    "        \n",
    "    return headers\n",
    "\n",
    "pdb_headers = pd.concat(get_pdb_headers(sifts_drop_dup, pdbstore),1).T\n",
    "pdb_headers.to_csv(os.path.join(REF_DIR, \"pdb_headers.txt\"), sep='\\t')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "db_match\n",
       "False     61\n",
       "True     492\n",
       "dtype: int64"
      ]
     },
     "execution_count": 19,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "pdb_headerpdb_headers = pd.read_csv(os.path.join(REF_DIR, \"pdb_headers.txt\"), sep='\\t', index_col=0)\n",
    "\n",
    "sifts_drop_dup = pd.merge(sifts_drop_dup.reset_index(), pdb_headers.reset_index().rename(columns={'index':'pdb_desc'}))\n",
    "sifts_drop_dup['db_match'] = sifts_drop_dup['SP_PRIMARY']==sifts_drop_dup['db_accession']\n",
    "sifts_drop_dup.groupby('db_match').size()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "['NP_000057.2' 'NP_000329.2' 'NP_001008240.1' 'NP_001129976.1'\n",
      " 'NP_001137407.1' 'NP_001158136.1' 'NP_001159737.1' 'NP_001171908.1'\n",
      " 'NP_001171915.1' 'NP_001254507.1' 'NP_001269115.1' 'NP_001284484.1'\n",
      " 'NP_001310876.1' 'NP_001336478.1' 'NP_001616.1' 'NP_001761.3'\n",
      " 'NP_003089.1' 'NP_005151.2' 'NP_005316.1' 'NP_005321.1' 'NP_056173.1'\n",
      " 'NP_066568.3' 'NP_542192.2' 'NP_940947.3' 'NP_942133.1']\n"
     ]
    }
   ],
   "source": [
    "accession_numbers_to_use = np.intersect1d(\n",
    "    mapped_acc_df[\n",
    "        (mapped_acc_df['blast']) & (mapped_acc_df['sifts'])\n",
    "    ].index,\n",
    "    ptm_comb_filt2_df['accession_number']\n",
    ")\n",
    "\n",
    "print(accession_numbers_to_use)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 7.1 PDB Matches\n",
    "---"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "CPU times: user 7.55 s, sys: 2.96 s, total: 10.5 s\n",
      "Wall time: 9.71 s\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "{None}"
      ]
     },
     "execution_count": 21,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "%%time\n",
    "\n",
    "clumpsptm.mp.get_pdb_matches(\n",
    "    accession_numbers_to_use,\n",
    "    ptm_comb_filt2_df,\n",
    "    sifts_drop_dup,\n",
    "    pdbstore,\n",
    "    os.path.join(REF_DIR, \"pdb_matches\"),\n",
    "    4\n",
    ")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 8. Select Best PDB Structure\n",
    "---\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      " 16 / 25 (64.00 %) structures with matching sites\n"
     ]
    }
   ],
   "source": [
    "pdb_mapped_df = pd.concat([\n",
    "    pd.read_parquet(x) for x in \n",
    "    glob.glob(os.path.join(REF_DIR, \"pdb_matches/acc_mapped_sifts/*\"))\n",
    "])\n",
    "\n",
    "pdb_mapped_df['db_match'] = pdb_mapped_df['db_match'].astype(bool)\n",
    "pdb_mapped_df = pdb_mapped_df[pdb_mapped_df['db_match']]\n",
    "pdb_mapped_df = pdb_mapped_df[pdb_mapped_df['db_ref']=='UniProt']\n",
    "pdb_mapped_df['percent_match'] = pdb_mapped_df['percent_match'].astype(float)\n",
    "pdb_mapped_df = pdb_mapped_df.sort_values(['percent_match','len','proteins'], ascending=False)\n",
    "\n",
    "# Select Top Structure by overlapping matched sites\n",
    "pdb_mapped_filt_df = pdb_mapped_df.drop_duplicates(subset=['proteins'])\n",
    "\n",
    "_prev = pdb_mapped_filt_df.shape[0]\n",
    "pdb_mapped_filt_df = pdb_mapped_filt_df[pdb_mapped_filt_df['percent_match']>0]\n",
    "pdb_mapped_filt_df['len'] = pdb_mapped_filt_df['len'].astype(int)\n",
    "print(\" {} / {} ({:.2f} %) structures with matching sites\".format(\n",
    "    pdb_mapped_filt_df.shape[0], _prev, 100*pdb_mapped_filt_df.shape[0]/_prev))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Text(0, 0.5, 'PDB Length')"
      ]
     },
     "execution_count": 25,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAVcAAADUCAYAAAAoYmDaAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+j8jraAAAgAElEQVR4nO3de5xdVX338c93ZnIHTDBBAwEDdqrFWC8dgdZSsfESfFQs1gvewAcNraCpxRcFtAJare2Dl6gUihXFl4qgjYUHaRAjFB6rSKBWBxCTAoGRkQRJMDHXmfk9f6x1yJ6Tc2bOXM6cS77v1+u8ztlrr7P32nuf8zvrrL3X2ooIzMxscnU0ugBmZu3IwdXMrA4cXM3M6sDB1cysDhxczczqwMHVzKwOuhpdgKkwf/78WLx4caOLYWZt5s4773wsIhZUmrdfBNfFixezdu3aRhfDzNqMpA3V5rlZwMysDhxczczqYL9oFmg1q3v7WblmHX2bd7Bo3ixWLO1m2ZKFjS6WmY2Bg2uTWd3bz/mrehmKYGZXB/1bdnL+ql4AB1izFuJmgSazcs06hiKY3tWBJKZ3dTAUwco16xpdNDMbAwfXJtO3eQfTOjUsbVqn6Nu8o0ElMrPxcHBtMovmzWLP4PBhIPcMBovmzWpQicxsPBxcm8yKpd10SOweGCIi2D0wRIfEiqXdjS6amY2BT2g1mdJJK18tYNbaHFyb0LIlCx1MzVqcmwXMzOrAwdXMrA4cXM3M6sDB1cysDhxczczqwMHVzKwOHFzNzOrAwdXMrA4cXM3M6sDB1cysDhxczczqwMHVzKwOHFzNzOrAwdXMrA6mLLhKOlzSzZLulXS3pBU5/WBJN0lal5/n5XRJ+qyk9ZJ+KumFhWWdmvOvk3TqVG2DmVmtprLmOgCcHRG/BxwHnCnpaOBcYE1EdANr8jTAiUB3fiwHLoUUjIELgGOBY4ALSgHZzKxZTFlwjYj+iLgrv94K3AscBpwEXJmzXQm8Lr8+CfhKJD8C5kpaCLwSuCkiHo+IzcBNwLKp2g4zs1o0pM1V0mLgBcDtwNMioh9SAAYOydkOAx4uvK0vp1VLL1/HcklrJa3dtGnTZG+CmdmIpjy4SjoA+FfgryLiNyNlrZAWI6QPT4i4PCJ6IqJnwYIF4yusmdk4TWlwlTSNFFi/FhGrcvKj+e8++XljTu8DDi+8fRHwyAjpZmZNYyqvFhDwReDeiPhUYdZ1QOmM/6nAtYX0d+SrBo4DnsjNBjcCr5A0L5/IekVOMzNrGlN599cXA28HfibpJzntfOATwDWSTgceAt6Q590AvApYD2wH3gkQEY9L+ihwR873kYh4fGo2wcysNorYp7my7fT09MTatWsbXQwzazOS7oyInkrz3EPLzKwOHFzNzOrAwdXMrA4cXM3M6sDB1cysDhxczczqwMHVzKwOHFzNzOqg5h5ako4FlpJGrRoWlCPifZNcLjOzllZTcJX0AeAfSV1RH2H4KFTt38XLzGyMaq25rgDeFxGfr2dhzMzaRa1trgeRBlIxM7Ma1Bpcr8K3UjEzq1nVZgFJf12YfBi4SNKLgZ8Ce4p5y8ZnNTPb743U5vresultwB/lR1EADq5mZgVVg2tEHDmVBTEzayc1tblK+rCk2RXSZ0n68OQXy8ystdV6QusC4IAK6bPzPDMzK6g1uIrKnQVeAPj+VWZmZUbsRCBpKymoBnC/pGKA7QRmApfVr3hmZq1ptB5aZ5FqrVcAHwSeKMzbDTwYET+sU9nMzFrWiME1Iq4EkPQA8J8RsWek/GZmltTa5voAsFDSEWWPwyUtqGUBkq6QtFFSbyHtQkm/lPST/HhVYd55ktZLuk/SKwvpy3Laeknn1rqhZmZTqdaBWx5khNGvJP0G+BJwTkQMVMn2ZeDzwFfK0j8dEReXLe9o4M3Ac4BDge9J+t08+xLg5UAfcIek6yLinhq3w8xsStQaXE8hDTl4GXB7TjsWWA5cCMwFPgRspcqlWRFxq6TFNa7vJOAbEbELeEDSeuCYPG99RNwPIOkbOa+Dq5k1lVqD618C74+IVYW070u6D1gRES+RtBG4iLFf93qWpHcAa4GzI2IzcBjwo0KevpwGaZyDYvqxY1yfmVnd1drmeizwswrpvcCL8usfAovGuP5LgWcCzwf6gU/mdFXIGyOk70PScklrJa3dtGnTGItlZjYxtQbXDaQmgHLvBh7Krxcwxg4FEfFoRAxGxBDwBfb+9e8DDi9kXUS6A0K19ErLvjwieiKiZ8GCms65mZlNmlqbBc4G/jWfzb+DVFt8EanW+fqc50XANWNZuaSFEdGfJ/+MVBMGuA74uqRPkU5odQM/JtVcuyUdCfySdNLrLWNZp5nZVKgpuEbEdyR1A+8BnkUKctcBl0XEQznPP420DElXAScA8yX1kdpmT5D0fFKwfhA4Iy/rbknXkE5UDQBnRsRgXs5ZwI2kHmJXRMTdY9lgM7OpoIj2v79gT09PrF27ttHFMLM2I+nOiOipNG8st9aeTTrxVOnW2qsqvsnMbD9V6621X0a6j9ZTK8wO0l90MzPLar1aYCXwHWBRRHSUPRxYzczK1NossBh4bURUvOzJzMyGq7Xm+gPSVQJmZlaDWmuulwEXSzqU1FOr/Nbad012wczMWlmtwfVb+fnyCvN8QsvMrEytwdW32TYzG4Nae2htqHdBzMzaSa0ntJB0oqTrJd0j6fCc9i5JS+tXPDOz1lRrJ4K3kk5q/QuwFJiWZ3UC5wBr6lI6M6uL1b39rFyzjr7NO1g0bxYrlnazbMnCRherrdRacz0HeHdEvJ80kErJj0hdYs2sRazu7ef8Vb30b9nJzK4O+rfs5PxVvazu7R/9zVazWoNrN2kw7HLbgIMmrzhmVm8r16xjKILpXR1IYnpXB0MRrFyzrtFFayu1BtdHgN+tkP4nwP9MXnHMrN76Nu9gWufwm3pM6xR9m3c0qETtqdbgejnwWUkvztOHSzqVdNPCS+tSMjOri0XzZrFncPhQo3sGg0XzZjWoRO2ppuAaEf8IrAJuAuYAN5NOcF0WEZfUr3hmNtlWLO2mQ2L3wBARwe6BITokViztbnTR2krN47lGxAclfQw4mhSU7wH2SDo6Inxra7MWUboqwFcL1FfNwRUgIraTboENgKTnAXfh7q9mLWXZkoUOpnVWcycCMzOrnYOrmVkdOLiamdXBiG2ukl44yvs9gLaZWQWjndBaSxqvVSPkaf97c5uZjdFozQJHAkfl52qPo2pZkaQrJG2U1FtIO1jSTZLW5ed5OV2SPitpvaSfFmvQkk7N+dfljgxmZk1nxOAaERtqedS4ri8Dy8rSzgXWREQ3aWStc3P6iaTxDLqB5eReYJIOBi4AjgWOAS4oBWQzs2YyZSe0IuJW4PGy5JOAK/PrK4HXFdK/EsmPgLmSFgKvBG6KiMcjYjOpx1h5wDYza7hGXy3wtIjoB8jPh+T0w4CHC/n6clq1dDOzptLo4FpNpRNo1U6sVTyhJmm5pLWS1m7atGlSC2dmNppGB9dH89998vPGnN4HHF7It4g07GG19H1ExOUR0RMRPQsWLJj0gpuZjaTRwfU6oHTG/1Tg2kL6O/JVA8cBT+RmgxuBV0ial09kvSKnmZk1lVGDq6RjJJ0h6QV5+nRJGyRtkvQFSTUNAinpKtLdDJ4lqU/S6cAngJdLWge8PE8D3ADcD6wHvgC8ByAiHgc+CtyRHx/JaWZmTWW0HlpvJ11C9TjwFEnnkILbV4Eh4G3Ar9l7CVVVEXFKlVn73D02IgI4s8pyrgCuGG19ZmaNNFrN9QPAeRGxAHgj6c4DZ0XEX0bEmcC7gD+vcxnNzFrOaMH1d4Br8utSe+iPC/NvZ/gJJjMzY/TgOgvYDk/+Vd+VHyW7GOOA22Zm+4PRAmMA8yQNkK4xDVJvqYPz/IOrvtPMbD82WnAV6V5Zxek7yqY9KpaZWZnRgutLp6QUZmZtZsTgGhH/MVUFMTNrJzWdjJJ0EPAy0titQbrA/3sRsbWOZTMza1mjBldJbwQuA+aWzdoi6YyI+GZdSmZm1sJG66H1XFJvrG8CF5NObgl4DqmDwVcl/Twiflbvgprtb1b39rNyzTr6Nu9g0bxZrFjazbIlCxtdLKvRaDXX9wH/HhFvLUu/EzhF0rU5z7vrUTiz/dXq3n7OX9XLUAQzuzro37KT81elOyQ5wLaG0ToRHE++xUoVl+Y8ZjaJVq5Zx1AE07s6kMT0rg6GIli5Zl2ji2Y1Gq3mehhw7wjzf04aU9XMJlHf5h0QweO/HWBwKOjsELOmdaR0awmj1VznACMdzR2kLrJmNokOmtnFtl2DDA6lPjqDQ8G2XYMcNNO9zVtFLUfq9yVVGzN1/mQWxswSoSfvayRBROk+R5XudGTNqJbgeiOV711V4u6vZpPsiZ17OGhmF9t3Dz7ZLDB7eidP7NzT6KJZjUYLrkdOSSnMbJhF82bRv2UnB8+Z/mTa7oEhFs6d2cBS2ViM1v11w1QVxMz2WrG0m/NX9bJ7YIhpnWLPYNAhsWJpd6OL1pbqcU3xiCe0JM2WdImkX0raKOnrktzOalZny5Ys5OMnL2Hh3JnszDXWj5+8xNe41kHpmuL+LTuHXVO8urd/QssdrVngIuA04GvATuAU0rWtb5jQWs1sVMuWLHQwnQLFa4oBpneJ3QNDrFyzbkL7f7TgejJwekR8A0DSV4EfSOqMiMFxr9XMrEn0bd7BzK7hf+KndWrC1xSPFlwPB24rTUTEj/NdCQ4FHp7Qms3MmsCiebPY8Nhv2T0YT16ZMb1TPGP+nAktd7ROBJ3A7rK0AXzfLDNrE8d3z2fHniEGcoeNgaFgx54hju+e2OmlWm7z8lVJxZsSzgS+IGl7KSEiXjuRQkh6ENgKDAIDEdGT79N1NbAYeBB4Y0RsliRgJfAq0s0TT4uIuyayfjPbf9227jFmT+9k18AQg0NBV4eY0dXBbesem9ByRwuuV1ZI++qE1ljdSyOiuDXnAmsi4hOSzs3TfwOcCHTnx7GkE2zH1qlMZtbm+jbvYPb0TubM2BsOI6K+ba4R8c4JLX1iTgJOyK+vBG4hBdeTgK/kW33/SNJcSQsjYmLXTZjZfqnUYWN6196OqHsGg0XzJjZsymhtrlMlgO9KulPS8pz2tFLAzM+H5PTDGH4yrS+nmZmN2Yql3XQoXX4VEeweGJqUDhvNcmLqxRHxiKRDgJsk/XyEvJXGOdhnfIMcpJcDHHHEEZNTSjNrO6VrWSe7h1ZTBNeIeCQ/b5T0beAY4NHS331JC4GNOXsf6RKxkkXAIxWWeTlwOUBPT48HlzGzqurRYaPhzQKS5kg6sPQaeAXQC1wHnJqznQpcm19fB7xDyXHAE25vNbNm0ww116cB305XWNEFfD0iVku6A7hG0unAQ+ztcnsD6TKs9aRLsRp50s3MrKKGB9eIuB94XoX0XwNLK6QHcOYUFM3MbNwa3ixgZtaOHFzNzOrAwdXMrA4cXM3M6qDhJ7T2Z/W4tYQ1no+rgYNrw5RuLTEUMezWEoC/iC3Mx9VK3CzQIMVbS0hielcHQxGsXLOu0UWzCfBxtRLXXBukXreWsMbycW1NU373V6ufRfNmsWdw+JAHkzHMmTWWj2vrWd3bz9nX/Df3/Wor23YOcN+vtnL2Nf894bu/OriWWd3bz4krb+W5F97IiStvnfAOrqZew5xZY/m4tp6PXn8P23cPEvk3MQK27x7ko9ffM6HlulmgYCpPRtRrmDNrLB/X1tP/xE4AVBjMNGJv+ng5uBbU6/7l1fi+9O3Jx7X1BMMHii6fHg83CxT0bd7BtM7hu9QnI8za26FPSe3hQ7H3UUwfLwfXAp+MMNv/vOr3nz6m9Fo5uBb4ZITZ/ue2dY8xZ3onXR1CQFeHmDO9s+631t6v+GSE2f6nIbfW3h/5ZITZ/qXdb61tZtYQK5Z2s2dwiMe27WLT1l08tm0XewaHJtwc6OBqZvs9KV16VboESxO9DgsHVzPbz61cs46ujg6eesAMDjlwBk89YAZdHR0THmzHba5mDeSxXxuvXoPtOLiaNYjHfm0Oi+bN4qFfb2fXwBCDQ0Fnh5jR1cERT509oeW6WcDa3lQNxjNWHvu1ORzfPZ/tuwcZyF2zBoaC7bsHOb57/oSW27LBVdIySfdJWi/p3EaXx5pTqXbYv2XnsNphMwRYd7duDrete4xZ0zro6kjHoqtDzJrWMeFOBC0ZXCV1ApcAJwJHA6dIOrqxpbJm1My1Q3e3bg59m3cwZ0YXB8+ZzoIDZ3DwnOnMmdE14R+5lgyuwDHA+oi4PyJ2A98ATmpwmawJNXPt0N2tm0O9fuRaNbgeBjxcmO7LaU+StFzSWklrN23aNKWFs+bRzLXDZUsW8vGTl7Bw7kx2DgyxcO5MPn7yEp/MmmL1+pFr1asFKl3iO+wbFBGXA5cD9PT0RIX8th9YsbSb81f1sntgiGmdYs9gNFXt0N2tG69eY4q0anDtAw4vTC8CHmlQWayJeTAeq0U9fuRaNbjeAXRLOhL4JfBm4C2NLZI1K9cOrRFaMrhGxICks4AbgU7gioi4u8HFMjN7UksGV4CIuAG4odHlMDOrpFWvFjAza2qKaP8T6ZI2ARsaXY4Jmg9MrMtIc2rX7YL23bZ23S4Y+7Y9IyIWVJqxXwTXdiBpbUT0NLock61dtwvad9vadbtgcrfNzQJmZnXg4GpmVgcOrq3j8kYXoE7adbugfbetXbcLJnHb3OZqZlYHrrmamdWBg2uTkXS4pJsl3SvpbkkrcvrBkm6StC4/z2t0WcdDUqek/5J0fZ4+UtLtebuuljS90WUcD0lzJX1L0s/zsfvDNjpm78+fxV5JV0ma2YrHTdIVkjZK6i2kVTxGSj6bB+P/qaQXjnV9Dq7NZwA4OyJ+DzgOODMPBH4usCYiuoE1eboVrQDuLUz/A/DpvF2bgdMbUqqJWwmsjohnA88jbWPLHzNJhwHvA3oiYgmpu/mbac3j9mVgWVlatWN0ItCdH8uBS8e8tojwo4kfwLXAy4H7gIU5bSFwX6PLNo5tWZQ/wH8KXE8aOvIxoCvP/0PgxkaXcxzbdRDwAPkcRiG9HY5Zaezkg0nd5a8HXtmqxw1YDPSOdoyAfwZOqZSv1odrrk1M0mLgBcDtwNMioh8gPx/SuJKN22eAc4ChPP1UYEtEDOTpfQY9bxFHAZuAL+Umj3+RNIc2OGYR8UvgYuAhoB94AriT9jhuUP0YjTog/2gcXJuUpAOAfwX+KiJ+0+jyTJSkVwMbI+LOYnKFrK14+UoX8ELg0oh4AfBbWrAJoJLcBnkScCRwKDCH9Je5XCset5FM+LPp4NqEJE0jBdavRcSqnPyopIV5/kJgY6PKN04vBl4r6UHSPc/+lFSTnSupNDpbqw563gf0RcTtefpbpGDb6scM4GXAAxGxKSL2AKuAP6I9jhtUP0YTHpDfwbXJSBLwReDeiPhUYdZ1wKn59amkttiWERHnRcSiiFhMOiHy/Yh4K3Az8Oc5W8ttF0BE/Ap4WNKzctJS4B5a/JhlDwHHSZqdP5ulbWv545ZVO0bXAe/IVw0cBzxRaj6olTsRNBlJfwzcBvyMvW2T55PaXa8BjiB94N8QEY83pJATJOkE4AMR8WpJR5FqsgcD/wW8LSJ2NbJ84yHp+cC/ANOB+4F3kiovLX/MJF0EvIl0Jct/Ae8itT+21HGTdBVwAmnkq0eBC4B/o8Ixyj8knyddXbAdeGdErB3T+hxczcwmn5sFzMzqwMHVzKwOHFzNzOrAwdXMrA4cXM3M6sDB1ZqWpFskfb4B6z1N0rapXm87k3RhcTSqKnk+L+mWKSpS3Tm4jpGkL0uK/Ngj6X5JF+e+5E2tlg94Ie9MSX+bh8/bKelxSddLOrbe5WwCV5PGC6grSYsLn6WQtFnSrZJeUshT/nnbqDQk5Zm5J19xebcU8g5J+pWkr5d6II1QjuI6BiQ9JOnSSR4i8WLgJaPmaiMOruPzPdIIOkcBHwLeQ/rwjEv5l6TR8tic3wX+Avg74FmknjkbgdskvWYK1t8wEbEjIqayq+oy0ufpJcBvgBskHVmYX/q8LQZeAfxf4CLSsSj/Uf9SzrsIOBk4GriihjIU1/Eu4DXAP41rayqIiG0R8evJWl5LaPQQYK32II0JeX1Z2heA/vxapJGf/gfYQepp9bZC3sWkASBOAb6f85yV552a8+8i9SD5cuF9TyHd32cjsBX4D9IYm6X5pwHbSEGwlzR4yM3AkYX5UfY4rco2lkauemGFef+WyzAb+N28nOeW5VlOGpJuWp4+GvhOLvdG4Crg6eX7FPgbUp/ujTn9FuDzhXxvA+4oLOebwGGF+Sfk8rwa+AmwkzSC0x9U2E+vAX6R89wMHFWepzB9Yd6nb87HdWveD/MLebqAT5PGNt2cX18K3DLCZ6n0WSgex8Ny2hnVPm85fQmwG7iokDZsf+W0s4BHx/GZ/iTw67K0d5K6vu7M++79QEdh/hmFfboJuJG9wxJeyPCh/jpJFZLS/vpM+f6qsj3Dysoo37ec58PABtL36lfAV6YiVrjmOjl2AKXa59+RBg4+kxRU/h74Z0n/q+w9f0+qGRwN/JukM0hjSH4J+H3gVcDd8OR4A98hffFeTRqG8Fbg+2V/+WYA5wH/mzTG5lzgsjzvatIX5j5SDWVhTqvkrcD3IuKuCvP+D7AAeHlE/AJYm/OXv//qiNiTy3crKTgdQxoI5ADgOknFz99L8nYvI/1AVDKd1GXxeXk/zCcF6nIXkwJ1D6kr6nckzS7Mn5GX807SfuoEvp33czWLSV1A/4xUe3wB8LHC/A+QgvK7SIOcdwBvGWF51WzPzyP+m4mIXmA18PpqeSQtAF5H6jpds9wleRmwp5D2buDjpED1e8DZpH38njy/B7iEVKN+Fuk4rx5hNWcD7yYF5NIxKP8c1WLE75uk15OOzXtIA1+/GvjxONYzdlMRwdvpwb6/nMeQamlXk4Zj2wEcX/aezwA35NeLSTWTs8vy9AGfqLLOPyXVtmaVpf8EOCe/Pi0v91mF+W8l1W468vSFFGoPI2zjDmBllXnz8npK611BqhWUulIfTqr1/mGe/ghppPdKyzimsE83ATPK8t1CWc2lbP6z83IW5ekT8vRbC3kOALYA7yrbTy8u5HkGMAi8rJCnvOa6E3hKIe2DwPrCdD9wbmFawM8ZQ801f34uI/Xhf26lz1vZ+z8BbC/bX7vzZ+W3edm9wKE1fKYH8vt2sPefzfsLeR4C3l72vr8C7smvTyaN9XpglXUM++yRRpj6YGG6g1TrvaVse6rWXKnt+/bXpArFtIl+98f6cM11fJZJ2iZpJ/BDUs3svaRfzpnA6jx/Wz7r/JfAM8uW8eQgEJIOIdVK11RZ3x+Q/oZvKlvukrLl7oqI+wrTj5BqQHPHsY2jDTpRmn8VaZzP4/P0W4D7I+KHhbL/SVm5S4MQF8veG6MM/CHphZKulbRB0lb27sMjyrKW1k1EbCP9VTy6MH+IQu0lIjaQ9lUxT7kNEfFEYfoR8sDKkp4CPL1smUFqwqjFrXm/bCU1V5wWET+r4X1i3+N0NfB8Uu3+j0k/fGvy+MAjliG/7xjgc8ANwGfhyRrw4aQaYfE4foK9x/CmvK4HJH1N0qmSDqxY6LS/FjL8OA0xxho2tX3fvpnzPCDpi5LeIGnGGNczLl2jZ7EKbiW1K+4BHok0ziWFkxCvIf3SF+0pm/5t4fVIf0ch/ao/yt4AVlQcSHugbF7pizfWH9FfAM+pMq8UgNYBRMRGSd8j1ZJvzc9fK+TvIDVpfKDCsh4tvP5thflPyidubiSdeHk7qc11PmkEsak4AVZ+/IJ99+t4R0F6C+kHYEuM7aTP0aRmj6InImJ9fr1e0umkWvWbSENZVrO98L73SboZ+FtSjbO0nX8B/GelN0fEVqWb+P0J6bZE5wEfl/SiiBjvWK9D7PvdKDaXlMpV9fsWEaWhIJeSmio+CVwg6diIGPEzN1GuuY7P9ohYHxEbSoE1u4fUaP6MPL/42FBtYRHxKPBLqrc13gU8DRiqsNyxnNXeTWrbGs3XgaWqfMfLc0jNIN8tpH0VeIOkPwCem6eLZX8OqeZXXvatYyj7s0nB9PyIuDUifk7126YcV3qRg/ISht8UsQN4USHPEaTadzFPzXKN9lekWl9pmSquYxR9EfE/YwmskpaQ2kW/NUrWwfw8e8Rc+7oI+BtJhxY+n8+scAxLAZmIGIiI70fEeaT28zmkNs5h8v7qZ/hxEoX9l20i1XCLnld4XdP3LSJ2RsR3IuL9pGPyHNLg7XXlmuskyr/eFwMX5w/LraQ2v+NIgfHyEd7+MeDTkh4l1fRmA0sj4pOk2toPgGslnUNqy3s66cv1vYi4rcYiPgg8IwfNh4CtVf6KfwZ4Lemk0zl53fNIdwE9EfjziNheyP9tUlvhF4EfR8S6wrxLSCcurpb0D6QvzFHAG0ntzrUG2IdIX6SzJF1COqny0Sp5PyRpE+mv+4dJPypfL8wfAD6jdNvyHaQz+3eT9vN4rQTOkfQL0pf+DFJgGNMAy1XMkPR00o/CAtKP8PmkKyHKLwGcnfNC+kH+EKm9+LuMQUTcIulu9l5qeCHwOUlbSE0G00h3WzgsIv5e6TY+zyR95h8HXgocSPUfrJXAeXl//Syvo3x/fZ90nF5Lajc9g9Q88WAu46jfN0mnkeLc7aQ25TeRarXFz2h9THUjb6s/GOEEQ54vUvtr6Vd1E6k96uV5/mLKLr8pvPf0/L7dpJrQFYV5B5I+kH15/sOkwYqfmeefRuEkTE47Ia9rfp6eQarpbGaES7Fy3lmkM+r35e3YQgr6x1XJ/5W8zPdWmNddWO+OvMzPAdNH2qfseynWm0iX3OwktW++Mq/zhLLtfS3w01zuu4AXFZZxGulLdhLpC7aLdFnb75TnKUxfSNmJwAp5ukg/Slvydn4qb9e/j7CPq34Wyj5vpRNMA6R/DbeQPmPTK+yvKDwez2kvHc9nmtRcUaoZQrp88K68/zcD/w94c573x6RL2n6dj3EvaYDpivuQvZeubcmPz7HvpVjTSD/Oj+XHR8rLyujft9eR2na3kF26TPAAAACKSURBVJqe7gBePRWxwoNlW9tQusPBzcCCiHisSp7TSAF7tBM8k1Geu4AfRMR7670uaz5uFjCbBJKeQapJ/wfpe7Wc1D64vJHlssZxcDWbHEPAO0idLDpIf1NPjDHed8nah5sFzMzqwJdimZnVgYOrmVkdOLiamdWBg6uZWR04uJqZ1YGDq5lZHfx/jRQ47+8n3sQAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 360x216 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "fig, ax = plt.subplots(figsize=(5,3))\n",
    "\n",
    "ax.scatter(pdb_mapped_filt_df['percent_match'], pdb_mapped_filt_df['len'], alpha=0.9)\n",
    "ax.set_xlabel(\"Percent Overlapping PDB Residues\", fontsize=14)\n",
    "ax.set_ylabel(\"PDB Length\", fontsize=14)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "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>index</th>\n",
       "      <th>PDB</th>\n",
       "      <th>CHAIN</th>\n",
       "      <th>SP_PRIMARY</th>\n",
       "      <th>RES_BEG</th>\n",
       "      <th>RES_END</th>\n",
       "      <th>PDB_BEG</th>\n",
       "      <th>PDB_END</th>\n",
       "      <th>SP_BEG</th>\n",
       "      <th>SP_END</th>\n",
       "      <th>len</th>\n",
       "      <th>pdb_desc</th>\n",
       "      <th>db_ref</th>\n",
       "      <th>db_accession</th>\n",
       "      <th>db_first_from</th>\n",
       "      <th>db_first_to</th>\n",
       "      <th>db_match</th>\n",
       "      <th>percent_match</th>\n",
       "      <th>proteins</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>66</th>\n",
       "      <td>58470</td>\n",
       "      <td>1xl7</td>\n",
       "      <td>A</td>\n",
       "      <td>Q9DC50</td>\n",
       "      <td>1</td>\n",
       "      <td>612</td>\n",
       "      <td>0</td>\n",
       "      <td>612</td>\n",
       "      <td>1</td>\n",
       "      <td>612</td>\n",
       "      <td>611</td>\n",
       "      <td>PEROXISOMAL CARNITINE O-OCTANOYLTRANSFERASE</td>\n",
       "      <td>UniProt</td>\n",
       "      <td>Q9DC50</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>True</td>\n",
       "      <td>100.000000</td>\n",
       "      <td>NP_001137407.1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>212</th>\n",
       "      <td>174536</td>\n",
       "      <td>3ojy</td>\n",
       "      <td>B</td>\n",
       "      <td>P07358</td>\n",
       "      <td>1</td>\n",
       "      <td>537</td>\n",
       "      <td>0</td>\n",
       "      <td>None</td>\n",
       "      <td>55</td>\n",
       "      <td>591</td>\n",
       "      <td>536</td>\n",
       "      <td>COMPLEMENT COMPONENT C8 BETA CHAIN</td>\n",
       "      <td>UniProt</td>\n",
       "      <td>P07358</td>\n",
       "      <td>1</td>\n",
       "      <td>55</td>\n",
       "      <td>True</td>\n",
       "      <td>100.000000</td>\n",
       "      <td>NP_000057.2</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>364</th>\n",
       "      <td>397057</td>\n",
       "      <td>5xtc</td>\n",
       "      <td>f</td>\n",
       "      <td>O43677</td>\n",
       "      <td>1</td>\n",
       "      <td>47</td>\n",
       "      <td>28</td>\n",
       "      <td>74</td>\n",
       "      <td>28</td>\n",
       "      <td>74</td>\n",
       "      <td>46</td>\n",
       "      <td>NADH DEHYDROGENASE [UBIQUINONE] 1 SUBUNIT C1, ...</td>\n",
       "      <td>UniProt</td>\n",
       "      <td>O43677</td>\n",
       "      <td>28</td>\n",
       "      <td>28</td>\n",
       "      <td>True</td>\n",
       "      <td>100.000000</td>\n",
       "      <td>NP_001171915.1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>74</th>\n",
       "      <td>63125</td>\n",
       "      <td>1z87</td>\n",
       "      <td>A</td>\n",
       "      <td>Q61234</td>\n",
       "      <td>1</td>\n",
       "      <td>263</td>\n",
       "      <td>2</td>\n",
       "      <td>264</td>\n",
       "      <td>2</td>\n",
       "      <td>264</td>\n",
       "      <td>262</td>\n",
       "      <td>ALPHA-1-SYNTROPHIN</td>\n",
       "      <td>UniProt</td>\n",
       "      <td>Q61234</td>\n",
       "      <td>2</td>\n",
       "      <td>2</td>\n",
       "      <td>True</td>\n",
       "      <td>100.000000</td>\n",
       "      <td>NP_003089.1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>84</th>\n",
       "      <td>71308</td>\n",
       "      <td>2c9y</td>\n",
       "      <td>A</td>\n",
       "      <td>P54819</td>\n",
       "      <td>2</td>\n",
       "      <td>240</td>\n",
       "      <td>0</td>\n",
       "      <td>None</td>\n",
       "      <td>1</td>\n",
       "      <td>239</td>\n",
       "      <td>238</td>\n",
       "      <td>ADENYLATE KINASE ISOENZYME 2, MITOCHONDRIAL</td>\n",
       "      <td>UniProt</td>\n",
       "      <td>P54819</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>True</td>\n",
       "      <td>100.000000</td>\n",
       "      <td>NP_001616.1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>277</th>\n",
       "      <td>288151</td>\n",
       "      <td>4wj7</td>\n",
       "      <td>A</td>\n",
       "      <td>Q9BSQ5</td>\n",
       "      <td>3</td>\n",
       "      <td>180</td>\n",
       "      <td>0</td>\n",
       "      <td>None</td>\n",
       "      <td>51</td>\n",
       "      <td>228</td>\n",
       "      <td>177</td>\n",
       "      <td>MALCAVERNIN</td>\n",
       "      <td>UniProt</td>\n",
       "      <td>Q9BSQ5</td>\n",
       "      <td>51</td>\n",
       "      <td>51</td>\n",
       "      <td>True</td>\n",
       "      <td>100.000000</td>\n",
       "      <td>NP_542192.2</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>365</th>\n",
       "      <td>397059</td>\n",
       "      <td>5xtc</td>\n",
       "      <td>h</td>\n",
       "      <td>O43920</td>\n",
       "      <td>1</td>\n",
       "      <td>104</td>\n",
       "      <td>2</td>\n",
       "      <td>105</td>\n",
       "      <td>2</td>\n",
       "      <td>105</td>\n",
       "      <td>103</td>\n",
       "      <td>NADH DEHYDROGENASE [UBIQUINONE] IRON-SULFUR PR...</td>\n",
       "      <td>UniProt</td>\n",
       "      <td>O43920</td>\n",
       "      <td>2</td>\n",
       "      <td>2</td>\n",
       "      <td>True</td>\n",
       "      <td>100.000000</td>\n",
       "      <td>NP_001171908.1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>547</th>\n",
       "      <td>625427</td>\n",
       "      <td>7s1x</td>\n",
       "      <td>A</td>\n",
       "      <td>P55011</td>\n",
       "      <td>6</td>\n",
       "      <td>1216</td>\n",
       "      <td>0</td>\n",
       "      <td>1212</td>\n",
       "      <td>2</td>\n",
       "      <td>1212</td>\n",
       "      <td>1210</td>\n",
       "      <td>SOLUTE CARRIER FAMILY 12 MEMBER 2</td>\n",
       "      <td>UniProt</td>\n",
       "      <td>P55011</td>\n",
       "      <td>2</td>\n",
       "      <td>2</td>\n",
       "      <td>True</td>\n",
       "      <td>75.000000</td>\n",
       "      <td>NP_000329.2</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>539</th>\n",
       "      <td>614613</td>\n",
       "      <td>7pmy</td>\n",
       "      <td>A</td>\n",
       "      <td>Q16348</td>\n",
       "      <td>1</td>\n",
       "      <td>729</td>\n",
       "      <td>0</td>\n",
       "      <td>None</td>\n",
       "      <td>1</td>\n",
       "      <td>729</td>\n",
       "      <td>728</td>\n",
       "      <td>SOLUTE CARRIER FAMILY 15 MEMBER 2</td>\n",
       "      <td>UniProt</td>\n",
       "      <td>Q16348</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>True</td>\n",
       "      <td>66.666667</td>\n",
       "      <td>NP_066568.3</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>175</th>\n",
       "      <td>130683</td>\n",
       "      <td>3cik</td>\n",
       "      <td>A</td>\n",
       "      <td>P25098</td>\n",
       "      <td>1</td>\n",
       "      <td>689</td>\n",
       "      <td>0</td>\n",
       "      <td>None</td>\n",
       "      <td>1</td>\n",
       "      <td>689</td>\n",
       "      <td>688</td>\n",
       "      <td>BETA-ADRENERGIC RECEPTOR KINASE 1</td>\n",
       "      <td>UniProt</td>\n",
       "      <td>P25098</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>True</td>\n",
       "      <td>66.666667</td>\n",
       "      <td>NP_005151.2</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>344</th>\n",
       "      <td>374558</td>\n",
       "      <td>5svy</td>\n",
       "      <td>A</td>\n",
       "      <td>Q14149</td>\n",
       "      <td>2</td>\n",
       "      <td>50</td>\n",
       "      <td>407</td>\n",
       "      <td>49</td>\n",
       "      <td>407</td>\n",
       "      <td>455</td>\n",
       "      <td>48</td>\n",
       "      <td>MORC FAMILY CW-TYPE ZINC FINGER PROTEIN 3</td>\n",
       "      <td>UniProt</td>\n",
       "      <td>Q14149</td>\n",
       "      <td>407</td>\n",
       "      <td>407</td>\n",
       "      <td>True</td>\n",
       "      <td>64.285714</td>\n",
       "      <td>NP_056173.1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>483</th>\n",
       "      <td>532608</td>\n",
       "      <td>6yfz</td>\n",
       "      <td>A</td>\n",
       "      <td>O14733</td>\n",
       "      <td>2</td>\n",
       "      <td>307</td>\n",
       "      <td>116</td>\n",
       "      <td>421</td>\n",
       "      <td>100</td>\n",
       "      <td>405</td>\n",
       "      <td>305</td>\n",
       "      <td>DUAL SPECIFICITY MITOGEN-ACTIVATED PROTEIN KIN...</td>\n",
       "      <td>UniProt</td>\n",
       "      <td>O14733</td>\n",
       "      <td>116</td>\n",
       "      <td>116</td>\n",
       "      <td>True</td>\n",
       "      <td>60.000000</td>\n",
       "      <td>NP_001284484.1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>221</th>\n",
       "      <td>191104</td>\n",
       "      <td>3tzd</td>\n",
       "      <td>T</td>\n",
       "      <td>P10412</td>\n",
       "      <td>2</td>\n",
       "      <td>19</td>\n",
       "      <td>19</td>\n",
       "      <td>None</td>\n",
       "      <td>19</td>\n",
       "      <td>36</td>\n",
       "      <td>17</td>\n",
       "      <td>HISTONE H1.4</td>\n",
       "      <td>UniProt</td>\n",
       "      <td>P10412</td>\n",
       "      <td>19</td>\n",
       "      <td>19</td>\n",
       "      <td>True</td>\n",
       "      <td>50.000000</td>\n",
       "      <td>NP_005316.1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>393</th>\n",
       "      <td>430676</td>\n",
       "      <td>6g2d</td>\n",
       "      <td>B</td>\n",
       "      <td>Q13085</td>\n",
       "      <td>62</td>\n",
       "      <td>2407</td>\n",
       "      <td>0</td>\n",
       "      <td>None</td>\n",
       "      <td>1</td>\n",
       "      <td>2346</td>\n",
       "      <td>2345</td>\n",
       "      <td>ACETYL-COA CARBOXYLASE 1</td>\n",
       "      <td>UniProt</td>\n",
       "      <td>Q13085</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>True</td>\n",
       "      <td>18.918919</td>\n",
       "      <td>NP_942133.1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>506</th>\n",
       "      <td>577723</td>\n",
       "      <td>7jic</td>\n",
       "      <td>A</td>\n",
       "      <td>P15391</td>\n",
       "      <td>16</td>\n",
       "      <td>323</td>\n",
       "      <td>0</td>\n",
       "      <td>None</td>\n",
       "      <td>20</td>\n",
       "      <td>327</td>\n",
       "      <td>307</td>\n",
       "      <td>B-LYMPHOCYTE ANTIGEN CD19</td>\n",
       "      <td>UniProt</td>\n",
       "      <td>P15391</td>\n",
       "      <td>20</td>\n",
       "      <td>20</td>\n",
       "      <td>True</td>\n",
       "      <td>16.666667</td>\n",
       "      <td>NP_001761.3</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>273</th>\n",
       "      <td>264693</td>\n",
       "      <td>4tvs</td>\n",
       "      <td>A</td>\n",
       "      <td>Q5JTV8</td>\n",
       "      <td>4</td>\n",
       "      <td>231</td>\n",
       "      <td>0</td>\n",
       "      <td>583</td>\n",
       "      <td>356</td>\n",
       "      <td>583</td>\n",
       "      <td>227</td>\n",
       "      <td>TORSIN-1A-INTERACTING PROTEIN 1</td>\n",
       "      <td>UniProt</td>\n",
       "      <td>Q5JTV8</td>\n",
       "      <td>356</td>\n",
       "      <td>356</td>\n",
       "      <td>True</td>\n",
       "      <td>6.976744</td>\n",
       "      <td>NP_001254507.1</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "      index   PDB CHAIN SP_PRIMARY RES_BEG RES_END PDB_BEG PDB_END SP_BEG  \\\n",
       "66    58470  1xl7     A     Q9DC50       1     612       0     612      1   \n",
       "212  174536  3ojy     B     P07358       1     537       0    None     55   \n",
       "364  397057  5xtc     f     O43677       1      47      28      74     28   \n",
       "74    63125  1z87     A     Q61234       1     263       2     264      2   \n",
       "84    71308  2c9y     A     P54819       2     240       0    None      1   \n",
       "277  288151  4wj7     A     Q9BSQ5       3     180       0    None     51   \n",
       "365  397059  5xtc     h     O43920       1     104       2     105      2   \n",
       "547  625427  7s1x     A     P55011       6    1216       0    1212      2   \n",
       "539  614613  7pmy     A     Q16348       1     729       0    None      1   \n",
       "175  130683  3cik     A     P25098       1     689       0    None      1   \n",
       "344  374558  5svy     A     Q14149       2      50     407      49    407   \n",
       "483  532608  6yfz     A     O14733       2     307     116     421    100   \n",
       "221  191104  3tzd     T     P10412       2      19      19    None     19   \n",
       "393  430676  6g2d     B     Q13085      62    2407       0    None      1   \n",
       "506  577723  7jic     A     P15391      16     323       0    None     20   \n",
       "273  264693  4tvs     A     Q5JTV8       4     231       0     583    356   \n",
       "\n",
       "    SP_END   len                                           pdb_desc   db_ref  \\\n",
       "66     612   611        PEROXISOMAL CARNITINE O-OCTANOYLTRANSFERASE  UniProt   \n",
       "212    591   536                 COMPLEMENT COMPONENT C8 BETA CHAIN  UniProt   \n",
       "364     74    46  NADH DEHYDROGENASE [UBIQUINONE] 1 SUBUNIT C1, ...  UniProt   \n",
       "74     264   262                                 ALPHA-1-SYNTROPHIN  UniProt   \n",
       "84     239   238        ADENYLATE KINASE ISOENZYME 2, MITOCHONDRIAL  UniProt   \n",
       "277    228   177                                        MALCAVERNIN  UniProt   \n",
       "365    105   103  NADH DEHYDROGENASE [UBIQUINONE] IRON-SULFUR PR...  UniProt   \n",
       "547   1212  1210                  SOLUTE CARRIER FAMILY 12 MEMBER 2  UniProt   \n",
       "539    729   728                  SOLUTE CARRIER FAMILY 15 MEMBER 2  UniProt   \n",
       "175    689   688                  BETA-ADRENERGIC RECEPTOR KINASE 1  UniProt   \n",
       "344    455    48          MORC FAMILY CW-TYPE ZINC FINGER PROTEIN 3  UniProt   \n",
       "483    405   305  DUAL SPECIFICITY MITOGEN-ACTIVATED PROTEIN KIN...  UniProt   \n",
       "221     36    17                                       HISTONE H1.4  UniProt   \n",
       "393   2346  2345                           ACETYL-COA CARBOXYLASE 1  UniProt   \n",
       "506    327   307                          B-LYMPHOCYTE ANTIGEN CD19  UniProt   \n",
       "273    583   227                    TORSIN-1A-INTERACTING PROTEIN 1  UniProt   \n",
       "\n",
       "    db_accession db_first_from db_first_to  db_match  percent_match  \\\n",
       "66        Q9DC50             1           1      True     100.000000   \n",
       "212       P07358             1          55      True     100.000000   \n",
       "364       O43677            28          28      True     100.000000   \n",
       "74        Q61234             2           2      True     100.000000   \n",
       "84        P54819             1           1      True     100.000000   \n",
       "277       Q9BSQ5            51          51      True     100.000000   \n",
       "365       O43920             2           2      True     100.000000   \n",
       "547       P55011             2           2      True      75.000000   \n",
       "539       Q16348             1           1      True      66.666667   \n",
       "175       P25098             1           1      True      66.666667   \n",
       "344       Q14149           407         407      True      64.285714   \n",
       "483       O14733           116         116      True      60.000000   \n",
       "221       P10412            19          19      True      50.000000   \n",
       "393       Q13085             1           1      True      18.918919   \n",
       "506       P15391            20          20      True      16.666667   \n",
       "273       Q5JTV8           356         356      True       6.976744   \n",
       "\n",
       "           proteins  \n",
       "66   NP_001137407.1  \n",
       "212     NP_000057.2  \n",
       "364  NP_001171915.1  \n",
       "74      NP_003089.1  \n",
       "84      NP_001616.1  \n",
       "277     NP_542192.2  \n",
       "365  NP_001171908.1  \n",
       "547     NP_000329.2  \n",
       "539     NP_066568.3  \n",
       "175     NP_005151.2  \n",
       "344     NP_056173.1  \n",
       "483  NP_001284484.1  \n",
       "221     NP_005316.1  \n",
       "393     NP_942133.1  \n",
       "506     NP_001761.3  \n",
       "273  NP_001254507.1  "
      ]
     },
     "execution_count": 26,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "pdb_mapped_filt_df"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 9. Finalize Sites\n",
    "---\n",
    "* Filter duplicates, etc."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "100%|██████████| 16/16 [00:00<00:00, 129.79it/s]\n"
     ]
    }
   ],
   "source": [
    "sites_mapped_df = list()\n",
    "\n",
    "for idx,row in tqdm(pdb_mapped_filt_df.iterrows(), total=pdb_mapped_filt_df.shape[0]):\n",
    "    _df = pd.read_parquet(\n",
    "        os.path.join(REF_DIR, \"pdb_matches/acc_mapped_sites\", \"{}.parquet\".format(row[\"proteins\"]))).loc[:,[\n",
    "        \"{}_{}_res_i\".format(row[\"PDB\"],row[\"CHAIN\"]),\n",
    "        \"{}_{}_res\".format(row[\"PDB\"],row[\"CHAIN\"]),\n",
    "        \"{}_{}_res_match\".format(row[\"PDB\"],row[\"CHAIN\"])\n",
    "    ]]\n",
    "    \n",
    "    _df.columns = ['pdb_res_i', 'pdb_res', 'pdb_res_match']\n",
    "    _df['pdb'] = row[\"PDB\"]\n",
    "    _df['chain'] = row[\"CHAIN\"]\n",
    "    sites_mapped_df.append(_df)\n",
    "    \n",
    "sites_mapped_df = pd.concat(sites_mapped_df)\n",
    "\n",
    "# Drop duplicates\n",
    "sites_mapped_df = sites_mapped_df.reset_index().drop_duplicates().set_index(\"id\")\n",
    "\n",
    "# Select for matched sites / rank multi sites by position\n",
    "sites_mapped_df = sites_mapped_df[sites_mapped_df['pdb_res_match']]\n",
    "sites_mapped_df['pdb_res_i_rank'] = sites_mapped_df.reset_index().groupby(\"id\")['pdb_res_i'].rank(\"dense\").astype(int).values"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "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>pdb_res_i</th>\n",
       "      <th>pdb_res</th>\n",
       "      <th>pdb_res_match</th>\n",
       "      <th>pdb</th>\n",
       "      <th>chain</th>\n",
       "      <th>pdb_res_i_rank</th>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>id</th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>NP_001137407.1_K386k_1_1_386_386</th>\n",
       "      <td>358</td>\n",
       "      <td>K</td>\n",
       "      <td>True</td>\n",
       "      <td>1xl7</td>\n",
       "      <td>A</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>NP_001137407.1_K523k_1_1_523_523</th>\n",
       "      <td>495</td>\n",
       "      <td>K</td>\n",
       "      <td>True</td>\n",
       "      <td>1xl7</td>\n",
       "      <td>A</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>NP_000057.2_S150s_1_1_150_150</th>\n",
       "      <td>96</td>\n",
       "      <td>S</td>\n",
       "      <td>True</td>\n",
       "      <td>3ojy</td>\n",
       "      <td>B</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>NP_000057.2_K259k_1_1_259_259</th>\n",
       "      <td>205</td>\n",
       "      <td>K</td>\n",
       "      <td>True</td>\n",
       "      <td>3ojy</td>\n",
       "      <td>B</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>NP_000057.2_K305k_1_1_305_305</th>\n",
       "      <td>251</td>\n",
       "      <td>K</td>\n",
       "      <td>True</td>\n",
       "      <td>3ojy</td>\n",
       "      <td>B</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",
       "      <td>...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>NP_942133.1_T1392t_1_1_1392_1392</th>\n",
       "      <td>1392</td>\n",
       "      <td>T</td>\n",
       "      <td>True</td>\n",
       "      <td>6g2d</td>\n",
       "      <td>B</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>NP_001761.3_T275t_1_1_275_275</th>\n",
       "      <td>275</td>\n",
       "      <td>T</td>\n",
       "      <td>True</td>\n",
       "      <td>7jic</td>\n",
       "      <td>A</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>NP_001254507.1_S467s_1_1_467_467</th>\n",
       "      <td>466</td>\n",
       "      <td>S</td>\n",
       "      <td>True</td>\n",
       "      <td>4tvs</td>\n",
       "      <td>A</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>NP_001254507.1_K471k_1_1_471_471</th>\n",
       "      <td>470</td>\n",
       "      <td>K</td>\n",
       "      <td>True</td>\n",
       "      <td>4tvs</td>\n",
       "      <td>A</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>NP_001254507.1_K505k_1_1_505_505</th>\n",
       "      <td>504</td>\n",
       "      <td>K</td>\n",
       "      <td>True</td>\n",
       "      <td>4tvs</td>\n",
       "      <td>A</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>72 rows × 6 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "                                  pdb_res_i pdb_res  pdb_res_match   pdb  \\\n",
       "id                                                                         \n",
       "NP_001137407.1_K386k_1_1_386_386        358       K           True  1xl7   \n",
       "NP_001137407.1_K523k_1_1_523_523        495       K           True  1xl7   \n",
       "NP_000057.2_S150s_1_1_150_150            96       S           True  3ojy   \n",
       "NP_000057.2_K259k_1_1_259_259           205       K           True  3ojy   \n",
       "NP_000057.2_K305k_1_1_305_305           251       K           True  3ojy   \n",
       "...                                     ...     ...            ...   ...   \n",
       "NP_942133.1_T1392t_1_1_1392_1392       1392       T           True  6g2d   \n",
       "NP_001761.3_T275t_1_1_275_275           275       T           True  7jic   \n",
       "NP_001254507.1_S467s_1_1_467_467        466       S           True  4tvs   \n",
       "NP_001254507.1_K471k_1_1_471_471        470       K           True  4tvs   \n",
       "NP_001254507.1_K505k_1_1_505_505        504       K           True  4tvs   \n",
       "\n",
       "                                 chain  pdb_res_i_rank  \n",
       "id                                                      \n",
       "NP_001137407.1_K386k_1_1_386_386     A               1  \n",
       "NP_001137407.1_K523k_1_1_523_523     A               1  \n",
       "NP_000057.2_S150s_1_1_150_150        B               1  \n",
       "NP_000057.2_K259k_1_1_259_259        B               1  \n",
       "NP_000057.2_K305k_1_1_305_305        B               1  \n",
       "...                                ...             ...  \n",
       "NP_942133.1_T1392t_1_1_1392_1392     B               1  \n",
       "NP_001761.3_T275t_1_1_275_275        A               1  \n",
       "NP_001254507.1_S467s_1_1_467_467     A               1  \n",
       "NP_001254507.1_K471k_1_1_471_471     A               1  \n",
       "NP_001254507.1_K505k_1_1_505_505     A               1  \n",
       "\n",
       "[72 rows x 6 columns]"
      ]
     },
     "execution_count": 29,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "sites_mapped_df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Drop duplicate ptmSites (i.e. NP_000005.2_K608k_1_1_608_608 - has two sites mapped from different cohorts)\n",
    "ptm_comb_filt_pdb_df = ptm_comb_filt_df.reset_index().drop_duplicates(['id','ptmSite']).set_index('id')\n",
    "\n",
    "# Rank multi-sites by position\n",
    "ptm_comb_filt_pdb_df['pdb_res_i_rank'] = ptm_comb_filt_pdb_df.reset_index().groupby(\"id\")['acc_res_i'].rank(\"dense\").astype(int).values\n",
    "\n",
    "# Combine pdb level annotations\n",
    "ptm_comb_filt_pdb_df = pd.merge(ptm_comb_filt_pdb_df.reset_index(), sites_mapped_df.reset_index(), how='left').set_index('id')\n",
    "ptm_comb_filt_pdb_df['pdb_res_match'] = ptm_comb_filt_pdb_df['pdb_res_match'].fillna(False)\n",
    "ptm_comb_filt_pdb_df['pdb_res_match'] = (ptm_comb_filt_pdb_df['pdb_res_match']) & (ptm_comb_filt_pdb_df['pdb_res'] == ptm_comb_filt_pdb_df['acc_res'])\n",
    "\n",
    "# Only matched sites\n",
    "ptm_comb_filt_match_pdb_df = ptm_comb_filt_pdb_df[ptm_comb_filt_pdb_df['pdb_res_match']==True].copy()\n",
    "ptm_comb_filt_match_pdb_df['pdb_res_i'] = ptm_comb_filt_match_pdb_df['pdb_res_i'].astype(int)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 10. Save Files + Summary\n",
    "---"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {},
   "outputs": [],
   "source": [
    "ptm_comb_filt_pdb_df.to_csv(os.path.join(REF_DIR, \"full_mapped_sites_to_pdbs.tsv\"), sep='\\t')\n",
    "ptm_comb_filt_match_pdb_df.to_csv(os.path.join(REF_DIR, \"mapped_sites_to_pdbs.tsv\"), sep='\\t')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      " 23 acetylsites mapped.\n",
      " 48 phosphosites mapped.\n",
      " 16 unique accession numbers.\n",
      " 16 unique uniprots.\n"
     ]
    }
   ],
   "source": [
    "print(\" {} acetylsites mapped.\".format(ptm_comb_filt_match_pdb_df.groupby('feature').size()['acetylome']))\n",
    "print(\" {} phosphosites mapped.\".format(ptm_comb_filt_match_pdb_df.groupby('feature').size()['phosphoproteome']))\n",
    "print(\" {} unique accession numbers.\".format(np.unique(ptm_comb_filt_match_pdb_df['accession_number']).shape[0]))\n",
    "print(\" {} unique uniprots.\".format(np.unique(ptm_comb_filt_match_pdb_df['uniprot']).shape[0]))"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}