notebooks/dataprep/01b-NCBIGeneProtein.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Load Viral Genome Data from NCBI\n",
    "**[Work in progress]**\n",
    "\n",
    "This notebook downloads and standardizes viral genome data from NCBI for ingestion into a Knowledge Graph.\n",
    "\n",
    "Data source: [NCBI](https://www.ncbi.nlm.nih.gov)\n",
    "\n",
    "Author: Peter Rose (pwrose@ucsd.edu)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import os\n",
    "import re\n",
    "import pandas as pd\n",
    "import dateutil\n",
    "import time\n",
    "from pathlib import Path\n",
    "import hashlib \n",
    "\n",
    "from Bio import SeqIO\n",
    "from Bio import Entrez\n",
    "Entrez.email = \"covid19@mybinder.org\""
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "pd.options.display.max_rows = None  # display all rows\n",
    "pd.options.display.max_columns = None  # display all columsns"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "/Users/peter/Library/Application Support/com.Neo4j.Relate/data/dbmss/dbms-8bf637fc-0d20-4d9f-9c6f-f7e72e92a4da/import\n"
     ]
    }
   ],
   "source": [
    "# Path will take care of handling operating system differences.\n",
    "NEO4J_IMPORT = Path(os.getenv('NEO4J_IMPORT'))\n",
    "print(NEO4J_IMPORT)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Download viral genome information"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "# TODO loop over a list of taxonomy ids. \n",
    "\n",
    "# for now SARS-CoV-2 is hardcoded.\n",
    "#tax_id = \"2697049\" # SARS-CoV-2, NC_045512\n",
    "#tax_id = \"227984\" # SARS, NC_004718\n",
    "#tax_id = \"1335626\" # MERS-CoV, NC_019843, NC_038294(UniProt)\n",
    "\n",
    "#tax_id = \"147711\" # Rhinovirus A NC_001617\n",
    "#tax_id = \"147712\" # Rhinovirus B (B14: NC_001490)\n",
    "#https://www.ncbi.nlm.nih.gov/genomes/GenomesGroup.cgi?taxid=12058&host=human\n",
    "\n",
    "#tax_id = \"186538\" # Zaire ebolavirus, NC_002549\n",
    "\n",
    "# https://www.ncbi.nlm.nih.gov/ipg/YP_009724389.1 identical protein groups\n",
    "#https://www.ebi.ac.uk/uniprot/api/covid-19/uniprotkb/accession/P0DTC1.txt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "organisms = pd.read_csv(\"../../reference_data/Organism.csv\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "pathogens = organisms.query(\"type == 'Pathogen'\")\n",
    "#refSeq = pathogens['genomeAccession'].values"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "genomes = pd.read_csv(\"../../reference_data/Genome.csv\", dtype=str)\n",
    "# restrict list to viral genomes\n",
    "genomes = genomes.query(\"chromosome == 'Viral Chromosome'\")\n",
    "genomes.dropna(inplace=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "refSeq = genomes['refSeq'].str.split(':', expand=True)[1]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0     NC_045512\n",
       "1     NC_038294\n",
       "2     NC_004718\n",
       "3     NC_006577\n",
       "4     NC_006213\n",
       "5     NC_002645\n",
       "6     NC_005831\n",
       "7     NC_001490\n",
       "9     NC_023760\n",
       "10    NC_009988\n",
       "Name: 1, dtype: object"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "refSeq"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Download reference sequence nucleotide accession numbers for the given taxonomy id"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "# TODO refactor function into a separate .py file\n",
    "def extract_genbank_features(gb_record, feature_type, qualifier_list) :\n",
    "    answers = list()\n",
    "    for (index, feature) in enumerate(gb_record.features) :\n",
    "        answer = dict()\n",
    "        if feature.type==feature_type :\n",
    "            for qualifier in qualifier_list:\n",
    "                if qualifier in feature.qualifiers :\n",
    "                    #print(feature, feature.qualifiers)\n",
    "                    # keep NCBI \"1\"-based indexing\n",
    "                    answer['start'] = min(feature.location) + 1\n",
    "                    answer['end'] = max(feature.location) + 1\n",
    "                    #There should only be one locus_tag per feature, but there\n",
    "                    #are usually several db_xref entries\n",
    "                    for value in feature.qualifiers[qualifier] :\n",
    "                        if value in answer :\n",
    "                                print(\"WARNING - Duplicate key %s for %s features %i and %i\" \\\n",
    "                                % (value, feature_type, answer[value], index))\n",
    "                        else :\n",
    "                            answer[qualifier] = value\n",
    "                else:\n",
    "                    answer[qualifier] = ''\n",
    "            answers.append(answer)\n",
    "    return answers"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [],
   "source": [
    "def process_protein_record(accession):\n",
    "    with Entrez.efetch(db=\"nucleotide\", rettype=\"gb\", retmode=\"text\", id=accession) as handle:\n",
    "        record = SeqIO.read(handle, \"gb\")\n",
    "        print(record)\n",
    "    \n",
    "    source = extract_genbank_features(record, \"source\", [\"db_xref\"])\n",
    "    taxonomy = source[0]['db_xref']\n",
    "    taxonomy = taxonomy.split(':')[1]\n",
    "    protein = extract_genbank_features(record, \"CDS\", [\"gene\",\"locus_tag\",\"db_xref\",\"product\",\"protein_id\"])\n",
    "    protein_df = pd.DataFrame(protein)\n",
    "    protein_df['fullLength'] = 'True'\n",
    "\n",
    "    mature_peptide = extract_genbank_features(record, \"mat_peptide\", [\"gene\",\"locus_tag\",\"db_xref\",\"product\",\"protein_id\"])\n",
    "    \n",
    "    # some reference sequences are missing the mature peptides (why???)\n",
    "    if len(mature_peptide) > 0:\n",
    "        mature_peptide_df = pd.DataFrame(mature_peptide)\n",
    "        mature_peptide_df['fullLength'] = 'False'\n",
    "\n",
    "        # assign db_xref to the mature proteins based on the precursor's gene id\n",
    "        mature_peptide_df['db_xref'] = mature_peptide_df['gene'].apply(lambda g: protein_df.query(f\"gene=='{g}'\")['db_xref'][0])\n",
    "        df = pd.concat([protein_df, mature_peptide_df])\n",
    "    else:\n",
    "        df = protein_df\n",
    "        \n",
    "    # TODO fix for taxonomy missmatch between NCBI and UniProt for SARS\n",
    "    if taxonomy == '227984':\n",
    "        taxonomy = '694009'\n",
    "        \n",
    "    df['taxonomyId'] = taxonomy\n",
    "    df['genomeName'] = record.description\n",
    "    df['genbank_accession'] = accession\n",
    "    \n",
    "    return df"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Concatenate all Genbank dataframes"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "ID: NC_045512.2\n",
      "Name: NC_045512\n",
      "Description: Severe acute respiratory syndrome coronavirus 2 isolate Wuhan-Hu-1, complete genome\n",
      "Database cross-references: BioProject:PRJNA485481\n",
      "Number of features: 57\n",
      "/molecule_type=ss-RNA\n",
      "/topology=linear\n",
      "/data_file_division=VRL\n",
      "/date=18-JUL-2020\n",
      "/accessions=['NC_045512']\n",
      "/sequence_version=2\n",
      "/keywords=['RefSeq']\n",
      "/source=Severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2)\n",
      "/organism=Severe acute respiratory syndrome coronavirus 2\n",
      "/taxonomy=['Viruses', 'Riboviria', 'Orthornavirae', 'Pisuviricota', 'Pisoniviricetes', 'Nidovirales', 'Cornidovirineae', 'Coronaviridae', 'Orthocoronavirinae', 'Betacoronavirus', 'Sarbecovirus']\n",
      "/references=[Reference(title='A new coronavirus associated with human respiratory disease in China', ...), Reference(title='Programmed ribosomal frameshifting in decoding the SARS-CoV genome', ...), Reference(title='The structure of a rigorously conserved RNA element within the SARS virus genome', ...), Reference(title=\"A phylogenetically conserved hairpin-type 3' untranslated region pseudoknot functions in coronavirus RNA replication\", ...), Reference(title='Direct Submission', ...), Reference(title='Direct Submission', ...)]\n",
      "/comment=REVIEWED REFSEQ: This record has been curated by NCBI staff. The\n",
      "reference sequence is identical to MN908947.\n",
      "On Jan 17, 2020 this sequence version replaced NC_045512.1.\n",
      "Annotation was added using homology to SARSr-CoV NC_004718.3. ###\n",
      "Formerly called 'Wuhan seafood market pneumonia virus.' If you have\n",
      "questions or suggestions, please email us at info@ncbi.nlm.nih.gov\n",
      "and include the accession number NC_045512.### Protein structures\n",
      "can be found at\n",
      "https://www.ncbi.nlm.nih.gov/structure/?term=sars-cov-2.### Find\n",
      "all other Severe acute respiratory syndrome coronavirus 2\n",
      "(SARS-CoV-2) sequences at\n",
      "https://www.ncbi.nlm.nih.gov/genbank/sars-cov-2-seqs/\n",
      "COMPLETENESS: full length.\n",
      "/structured_comment=OrderedDict([('Assembly-Data', OrderedDict([('Assembly Method', 'Megahit v. V1.1.3'), ('Sequencing Technology', 'Illumina')]))])\n",
      "Seq('ATTAAAGGTTTATACCTTCCCAGGTAACAAACCAACCAACTTTCGATCTCTTGT...AAA')\n",
      "ID: NC_038294.1\n",
      "Name: NC_038294\n",
      "Description: Betacoronavirus England 1 isolate H123990006, complete genome\n",
      "Database cross-references: BioProject:PRJNA485481\n",
      "Number of features: 50\n",
      "/molecule_type=RNA\n",
      "/topology=linear\n",
      "/data_file_division=VRL\n",
      "/date=22-NOV-2020\n",
      "/accessions=['NC_038294']\n",
      "/sequence_version=1\n",
      "/keywords=['RefSeq']\n",
      "/source=Betacoronavirus England 1\n",
      "/organism=Betacoronavirus England 1\n",
      "/taxonomy=['Viruses', 'Riboviria', 'Orthornavirae', 'Pisuviricota', 'Pisoniviricetes', 'Nidovirales', 'Cornidovirineae', 'Coronaviridae', 'Orthocoronavirinae', 'Betacoronavirus', 'Merbecovirus']\n",
      "/references=[Reference(title='Middle East respiratory syndrome coronavirus (MERS-CoV): announcement of the Coronavirus Study Group', ...), Reference(title='Full-genome deep sequencing and phylogenetic analysis of novel human betacoronavirus', ...), Reference(title='Direct Submission', ...), Reference(title='Direct Submission', ...), Reference(title='Direct Submission', ...)]\n",
      "/comment=REVIEWED REFSEQ: This record has been curated by NCBI staff. The\n",
      "reference sequence is identical to KC164505.\n",
      "Annotation was based on information found in PMID: 23170002 and the\n",
      "annotation of Middle East respiratory syndrome coronavirus strain\n",
      "HCoV-EMC, NC_019843.3. The sequence of NC_019843.3 is 99.67%\n",
      "similar to this sequence.\n",
      "COMPLETENESS: full length.\n",
      "Seq('ATTTAAGTGAATAGCTTGGCTATCTCACTTCCCCTCGTTCTCTTGCAGAACTTT...AAA')\n",
      "ID: NC_004718.3\n",
      "Name: NC_004718\n",
      "Description: SARS coronavirus Tor2, complete genome\n",
      "Database cross-references: BioProject:PRJNA485481\n",
      "Number of features: 57\n",
      "/molecule_type=RNA\n",
      "/topology=linear\n",
      "/data_file_division=VRL\n",
      "/date=20-NOV-2020\n",
      "/accessions=['NC_004718']\n",
      "/sequence_version=3\n",
      "/keywords=['RefSeq']\n",
      "/source=SARS coronavirus Tor2\n",
      "/organism=SARS coronavirus Tor2\n",
      "/taxonomy=['Viruses', 'Riboviria', 'Orthornavirae', 'Pisuviricota', 'Pisoniviricetes', 'Nidovirales', 'Cornidovirineae', 'Coronaviridae', 'Orthocoronavirinae', 'Betacoronavirus', 'Sarbecovirus']\n",
      "/references=[Reference(title='Analysis of multimerization of the SARS coronavirus nucleocapsid protein', ...), Reference(title='Direct Submission', ...), Reference(title='Direct Submission', ...), Reference(title='Direct Submission', ...), Reference(title='Direct Submission', ...)]\n",
      "/comment=REVIEWED REFSEQ: This record has been curated by NCBI staff. The\n",
      "reference sequence is identical to AY274119.\n",
      "On or before Mar 28, 2016 this sequence version replaced\n",
      "NC_028858.1, NC_028866.1, NC_028884.1, NC_028893.1, NC_028873.1,\n",
      "NC_028845.1, NC_009696.1, NC_009695.1, NC_013664.1, NC_009693.1,\n",
      "NC_009694.1, NC_004718.2.\n",
      "Annotation based on that found in PMID: 31987001, PMID: 31967327\n",
      "and the annotation of P0C6X7.1.\n",
      "COMPLETENESS: full length.\n",
      "Seq('ATATTAGGTTTTTACCTACCCAGGAAAAGCCAACCAACCTCGATCTCTTGTAGA...AAA')\n",
      "ID: NC_006577.2\n",
      "Name: NC_006577\n",
      "Description: Human coronavirus HKU1, complete genome\n",
      "Database cross-references: BioProject:PRJNA485481\n",
      "Number of features: 48\n",
      "/molecule_type=ss-RNA\n",
      "/topology=linear\n",
      "/data_file_division=VRL\n",
      "/date=20-NOV-2020\n",
      "/accessions=['NC_006577']\n",
      "/sequence_version=2\n",
      "/keywords=['RefSeq']\n",
      "/source=Human coronavirus HKU1 (HCoV-HKU1)\n",
      "/organism=Human coronavirus HKU1\n",
      "/taxonomy=['Viruses', 'Riboviria', 'Orthornavirae', 'Pisuviricota', 'Pisoniviricetes', 'Nidovirales', 'Cornidovirineae', 'Coronaviridae', 'Orthocoronavirinae', 'Betacoronavirus', 'Embecovirus']\n",
      "/references=[Reference(title='Characterization and complete genome sequence of a novel coronavirus, coronavirus HKU1, from patients with pneumonia', ...), Reference(title='Direct Submission', ...), Reference(title='Direct Submission', ...), Reference(title='Direct Submission', ...)]\n",
      "/comment=REVIEWED REFSEQ: This record has been curated by NCBI staff. The\n",
      "reference sequence was derived from AY597011.\n",
      "On Jan 23, 2006 this sequence version replaced NC_006577.1.\n",
      "The mature peptides were designated nsp1 through nsp16 following\n",
      "the nomenclature used in the MHV-A59 and SARS-CoV RefSeqs\n",
      "(NC_001846 and NC_004718, respectively).\n",
      "COMPLETENESS: full length.\n",
      "Seq('GAGTTTGAGCGATTGACGTTCGTACCGTCTATCAGCTTACGATCTCTTGTCAGA...TAC')\n",
      "ID: NC_006213.1\n",
      "Name: NC_006213\n",
      "Description: Human coronavirus OC43 strain ATCC VR-759, complete genome\n",
      "Database cross-references: BioProject:PRJNA485481\n",
      "Number of features: 47\n",
      "/molecule_type=RNA\n",
      "/topology=linear\n",
      "/data_file_division=VRL\n",
      "/date=16-SEP-2020\n",
      "/accessions=['NC_006213']\n",
      "/sequence_version=1\n",
      "/keywords=['RefSeq']\n",
      "/source=Human coronavirus OC43 (HCoV-OC43)\n",
      "/organism=Human coronavirus OC43\n",
      "/taxonomy=['Viruses', 'Riboviria', 'Orthornavirae', 'Pisuviricota', 'Pisoniviricetes', 'Nidovirales', 'Cornidovirineae', 'Coronaviridae', 'Orthocoronavirinae', 'Betacoronavirus', 'Embecovirus']\n",
      "/references=[Reference(title='Human Respiratory Coronavirus OC43: Genetic Stability and Neuroinvasion', ...), Reference(title='Direct Submission', ...), Reference(title='Direct Submission', ...)]\n",
      "/comment=REVIEWED REFSEQ: This record has been curated by NCBI staff. The\n",
      "reference sequence is identical to AY585228.\n",
      "Annotation was based on information found in PMID: 16537637, PMID:\n",
      "27381385 and the annotation of KF530089.1.\n",
      "COMPLETENESS: full length.\n",
      "Seq('ATTGTGAGCGATTTGCGTGCGTGCATCCCGCTTCACTGATCTCTTGTTAGATCT...AAA')\n",
      "ID: NC_002645.1\n",
      "Name: NC_002645\n",
      "Description: Human coronavirus 229E, complete genome\n",
      "Database cross-references: BioProject:PRJNA485481\n",
      "Number of features: 16\n",
      "/molecule_type=RNA\n",
      "/topology=linear\n",
      "/data_file_division=VRL\n",
      "/date=13-AUG-2018\n",
      "/accessions=['NC_002645']\n",
      "/sequence_version=1\n",
      "/keywords=['RefSeq']\n",
      "/source=Human coronavirus 229E\n",
      "/organism=Human coronavirus 229E\n",
      "/taxonomy=['Viruses', 'Riboviria', 'Orthornavirae', 'Pisuviricota', 'Pisoniviricetes', 'Nidovirales', 'Cornidovirineae', 'Coronaviridae', 'Orthocoronavirinae', 'Alphacoronavirus', 'Duvinacovirus']\n",
      "/references=[Reference(title='Infectious RNA transcribed in vitro from a cDNA copy of the human coronavirus genome cloned in vaccinia virus', ...), Reference(title='Direct Submission', ...), Reference(title='Direct Submission', ...)]\n",
      "/comment=PROVISIONAL REFSEQ: This record has not yet been subject to final\n",
      "NCBI review. The reference sequence is identical to AF304460.\n",
      "COMPLETENESS: full length.\n",
      "Seq('ACTTAAGTACCTTATCTATCTACAGATAGAAAAGTTGCTTTTTAGACTTTGTGT...AAA')\n",
      "ID: NC_005831.2\n",
      "Name: NC_005831\n",
      "Description: Human Coronavirus NL63, complete genome\n",
      "Database cross-references: BioProject:PRJNA485481\n",
      "Number of features: 15\n",
      "/molecule_type=ss-RNA\n",
      "/topology=linear\n",
      "/data_file_division=VRL\n",
      "/date=13-AUG-2018\n",
      "/accessions=['NC_005831']\n",
      "/sequence_version=2\n",
      "/keywords=['RefSeq']\n",
      "/source=Human coronavirus NL63 (HCoV-NL63)\n",
      "/organism=Human coronavirus NL63\n",
      "/taxonomy=['Viruses', 'Riboviria', 'Orthornavirae', 'Pisuviricota', 'Pisoniviricetes', 'Nidovirales', 'Cornidovirineae', 'Coronaviridae', 'Orthocoronavirinae', 'Alphacoronavirus', 'Setracovirus']\n",
      "/references=[Reference(title='Identification of a new human coronavirus', ...), Reference(title='Direct Submission', ...), Reference(title='Direct Submission', ...), Reference(title='Direct Submission', ...)]\n",
      "/comment=PROVISIONAL REFSEQ: This record has not yet been subject to final\n",
      "NCBI review. The reference sequence is identical to AY567487.\n",
      "On Jun 24, 2004 this sequence version replaced NC_005831.1.\n",
      "COMPLETENESS: full length.\n",
      "Seq('CTTAAAGAATTTTTCTATCTATAGATAGAGAATTTTCTTATTTAGACTTTGTGT...CAC')\n",
      "ID: NC_001490.1\n",
      "Name: NC_001490\n",
      "Description: Human rhinovirus 14, complete genome\n",
      "Database cross-references: BioProject:PRJNA485481\n",
      "Number of features: 14\n",
      "/molecule_type=ss-RNA\n",
      "/topology=linear\n",
      "/data_file_division=VRL\n",
      "/date=13-AUG-2018\n",
      "/accessions=['NC_001490']\n",
      "/sequence_version=1\n",
      "/keywords=['RefSeq']\n",
      "/source=Rhinovirus B14\n",
      "/organism=Rhinovirus B14\n",
      "/taxonomy=['Viruses', 'Riboviria', 'Orthornavirae', 'Pisuviricota', 'Pisoniviricetes', 'Picornavirales', 'Picornaviridae', 'Enterovirus']\n",
      "/references=[Reference(title='Molecular cloning and complete sequence determination of RNA genome of human rhinovirus type 14', ...), Reference(title='The complete nucleotide sequence of a common cold virus: human rhinovirus 14', ...), Reference(title='Direct Submission', ...)]\n",
      "/comment=PROVISIONAL REFSEQ: This record has not yet been subject to final\n",
      "NCBI review. The reference sequence is identical to K02121.\n",
      "The HRV-14 virion contains a ss-RNA genome that serves as a\n",
      "monocistronic mRNA for the synthesis of all viral structural and\n",
      "non-structural proteins.  The genome is polyadenylated on the 3'\n",
      "end and may be capped with a protein (VPg) on the 5' end (by\n",
      "analogy with other picornaviruses).\n",
      "Based on homology, HRV14 appears most closely related to poliovirus\n",
      "type 1 (PV-1); on these grounds, both [1] and [2] tentatively\n",
      "predict a number of mature products resulting from cleavage of the\n",
      "polyprotein precursor. The only major differences in these\n",
      "identifications concern the VP3 and VP1 coding sequences, about\n",
      "which [1] is uncertain.\n",
      "The sequence differences annotated as conflicts herein may in fact\n",
      "be natural variations.\n",
      "COMPLETENESS: full length.\n",
      "Seq('TTAAAACAGCGGATGGGTATCCCACCATTCGACCCATTGGGTGTAGTACTCTGG...TAT')\n",
      "ID: NC_023760.1\n",
      "Name: NC_023760\n",
      "Description: Mink coronavirus strain WD1127, complete genome\n",
      "Database cross-references: BioProject:PRJNA485481\n",
      "Number of features: 21\n",
      "/molecule_type=ss-RNA\n",
      "/topology=linear\n",
      "/data_file_division=VRL\n",
      "/date=13-AUG-2018\n",
      "/accessions=['NC_023760']\n",
      "/sequence_version=1\n",
      "/keywords=['RefSeq']\n",
      "/source=Mink coronavirus strain WD1127\n",
      "/organism=Mink coronavirus strain WD1127\n",
      "/taxonomy=['Viruses', 'Riboviria', 'Orthornavirae', 'Pisuviricota', 'Pisoniviricetes', 'Nidovirales', 'Cornidovirineae', 'Coronaviridae', 'Orthocoronavirinae', 'Alphacoronavirus', 'Minacovirus']\n",
      "/references=[Reference(title='Molecular characterization of a new species in the genus Alphacoronavirus associated with mink epizootic catarrhal gastroenteritis', ...), Reference(title='Direct Submission', ...), Reference(title='Direct Submission', ...)]\n",
      "/comment=PROVISIONAL REFSEQ: This record has not yet been subject to final\n",
      "NCBI review. The reference sequence is identical to HM245925.\n",
      "COMPLETENESS: full length.\n",
      "Seq('ACTTTTAAAGATAAGTGCGTGTAGCGTAGCTGCCCACCTTTCTTTAACTTGACT...AAA')\n",
      "ID: NC_009988.1\n",
      "Name: NC_009988\n",
      "Description: Bat coronavirus HKU2, complete genome\n",
      "Database cross-references: BioProject:PRJNA485481\n",
      "Number of features: 17\n",
      "/molecule_type=ss-RNA\n",
      "/topology=linear\n",
      "/data_file_division=VRL\n",
      "/date=13-AUG-2018\n",
      "/accessions=['NC_009988']\n",
      "/sequence_version=1\n",
      "/keywords=['RefSeq']\n",
      "/source=Rhinolophus bat coronavirus HKU2\n",
      "/organism=Rhinolophus bat coronavirus HKU2\n",
      "/taxonomy=['Viruses', 'Riboviria', 'Orthornavirae', 'Pisuviricota', 'Pisoniviricetes', 'Nidovirales', 'Cornidovirineae', 'Coronaviridae', 'Orthocoronavirinae', 'Alphacoronavirus', 'Rhinacovirus']\n",
      "/references=[Reference(title='Complete genome sequence of bat coronavirus HKU2 from Chinese horseshoe bats revealed a much smaller spike gene with a different evolutionary lineage from the rest of the genome', ...), Reference(title='Direct Submission', ...), Reference(title='Direct Submission', ...)]\n",
      "/comment=PROVISIONAL REFSEQ: This record has not yet been subject to final\n",
      "NCBI review. The reference sequence was derived from EF203064.\n",
      "COMPLETENESS: full length.\n",
      "Seq('GACTTAAAGATATAATCTATCTGTCGATAGAGTCCTTATCTTTTTAGACTTTCC...AAA')\n"
     ]
    }
   ],
   "source": [
    "df = pd.concat((process_protein_record(accession)) for accession in refSeq)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Get protein sequence and calculate md5 hash code as a unique id"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_sequence(accession):\n",
    "    with Entrez.efetch(db=\"protein\", rettype=\"gb\", retmode=\"text\", id=accession) as handle:\n",
    "        record = SeqIO.read(handle, \"gb\")\n",
    "        time.sleep(0.5)\n",
    "        return str(record.seq)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "df['sequence'] = df['protein_id'].apply(get_sequence)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [],
   "source": [
    "df['id'] = df['sequence'].apply(lambda seq: 'md5:' + hashlib.md5(seq.encode()).hexdigest())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [],
   "source": [
    "# disambiguate id by taxonomyId (same sequence for different organisms)\n",
    "df['id'] = df['id'] + '-' + df['taxonomyId']"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Create unique and interoperable identifiers"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**proteinAccession**: CURIE: [ncbiprotein](https://registry.identifiers.org/registry/ncbiprotein) (NCBI Reference Sequences, Refseq)\n",
    "\n",
    "**geneAccession**: CURIE: [ncbigene](https://registry.identifiers.org/registry/ncbigene) (NCBI Gene)\n",
    "\n",
    "A [CURIE](https://en.wikipedia.org/wiki/CURIE) (Compact URI) is a compact abbreviation for Uniform Resource Identifiers (URIs) that can be resolved by [Identifiers.org](https://identifiers.org/)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [],
   "source": [
    "df['taxonomyId'] = 'taxonomy:' + df['taxonomyId']\n",
    "df['proteinAccession'] = 'ncbiprotein:' + df['protein_id'].apply(lambda s: s.split('.')[0])\n",
    "df['geneAccession'] = 'ncbigene:' + df['db_xref'].str.split(':', expand=True)[1]\n",
    "df.rename(columns={'gene': 'geneName', 'product': 'proteinName', 'start': 'geneStart', 'end': 'geneEnd'}, inplace=True)\n",
    "df.fillna('', inplace=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "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>geneStart</th>\n",
       "      <th>geneEnd</th>\n",
       "      <th>geneName</th>\n",
       "      <th>locus_tag</th>\n",
       "      <th>db_xref</th>\n",
       "      <th>proteinName</th>\n",
       "      <th>protein_id</th>\n",
       "      <th>fullLength</th>\n",
       "      <th>taxonomyId</th>\n",
       "      <th>genomeName</th>\n",
       "      <th>genbank_accession</th>\n",
       "      <th>sequence</th>\n",
       "      <th>id</th>\n",
       "      <th>proteinAccession</th>\n",
       "      <th>geneAccession</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>266</td>\n",
       "      <td>21555</td>\n",
       "      <td>ORF1ab</td>\n",
       "      <td>GU280_gp01</td>\n",
       "      <td>GeneID:43740578</td>\n",
       "      <td>ORF1ab polyprotein</td>\n",
       "      <td>YP_009724389.1</td>\n",
       "      <td>True</td>\n",
       "      <td>taxonomy:2697049</td>\n",
       "      <td>Severe acute respiratory syndrome coronavirus ...</td>\n",
       "      <td>NC_045512</td>\n",
       "      <td>MESLVPGFNEKTHVQLSLPVLQVRDVLVRGFGDSVEEVLSEARQHL...</td>\n",
       "      <td>md5:e6608b50fcd6e004708a875615ddf2d9-2697049</td>\n",
       "      <td>ncbiprotein:YP_009724389</td>\n",
       "      <td>ncbigene:43740578</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>266</td>\n",
       "      <td>13483</td>\n",
       "      <td>ORF1ab</td>\n",
       "      <td>GU280_gp01</td>\n",
       "      <td>GeneID:43740578</td>\n",
       "      <td>ORF1a polyprotein</td>\n",
       "      <td>YP_009725295.1</td>\n",
       "      <td>True</td>\n",
       "      <td>taxonomy:2697049</td>\n",
       "      <td>Severe acute respiratory syndrome coronavirus ...</td>\n",
       "      <td>NC_045512</td>\n",
       "      <td>MESLVPGFNEKTHVQLSLPVLQVRDVLVRGFGDSVEEVLSEARQHL...</td>\n",
       "      <td>md5:e781b58591b8dbdd15f84dcbdec82105-2697049</td>\n",
       "      <td>ncbiprotein:YP_009725295</td>\n",
       "      <td>ncbigene:43740578</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>21563</td>\n",
       "      <td>25384</td>\n",
       "      <td>S</td>\n",
       "      <td>GU280_gp02</td>\n",
       "      <td>GeneID:43740568</td>\n",
       "      <td>surface glycoprotein</td>\n",
       "      <td>YP_009724390.1</td>\n",
       "      <td>True</td>\n",
       "      <td>taxonomy:2697049</td>\n",
       "      <td>Severe acute respiratory syndrome coronavirus ...</td>\n",
       "      <td>NC_045512</td>\n",
       "      <td>MFVFLVLLPLVSSQCVNLTTRTQLPPAYTNSFTRGVYYPDKVFRSS...</td>\n",
       "      <td>md5:4c35f09aac2f7be4f3cffd30c6aecac8-2697049</td>\n",
       "      <td>ncbiprotein:YP_009724390</td>\n",
       "      <td>ncbigene:43740568</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>25393</td>\n",
       "      <td>26220</td>\n",
       "      <td>ORF3a</td>\n",
       "      <td>GU280_gp03</td>\n",
       "      <td>GeneID:43740569</td>\n",
       "      <td>ORF3a protein</td>\n",
       "      <td>YP_009724391.1</td>\n",
       "      <td>True</td>\n",
       "      <td>taxonomy:2697049</td>\n",
       "      <td>Severe acute respiratory syndrome coronavirus ...</td>\n",
       "      <td>NC_045512</td>\n",
       "      <td>MDLFMRIFTIGTVTLKQGEIKDATPSDFVRATATIPIQASLPFGWL...</td>\n",
       "      <td>md5:f5c8b89ceac3f14e456577557df1ef40-2697049</td>\n",
       "      <td>ncbiprotein:YP_009724391</td>\n",
       "      <td>ncbigene:43740569</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>26245</td>\n",
       "      <td>26472</td>\n",
       "      <td>E</td>\n",
       "      <td>GU280_gp04</td>\n",
       "      <td>GeneID:43740570</td>\n",
       "      <td>envelope protein</td>\n",
       "      <td>YP_009724392.1</td>\n",
       "      <td>True</td>\n",
       "      <td>taxonomy:2697049</td>\n",
       "      <td>Severe acute respiratory syndrome coronavirus ...</td>\n",
       "      <td>NC_045512</td>\n",
       "      <td>MYSFVSEETGTLIVNSVLLFLAFVVFLLVTLAILTALRLCAYCCNI...</td>\n",
       "      <td>md5:375e0f905c315e06a99c80b736c125d2-2697049</td>\n",
       "      <td>ncbiprotein:YP_009724392</td>\n",
       "      <td>ncbigene:43740570</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   geneStart  geneEnd geneName   locus_tag          db_xref  \\\n",
       "0        266    21555   ORF1ab  GU280_gp01  GeneID:43740578   \n",
       "1        266    13483   ORF1ab  GU280_gp01  GeneID:43740578   \n",
       "2      21563    25384        S  GU280_gp02  GeneID:43740568   \n",
       "3      25393    26220    ORF3a  GU280_gp03  GeneID:43740569   \n",
       "4      26245    26472        E  GU280_gp04  GeneID:43740570   \n",
       "\n",
       "            proteinName      protein_id fullLength        taxonomyId  \\\n",
       "0    ORF1ab polyprotein  YP_009724389.1       True  taxonomy:2697049   \n",
       "1     ORF1a polyprotein  YP_009725295.1       True  taxonomy:2697049   \n",
       "2  surface glycoprotein  YP_009724390.1       True  taxonomy:2697049   \n",
       "3         ORF3a protein  YP_009724391.1       True  taxonomy:2697049   \n",
       "4      envelope protein  YP_009724392.1       True  taxonomy:2697049   \n",
       "\n",
       "                                          genomeName genbank_accession  \\\n",
       "0  Severe acute respiratory syndrome coronavirus ...         NC_045512   \n",
       "1  Severe acute respiratory syndrome coronavirus ...         NC_045512   \n",
       "2  Severe acute respiratory syndrome coronavirus ...         NC_045512   \n",
       "3  Severe acute respiratory syndrome coronavirus ...         NC_045512   \n",
       "4  Severe acute respiratory syndrome coronavirus ...         NC_045512   \n",
       "\n",
       "                                            sequence  \\\n",
       "0  MESLVPGFNEKTHVQLSLPVLQVRDVLVRGFGDSVEEVLSEARQHL...   \n",
       "1  MESLVPGFNEKTHVQLSLPVLQVRDVLVRGFGDSVEEVLSEARQHL...   \n",
       "2  MFVFLVLLPLVSSQCVNLTTRTQLPPAYTNSFTRGVYYPDKVFRSS...   \n",
       "3  MDLFMRIFTIGTVTLKQGEIKDATPSDFVRATATIPIQASLPFGWL...   \n",
       "4  MYSFVSEETGTLIVNSVLLFLAFVVFLLVTLAILTALRLCAYCCNI...   \n",
       "\n",
       "                                             id          proteinAccession  \\\n",
       "0  md5:e6608b50fcd6e004708a875615ddf2d9-2697049  ncbiprotein:YP_009724389   \n",
       "1  md5:e781b58591b8dbdd15f84dcbdec82105-2697049  ncbiprotein:YP_009725295   \n",
       "2  md5:4c35f09aac2f7be4f3cffd30c6aecac8-2697049  ncbiprotein:YP_009724390   \n",
       "3  md5:f5c8b89ceac3f14e456577557df1ef40-2697049  ncbiprotein:YP_009724391   \n",
       "4  md5:375e0f905c315e06a99c80b736c125d2-2697049  ncbiprotein:YP_009724392   \n",
       "\n",
       "       geneAccession  \n",
       "0  ncbigene:43740578  \n",
       "1  ncbigene:43740578  \n",
       "2  ncbigene:43740568  \n",
       "3  ncbigene:43740569  \n",
       "4  ncbigene:43740570  "
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "df.head()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Save data for Knowledge Graph Import"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "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>id</th>\n",
       "      <th>geneName</th>\n",
       "      <th>geneStart</th>\n",
       "      <th>geneEnd</th>\n",
       "      <th>geneAccession</th>\n",
       "      <th>proteinAccession</th>\n",
       "      <th>proteinName</th>\n",
       "      <th>sequence</th>\n",
       "      <th>fullLength</th>\n",
       "      <th>taxonomyId</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>md5:04b37f827f46ae57d50f68cdc0544990-693998</td>\n",
       "      <td>NS3</td>\n",
       "      <td>23862</td>\n",
       "      <td>24551</td>\n",
       "      <td>ncbigene:5741385</td>\n",
       "      <td>ncbiprotein:YP_001552237</td>\n",
       "      <td>hypothetical protein</td>\n",
       "      <td>MFGGLFQLTLEGVVNSTIRAANLGPKDEAIIREQVQPAINGANMLG...</td>\n",
       "      <td>True</td>\n",
       "      <td>taxonomy:693998</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>md5:cd33bede1b54e7286126465bf03d5a06-693998</td>\n",
       "      <td>E</td>\n",
       "      <td>24532</td>\n",
       "      <td>24759</td>\n",
       "      <td>ncbigene:5741386</td>\n",
       "      <td>ncbiprotein:YP_001552238</td>\n",
       "      <td>small membrane protein</td>\n",
       "      <td>MFLKIVEDDGLFINTVLWLLVLILVLLVAITVIKLIQLCFSCHRLM...</td>\n",
       "      <td>True</td>\n",
       "      <td>taxonomy:693998</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5</th>\n",
       "      <td>md5:ec60b60cded44a243c07ca514f438faf-693998</td>\n",
       "      <td>M</td>\n",
       "      <td>24768</td>\n",
       "      <td>25457</td>\n",
       "      <td>ncbigene:5741387</td>\n",
       "      <td>ncbiprotein:YP_001552239</td>\n",
       "      <td>membrane glycoprotein</td>\n",
       "      <td>MSSNSTESVPVTEVLEHLRNWNFSWNIILTVFIAVLQYGNMKYSFF...</td>\n",
       "      <td>True</td>\n",
       "      <td>taxonomy:693998</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6</th>\n",
       "      <td>md5:eef9118db3bed6fa64197cf46bf2d8b9-693998</td>\n",
       "      <td>N</td>\n",
       "      <td>25469</td>\n",
       "      <td>26596</td>\n",
       "      <td>ncbigene:5741383</td>\n",
       "      <td>ncbiprotein:YP_001552240</td>\n",
       "      <td>nucleocapsid phosphoprotein</td>\n",
       "      <td>MATVNWGDAVEQAESRGRKRIPLSLFAPLRVTDGKNFWNVMPRNGV...</td>\n",
       "      <td>True</td>\n",
       "      <td>taxonomy:693998</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7</th>\n",
       "      <td>md5:7d6ed931572d7829c1a73a26e5c637fe-693998</td>\n",
       "      <td>NS7a</td>\n",
       "      <td>26608</td>\n",
       "      <td>26907</td>\n",
       "      <td>ncbigene:5741382</td>\n",
       "      <td>ncbiprotein:YP_001552241</td>\n",
       "      <td>hypothetical protein</td>\n",
       "      <td>MNQLIFFLCMMCCCAILFDWLFNLFFYACQVDTWQEFAFSCNWSWS...</td>\n",
       "      <td>True</td>\n",
       "      <td>taxonomy:693998</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "                                            id geneName  geneStart  geneEnd  \\\n",
       "3  md5:04b37f827f46ae57d50f68cdc0544990-693998      NS3      23862    24551   \n",
       "4  md5:cd33bede1b54e7286126465bf03d5a06-693998        E      24532    24759   \n",
       "5  md5:ec60b60cded44a243c07ca514f438faf-693998        M      24768    25457   \n",
       "6  md5:eef9118db3bed6fa64197cf46bf2d8b9-693998        N      25469    26596   \n",
       "7  md5:7d6ed931572d7829c1a73a26e5c637fe-693998     NS7a      26608    26907   \n",
       "\n",
       "      geneAccession          proteinAccession                  proteinName  \\\n",
       "3  ncbigene:5741385  ncbiprotein:YP_001552237         hypothetical protein   \n",
       "4  ncbigene:5741386  ncbiprotein:YP_001552238       small membrane protein   \n",
       "5  ncbigene:5741387  ncbiprotein:YP_001552239        membrane glycoprotein   \n",
       "6  ncbigene:5741383  ncbiprotein:YP_001552240  nucleocapsid phosphoprotein   \n",
       "7  ncbigene:5741382  ncbiprotein:YP_001552241         hypothetical protein   \n",
       "\n",
       "                                            sequence fullLength  \\\n",
       "3  MFGGLFQLTLEGVVNSTIRAANLGPKDEAIIREQVQPAINGANMLG...       True   \n",
       "4  MFLKIVEDDGLFINTVLWLLVLILVLLVAITVIKLIQLCFSCHRLM...       True   \n",
       "5  MSSNSTESVPVTEVLEHLRNWNFSWNIILTVFIAVLQYGNMKYSFF...       True   \n",
       "6  MATVNWGDAVEQAESRGRKRIPLSLFAPLRVTDGKNFWNVMPRNGV...       True   \n",
       "7  MNQLIFFLCMMCCCAILFDWLFNLFFYACQVDTWQEFAFSCNWSWS...       True   \n",
       "\n",
       "        taxonomyId  \n",
       "3  taxonomy:693998  \n",
       "4  taxonomy:693998  \n",
       "5  taxonomy:693998  \n",
       "6  taxonomy:693998  \n",
       "7  taxonomy:693998  "
      ]
     },
     "execution_count": 19,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "df = df[['id', 'geneName', 'geneStart', 'geneEnd', 'geneAccession', 'proteinAccession', 'proteinName', 'sequence', 'fullLength', 'taxonomyId']]\n",
    "df.tail()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### TODO find a premanent fix for missing data in GenBank records"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "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>id</th>\n",
       "      <th>geneName</th>\n",
       "      <th>geneStart</th>\n",
       "      <th>geneEnd</th>\n",
       "      <th>geneAccession</th>\n",
       "      <th>proteinAccession</th>\n",
       "      <th>proteinName</th>\n",
       "      <th>sequence</th>\n",
       "      <th>fullLength</th>\n",
       "      <th>taxonomyId</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>6</th>\n",
       "      <td>md5:6fc05e513ec08d8b439f73b2851bfba0-31631</td>\n",
       "      <td></td>\n",
       "      <td>28108</td>\n",
       "      <td>28362</td>\n",
       "      <td></td>\n",
       "      <td>ncbiprotein:YP_009555243</td>\n",
       "      <td>envelope protein</td>\n",
       "      <td>MFMADAYLADTVWYVGQIIFIVAICLLVTIVVVAFLATFKLCIQLC...</td>\n",
       "      <td>True</td>\n",
       "      <td>taxonomy:31631</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>md5:cbdd4c2307a04becd1132db6820e61d8-693998</td>\n",
       "      <td></td>\n",
       "      <td>297</td>\n",
       "      <td>20479</td>\n",
       "      <td></td>\n",
       "      <td>ncbiprotein:YP_001552234</td>\n",
       "      <td>orf1ab polyprotein</td>\n",
       "      <td>MSINQLTLAVASDQEISAHGYPTMSDAVEHFSSSASHGFKDCRFVA...</td>\n",
       "      <td>True</td>\n",
       "      <td>taxonomy:693998</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "                                            id geneName  geneStart  geneEnd  \\\n",
       "6   md5:6fc05e513ec08d8b439f73b2851bfba0-31631               28108    28362   \n",
       "0  md5:cbdd4c2307a04becd1132db6820e61d8-693998                 297    20479   \n",
       "\n",
       "  geneAccession          proteinAccession         proteinName  \\\n",
       "6                ncbiprotein:YP_009555243    envelope protein   \n",
       "0                ncbiprotein:YP_001552234  orf1ab polyprotein   \n",
       "\n",
       "                                            sequence fullLength  \\\n",
       "6  MFMADAYLADTVWYVGQIIFIVAICLLVTIVVVAFLATFKLCIQLC...       True   \n",
       "0  MSINQLTLAVASDQEISAHGYPTMSDAVEHFSSSASHGFKDCRFVA...       True   \n",
       "\n",
       "        taxonomyId  \n",
       "6   taxonomy:31631  \n",
       "0  taxonomy:693998  "
      ]
     },
     "execution_count": 20,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "df.query(\"geneAccession == ''\")"
   ]
  },
  {
   "cell_type": "raw",
   "metadata": {},
   "source": [
    "     gene            28108..28362\n",
    "                     /gene=\"E\"\n",
    "                     /locus_tag=\"EYW02_gp6\"\n",
    "                     /db_xref=\"GeneID:39105217\"\n",
    "     CDS             28108..28362\n",
    "                     /codon_start=1\n",
    "                     /product=\"envelope protein\"\n",
    "                     /protein_id=\"YP_009555243.1\"\n",
    "                     /translation=\"MFMADAYLADTVWYVGQIIFIVAICLLVTIVVVAFLATFKLCIQ\n",
    "                     LCGMCNTLVLSPSIYVFNRGRQFYEFYNDVKPPVLDVDDV\""
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [],
   "source": [
    "df.reset_index(drop=True,inplace=True)\n",
    "index = df.query(\"proteinAccession == 'ncbiprotein:YP_009555243'\").index.values[0]\n",
    "df.at[index, 'geneName'] = 'E'\n",
    "df.at[index, 'geneAccession'] = 'ncbigene:39105217'"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "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>id</th>\n",
       "      <th>geneName</th>\n",
       "      <th>geneStart</th>\n",
       "      <th>geneEnd</th>\n",
       "      <th>geneAccession</th>\n",
       "      <th>proteinAccession</th>\n",
       "      <th>proteinName</th>\n",
       "      <th>sequence</th>\n",
       "      <th>fullLength</th>\n",
       "      <th>taxonomyId</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "Empty DataFrame\n",
       "Columns: [id, geneName, geneStart, geneEnd, geneAccession, proteinAccession, proteinName, sequence, fullLength, taxonomyId]\n",
       "Index: []"
      ]
     },
     "execution_count": 22,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "df.query(\"sequence == ''\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "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>id</th>\n",
       "      <th>geneName</th>\n",
       "      <th>geneStart</th>\n",
       "      <th>geneEnd</th>\n",
       "      <th>geneAccession</th>\n",
       "      <th>proteinAccession</th>\n",
       "      <th>proteinName</th>\n",
       "      <th>sequence</th>\n",
       "      <th>fullLength</th>\n",
       "      <th>taxonomyId</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>157</th>\n",
       "      <td>md5:6fc05e513ec08d8b439f73b2851bfba0-31631</td>\n",
       "      <td>E</td>\n",
       "      <td>28108</td>\n",
       "      <td>28362</td>\n",
       "      <td>ncbigene:39105217</td>\n",
       "      <td>ncbiprotein:YP_009555243</td>\n",
       "      <td>envelope protein</td>\n",
       "      <td>MFMADAYLADTVWYVGQIIFIVAICLLVTIVVVAFLATFKLCIQLC...</td>\n",
       "      <td>True</td>\n",
       "      <td>taxonomy:31631</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "                                             id geneName  geneStart  geneEnd  \\\n",
       "157  md5:6fc05e513ec08d8b439f73b2851bfba0-31631        E      28108    28362   \n",
       "\n",
       "         geneAccession          proteinAccession       proteinName  \\\n",
       "157  ncbigene:39105217  ncbiprotein:YP_009555243  envelope protein   \n",
       "\n",
       "                                              sequence fullLength  \\\n",
       "157  MFMADAYLADTVWYVGQIIFIVAICLLVTIVVVAFLATFKLCIQLC...       True   \n",
       "\n",
       "         taxonomyId  \n",
       "157  taxonomy:31631  "
      ]
     },
     "execution_count": 23,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "df.query(\"proteinAccession == 'ncbiprotein:YP_009555243'\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "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>id</th>\n",
       "      <th>geneName</th>\n",
       "      <th>geneStart</th>\n",
       "      <th>geneEnd</th>\n",
       "      <th>geneAccession</th>\n",
       "      <th>proteinAccession</th>\n",
       "      <th>proteinName</th>\n",
       "      <th>sequence</th>\n",
       "      <th>fullLength</th>\n",
       "      <th>taxonomyId</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>77</th>\n",
       "      <td>md5:5b39d604648221117b0614889d3d6067-694009</td>\n",
       "      <td>S</td>\n",
       "      <td>21492</td>\n",
       "      <td>25259</td>\n",
       "      <td>ncbigene:1489668</td>\n",
       "      <td>ncbiprotein:YP_009825051</td>\n",
       "      <td>spike glycoprotein</td>\n",
       "      <td>MFIFLLFLTLTSGSDLDRCTTFDDVQAPNYTQHTSSMRGVYYPDEI...</td>\n",
       "      <td>True</td>\n",
       "      <td>taxonomy:694009</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "                                             id geneName  geneStart  geneEnd  \\\n",
       "77  md5:5b39d604648221117b0614889d3d6067-694009        S      21492    25259   \n",
       "\n",
       "       geneAccession          proteinAccession         proteinName  \\\n",
       "77  ncbigene:1489668  ncbiprotein:YP_009825051  spike glycoprotein   \n",
       "\n",
       "                                             sequence fullLength  \\\n",
       "77  MFIFLLFLTLTSGSDLDRCTTFDDVQAPNYTQHTSSMRGVYYPDEI...       True   \n",
       "\n",
       "         taxonomyId  \n",
       "77  taxonomy:694009  "
      ]
     },
     "execution_count": 24,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "df.query(\"proteinAccession == 'ncbiprotein:YP_009825051'\") # protein name not connected to protein"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "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>id</th>\n",
       "      <th>geneName</th>\n",
       "      <th>geneStart</th>\n",
       "      <th>geneEnd</th>\n",
       "      <th>geneAccession</th>\n",
       "      <th>proteinAccession</th>\n",
       "      <th>proteinName</th>\n",
       "      <th>sequence</th>\n",
       "      <th>fullLength</th>\n",
       "      <th>taxonomyId</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>155</th>\n",
       "      <td>md5:a9fc2c01d20938618b58f690ae776309-31631</td>\n",
       "      <td>S</td>\n",
       "      <td>23643</td>\n",
       "      <td>27704</td>\n",
       "      <td>ncbigene:39105218</td>\n",
       "      <td>ncbiprotein:YP_009555241</td>\n",
       "      <td>spike surface glycoprotein</td>\n",
       "      <td>MFLILLISLPTAFAVIGDLKCTSDNINDKDTGPPPISTDTVDVTNG...</td>\n",
       "      <td>True</td>\n",
       "      <td>taxonomy:31631</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "                                             id geneName  geneStart  geneEnd  \\\n",
       "155  md5:a9fc2c01d20938618b58f690ae776309-31631        S      23643    27704   \n",
       "\n",
       "         geneAccession          proteinAccession                 proteinName  \\\n",
       "155  ncbigene:39105218  ncbiprotein:YP_009555241  spike surface glycoprotein   \n",
       "\n",
       "                                              sequence fullLength  \\\n",
       "155  MFLILLISLPTAFAVIGDLKCTSDNINDKDTGPPPISTDTVDVTNG...       True   \n",
       "\n",
       "         taxonomyId  \n",
       "155  taxonomy:31631  "
      ]
     },
     "execution_count": 25,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "df.query(\"proteinAccession == 'ncbiprotein:YP_009555241'\") # protein name not connected to protein"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [],
   "source": [
    "df.to_csv(NEO4J_IMPORT / \"01b-NCBIGeneProtein.csv\", index=False)"
   ]
  }
 ],
 "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.8"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}