08 - Introdução Biopython, Seq, Seqrecord.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Introdução\n",
    "\n",
    "### O que é o Biopython?\n",
    "\n",
    "O Biopython é uma biblioteca que facilita o uso da linguagem python para sua aplicação na bioinformática. O site do [Biopython](http://www.biopython.org/) fornece toda a documentação necessária para entender e exemplos para aplicar os módulos e classes da biblioteca.\n",
    "\n",
    "### O que tem no pacote Biopyhon\n",
    "\n",
    "O Biopython tem diversas funcionalidades, entre elas temos:\n",
    "\n",
    "1. Parsers para os principais tipos de arquivos utilizados em bioinformática\n",
    "   - Blast output\n",
    "   - Clustalw\n",
    "   - FASTA\n",
    "   - GenBank\n",
    "   - PubMed and Medline\n",
    "   - ExPASy\n",
    "   - SCOP\n",
    "   - UniGene\n",
    "   - SwissProt\n",
    "\n",
    "2. Suporte para executar ferramentas online e standalone de bioinformática:\n",
    "   - Blast (online, standalone)\n",
    "   - ExPASy (online)\n",
    "   - Clustalw\n",
    "   - EMBOSS\n",
    "\n",
    "3. Operações com arquivos de sequências (transcrição, tradução, cálculos de pesos)\n",
    "\n",
    "4. Algoritmos de classificação de dados (K-Nearest, Naive Bayes ou Support Vector Machines)\n",
    "\n",
    "E estas são apenas algumas das funcionalidades do Biopython.\n",
    "\n",
    "### Instalação do Biopython\n",
    "\n",
    "A forma mais fácil de instalar o Biopython é utilizando o gerenciador de pacotes Anaconda:\n",
    "\n",
    "``` conda install -c anaconda biopython```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1. Classe Seq\n",
    "\n",
    "Neste etapa vamos aprender como trabalhar com as sequências, que são um dos principais tipo de dados que usamos em bioinformática. A primeira pergunta que pode ocorrer é por que ter uma classe Seq e não trabalhar com uma string do python? E a resposta é simples, é porque elas são classes diferentes. Apesar da classe `Seq` ter quase todos os métodos da classe string, ela possuí alguns métodos diferentes, como o translate() que faz a tradução biológica e o método reverse_complement() que também é bastante útil. Outra diferença importante é que a classe `Seq` possuí um atributo Alphabet() que indica se a sequência é de DNA, RNA ou aminoácios.\n",
    "\n",
    "### 1.1 Sequências e Alfabetos\n",
    "\n",
    "Vamos a alguns exemplos de como criar sequências usando biopython:\n",
    "\n",
    "Inicialmente vamos importar a classe `Seq` do Biopython e em seguida definir uma sequência."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Help on class Seq in module Bio.Seq:\n",
      "\n",
      "class Seq(builtins.object)\n",
      " |  Seq(data, alphabet=Alphabet())\n",
      " |  \n",
      " |  Read-only sequence object (essentially a string with an alphabet).\n",
      " |  \n",
      " |  Like normal python strings, our basic sequence object is immutable.\n",
      " |  This prevents you from doing my_seq[5] = \"A\" for example, but does allow\n",
      " |  Seq objects to be used as dictionary keys.\n",
      " |  \n",
      " |  The Seq object provides a number of string like methods (such as count,\n",
      " |  find, split and strip), which are alphabet aware where appropriate.\n",
      " |  \n",
      " |  In addition to the string like sequence, the Seq object has an alphabet\n",
      " |  property. This is an instance of an Alphabet class from Bio.Alphabet,\n",
      " |  for example generic DNA, or IUPAC DNA. This describes the type of molecule\n",
      " |  (e.g. RNA, DNA, protein) and may also indicate the expected symbols\n",
      " |  (letters).\n",
      " |  \n",
      " |  The Seq object also provides some biological methods, such as complement,\n",
      " |  reverse_complement, transcribe, back_transcribe and translate (which are\n",
      " |  not applicable to sequences with a protein alphabet).\n",
      " |  \n",
      " |  Methods defined here:\n",
      " |  \n",
      " |  __add__(self, other)\n",
      " |      Add another sequence or string to this sequence.\n",
      " |      \n",
      " |      If adding a string to a Seq, the alphabet is preserved:\n",
      " |      \n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> from Bio.Alphabet import generic_protein\n",
      " |      >>> Seq(\"MELKI\", generic_protein) + \"LV\"\n",
      " |      Seq('MELKILV', ProteinAlphabet())\n",
      " |      \n",
      " |      When adding two Seq (like) objects, the alphabets are important.\n",
      " |      Consider this example:\n",
      " |      \n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> from Bio.Alphabet.IUPAC import unambiguous_dna, ambiguous_dna\n",
      " |      >>> unamb_dna_seq = Seq(\"ACGT\", unambiguous_dna)\n",
      " |      >>> ambig_dna_seq = Seq(\"ACRGT\", ambiguous_dna)\n",
      " |      >>> unamb_dna_seq\n",
      " |      Seq('ACGT', IUPACUnambiguousDNA())\n",
      " |      >>> ambig_dna_seq\n",
      " |      Seq('ACRGT', IUPACAmbiguousDNA())\n",
      " |      \n",
      " |      If we add the ambiguous and unambiguous IUPAC DNA alphabets, we get\n",
      " |      the more general ambiguous IUPAC DNA alphabet:\n",
      " |      \n",
      " |      >>> unamb_dna_seq + ambig_dna_seq\n",
      " |      Seq('ACGTACRGT', IUPACAmbiguousDNA())\n",
      " |      \n",
      " |      However, if the default generic alphabet is included, the result is\n",
      " |      a generic alphabet:\n",
      " |      \n",
      " |      >>> Seq(\"\") + ambig_dna_seq\n",
      " |      Seq('ACRGT')\n",
      " |      \n",
      " |      You can't add RNA and DNA sequences:\n",
      " |      \n",
      " |      >>> from Bio.Alphabet import generic_dna, generic_rna\n",
      " |      >>> Seq(\"ACGT\", generic_dna) + Seq(\"ACGU\", generic_rna)\n",
      " |      Traceback (most recent call last):\n",
      " |         ...\n",
      " |      TypeError: Incompatible alphabets DNAAlphabet() and RNAAlphabet()\n",
      " |      \n",
      " |      You can't add nucleotide and protein sequences:\n",
      " |      \n",
      " |      >>> from Bio.Alphabet import generic_dna, generic_protein\n",
      " |      >>> Seq(\"ACGT\", generic_dna) + Seq(\"MELKI\", generic_protein)\n",
      " |      Traceback (most recent call last):\n",
      " |         ...\n",
      " |      TypeError: Incompatible alphabets DNAAlphabet() and ProteinAlphabet()\n",
      " |  \n",
      " |  __contains__(self, char)\n",
      " |      Implement the 'in' keyword, like a python string.\n",
      " |      \n",
      " |      e.g.\n",
      " |      \n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> from Bio.Alphabet import generic_dna, generic_rna, generic_protein\n",
      " |      >>> my_dna = Seq(\"ATATGAAATTTGAAAA\", generic_dna)\n",
      " |      >>> \"AAA\" in my_dna\n",
      " |      True\n",
      " |      >>> Seq(\"AAA\") in my_dna\n",
      " |      True\n",
      " |      >>> Seq(\"AAA\", generic_dna) in my_dna\n",
      " |      True\n",
      " |      \n",
      " |      Like other Seq methods, this will raise a type error if another Seq\n",
      " |      (or Seq like) object with an incompatible alphabet is used:\n",
      " |      \n",
      " |      >>> Seq(\"AAA\", generic_rna) in my_dna\n",
      " |      Traceback (most recent call last):\n",
      " |         ...\n",
      " |      TypeError: Incompatible alphabets DNAAlphabet() and RNAAlphabet()\n",
      " |      >>> Seq(\"AAA\", generic_protein) in my_dna\n",
      " |      Traceback (most recent call last):\n",
      " |         ...\n",
      " |      TypeError: Incompatible alphabets DNAAlphabet() and ProteinAlphabet()\n",
      " |  \n",
      " |  __eq__(self, other)\n",
      " |      Compare the sequence to another sequence or a string (README).\n",
      " |      \n",
      " |      Historically comparing Seq objects has done Python object comparison.\n",
      " |      After considerable discussion (keeping in mind constraints of the\n",
      " |      Python language, hashes and dictionary support), Biopython now uses\n",
      " |      simple string comparison (with a warning about the change).\n",
      " |      \n",
      " |      Note that incompatible alphabets (e.g. DNA to RNA) will trigger a\n",
      " |      warning.\n",
      " |      \n",
      " |      During this transition period, please just do explicit comparisons:\n",
      " |      \n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> from Bio.Alphabet import generic_dna\n",
      " |      >>> seq1 = Seq(\"ACGT\")\n",
      " |      >>> seq2 = Seq(\"ACGT\")\n",
      " |      >>> id(seq1) == id(seq2)\n",
      " |      False\n",
      " |      >>> str(seq1) == str(seq2)\n",
      " |      True\n",
      " |      \n",
      " |      The new behaviour is to use string-like equality:\n",
      " |      \n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> from Bio.Alphabet import generic_dna\n",
      " |      >>> seq1 == seq2\n",
      " |      True\n",
      " |      >>> seq1 == \"ACGT\"\n",
      " |      True\n",
      " |      >>> seq1 == Seq(\"ACGT\", generic_dna)\n",
      " |      True\n",
      " |  \n",
      " |  __ge__(self, other)\n",
      " |      Implement the greater-than or equal operand.\n",
      " |  \n",
      " |  __getitem__(self, index)\n",
      " |      Return a subsequence of single letter, use my_seq[index].\n",
      " |      \n",
      " |      >>> my_seq = Seq('ACTCGACGTCG')\n",
      " |      >>> my_seq[5]\n",
      " |      'A'\n",
      " |  \n",
      " |  __gt__(self, other)\n",
      " |      Implement the greater-than operand.\n",
      " |  \n",
      " |  __hash__(self)\n",
      " |      Hash for comparison.\n",
      " |      \n",
      " |      See the __cmp__ documentation - this has changed from past\n",
      " |      versions of Biopython!\n",
      " |  \n",
      " |  __imul__(self, other)\n",
      " |      Multiply Seq in-place.\n",
      " |      \n",
      " |      Note although Seq is immutable, the in-place method is\n",
      " |      included to match the behaviour for regular Python strings.\n",
      " |      \n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> from Bio.Alphabet import generic_dna\n",
      " |      >>> seq = Seq('ATG', generic_dna)\n",
      " |      >>> seq *= 2\n",
      " |      >>> seq\n",
      " |      Seq('ATGATG', DNAAlphabet())\n",
      " |  \n",
      " |  __init__(self, data, alphabet=Alphabet())\n",
      " |      Create a Seq object.\n",
      " |      \n",
      " |      Arguments:\n",
      " |       - seq - Sequence, required (string)\n",
      " |       - alphabet - Optional argument, an Alphabet object from\n",
      " |         Bio.Alphabet\n",
      " |      \n",
      " |      You will typically use Bio.SeqIO to read in sequences from files as\n",
      " |      SeqRecord objects, whose sequence will be exposed as a Seq object via\n",
      " |      the seq property.\n",
      " |      \n",
      " |      However, will often want to create your own Seq objects directly:\n",
      " |      \n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> from Bio.Alphabet import IUPAC\n",
      " |      >>> my_seq = Seq(\"MKQHKAMIVALIVICITAVVAALVTRKDLCEVHIRTGQTEVAVF\",\n",
      " |      ...              IUPAC.protein)\n",
      " |      >>> my_seq\n",
      " |      Seq('MKQHKAMIVALIVICITAVVAALVTRKDLCEVHIRTGQTEVAVF', IUPACProtein())\n",
      " |      >>> print(my_seq)\n",
      " |      MKQHKAMIVALIVICITAVVAALVTRKDLCEVHIRTGQTEVAVF\n",
      " |      >>> my_seq.alphabet\n",
      " |      IUPACProtein()\n",
      " |  \n",
      " |  __le__(self, other)\n",
      " |      Implement the less-than or equal operand.\n",
      " |  \n",
      " |  __len__(self)\n",
      " |      Return the length of the sequence, use len(my_seq).\n",
      " |  \n",
      " |  __lt__(self, other)\n",
      " |      Implement the less-than operand.\n",
      " |  \n",
      " |  __mul__(self, other)\n",
      " |      Multiply Seq by integer.\n",
      " |      \n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> from Bio.Alphabet import generic_dna\n",
      " |      >>> Seq('ATG') * 2\n",
      " |      Seq('ATGATG')\n",
      " |      >>> Seq('ATG', generic_dna) * 2\n",
      " |      Seq('ATGATG', DNAAlphabet())\n",
      " |  \n",
      " |  __ne__(self, other)\n",
      " |      Implement the not-equal operand.\n",
      " |  \n",
      " |  __radd__(self, other)\n",
      " |      Add a sequence on the left.\n",
      " |      \n",
      " |      If adding a string to a Seq, the alphabet is preserved:\n",
      " |      \n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> from Bio.Alphabet import generic_protein\n",
      " |      >>> \"LV\" + Seq(\"MELKI\", generic_protein)\n",
      " |      Seq('LVMELKI', ProteinAlphabet())\n",
      " |      \n",
      " |      Adding two Seq (like) objects is handled via the __add__ method.\n",
      " |  \n",
      " |  __repr__(self)\n",
      " |      Return (truncated) representation of the sequence for debugging.\n",
      " |  \n",
      " |  __rmul__(self, other)\n",
      " |      Multiply integer by Seq.\n",
      " |      \n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> from Bio.Alphabet import generic_dna\n",
      " |      >>> 2 * Seq('ATG')\n",
      " |      Seq('ATGATG')\n",
      " |      >>> 2 * Seq('ATG', generic_dna)\n",
      " |      Seq('ATGATG', DNAAlphabet())\n",
      " |  \n",
      " |  __str__(self)\n",
      " |      Return the full sequence as a python string, use str(my_seq).\n",
      " |      \n",
      " |      Note that Biopython 1.44 and earlier would give a truncated\n",
      " |      version of repr(my_seq) for str(my_seq).  If you are writing code\n",
      " |      which need to be backwards compatible with old Biopython, you\n",
      " |      should continue to use my_seq.tostring() rather than str(my_seq).\n",
      " |  \n",
      " |  back_transcribe(self)\n",
      " |      Return the DNA sequence from an RNA sequence by creating a new Seq object.\n",
      " |      \n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> from Bio.Alphabet import IUPAC\n",
      " |      >>> messenger_rna = Seq(\"AUGGCCAUUGUAAUGGGCCGCUGAAAGGGUGCCCGAUAG\",\n",
      " |      ...                     IUPAC.unambiguous_rna)\n",
      " |      >>> messenger_rna\n",
      " |      Seq('AUGGCCAUUGUAAUGGGCCGCUGAAAGGGUGCCCGAUAG', IUPACUnambiguousRNA())\n",
      " |      >>> messenger_rna.back_transcribe()\n",
      " |      Seq('ATGGCCATTGTAATGGGCCGCTGAAAGGGTGCCCGATAG', IUPACUnambiguousDNA())\n",
      " |      \n",
      " |      Trying to back-transcribe a protein or DNA sequence raises an\n",
      " |      exception:\n",
      " |      \n",
      " |      >>> my_protein = Seq(\"MAIVMGR\", IUPAC.protein)\n",
      " |      >>> my_protein.back_transcribe()\n",
      " |      Traceback (most recent call last):\n",
      " |         ...\n",
      " |      ValueError: Proteins cannot be back transcribed!\n",
      " |  \n",
      " |  complement(self)\n",
      " |      Return the complement sequence by creating a new Seq object.\n",
      " |      \n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> from Bio.Alphabet import IUPAC\n",
      " |      >>> my_dna = Seq(\"CCCCCGATAG\", IUPAC.unambiguous_dna)\n",
      " |      >>> my_dna\n",
      " |      Seq('CCCCCGATAG', IUPACUnambiguousDNA())\n",
      " |      >>> my_dna.complement()\n",
      " |      Seq('GGGGGCTATC', IUPACUnambiguousDNA())\n",
      " |      \n",
      " |      You can of course used mixed case sequences,\n",
      " |      \n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> from Bio.Alphabet import generic_dna\n",
      " |      >>> my_dna = Seq(\"CCCCCgatA-GD\", generic_dna)\n",
      " |      >>> my_dna\n",
      " |      Seq('CCCCCgatA-GD', DNAAlphabet())\n",
      " |      >>> my_dna.complement()\n",
      " |      Seq('GGGGGctaT-CH', DNAAlphabet())\n",
      " |      \n",
      " |      Note in the above example, ambiguous character D denotes\n",
      " |      G, A or T so its complement is H (for C, T or A).\n",
      " |      \n",
      " |      Trying to complement a protein sequence raises an exception.\n",
      " |      \n",
      " |      >>> my_protein = Seq(\"MAIVMGR\", IUPAC.protein)\n",
      " |      >>> my_protein.complement()\n",
      " |      Traceback (most recent call last):\n",
      " |         ...\n",
      " |      ValueError: Proteins do not have complements!\n",
      " |  \n",
      " |  count(self, sub, start=0, end=9223372036854775807)\n",
      " |      Return a non-overlapping count, like that of a python string.\n",
      " |      \n",
      " |      This behaves like the python string method of the same name,\n",
      " |      which does a non-overlapping count!\n",
      " |      \n",
      " |      For an overlapping search use the newer count_overlap() method.\n",
      " |      \n",
      " |      Returns an integer, the number of occurrences of substring\n",
      " |      argument sub in the (sub)sequence given by [start:end].\n",
      " |      Optional arguments start and end are interpreted as in slice\n",
      " |      notation.\n",
      " |      \n",
      " |      Arguments:\n",
      " |       - sub - a string or another Seq object to look for\n",
      " |       - start - optional integer, slice start\n",
      " |       - end - optional integer, slice end\n",
      " |      \n",
      " |      e.g.\n",
      " |      \n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> my_seq = Seq(\"AAAATGA\")\n",
      " |      >>> print(my_seq.count(\"A\"))\n",
      " |      5\n",
      " |      >>> print(my_seq.count(\"ATG\"))\n",
      " |      1\n",
      " |      >>> print(my_seq.count(Seq(\"AT\")))\n",
      " |      1\n",
      " |      >>> print(my_seq.count(\"AT\", 2, -1))\n",
      " |      1\n",
      " |      \n",
      " |      HOWEVER, please note because python strings and Seq objects (and\n",
      " |      MutableSeq objects) do a non-overlapping search, this may not give\n",
      " |      the answer you expect:\n",
      " |      \n",
      " |      >>> \"AAAA\".count(\"AA\")\n",
      " |      2\n",
      " |      >>> print(Seq(\"AAAA\").count(\"AA\"))\n",
      " |      2\n",
      " |      \n",
      " |      An overlapping search, as implemented in .count_overlap(),\n",
      " |      would give the answer as three!\n",
      " |  \n",
      " |  count_overlap(self, sub, start=0, end=9223372036854775807)\n",
      " |      Return an overlapping count.\n",
      " |      \n",
      " |      For a non-overlapping search use the count() method.\n",
      " |      \n",
      " |      Returns an integer, the number of occurrences of substring\n",
      " |      argument sub in the (sub)sequence given by [start:end].\n",
      " |      Optional arguments start and end are interpreted as in slice\n",
      " |      notation.\n",
      " |      \n",
      " |      Arguments:\n",
      " |       - sub - a string or another Seq object to look for\n",
      " |       - start - optional integer, slice start\n",
      " |       - end - optional integer, slice end\n",
      " |      \n",
      " |      e.g.\n",
      " |      \n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> print(Seq(\"AAAA\").count_overlap(\"AA\"))\n",
      " |      3\n",
      " |      >>> print(Seq(\"ATATATATA\").count_overlap(\"ATA\"))\n",
      " |      4\n",
      " |      >>> print(Seq(\"ATATATATA\").count_overlap(\"ATA\", 3, -1))\n",
      " |      1\n",
      " |      \n",
      " |      Where substrings do not overlap, should behave the same as\n",
      " |      the count() method:\n",
      " |      \n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> my_seq = Seq(\"AAAATGA\")\n",
      " |      >>> print(my_seq.count_overlap(\"A\"))\n",
      " |      5\n",
      " |      >>> my_seq.count_overlap(\"A\") == my_seq.count(\"A\")\n",
      " |      True\n",
      " |      >>> print(my_seq.count_overlap(\"ATG\"))\n",
      " |      1\n",
      " |      >>> my_seq.count_overlap(\"ATG\") == my_seq.count(\"ATG\")\n",
      " |      True\n",
      " |      >>> print(my_seq.count_overlap(Seq(\"AT\")))\n",
      " |      1\n",
      " |      >>> my_seq.count_overlap(Seq(\"AT\")) == my_seq.count(Seq(\"AT\"))\n",
      " |      True\n",
      " |      >>> print(my_seq.count_overlap(\"AT\", 2, -1))\n",
      " |      1\n",
      " |      >>> my_seq.count_overlap(\"AT\", 2, -1) == my_seq.count(\"AT\", 2, -1)\n",
      " |      True\n",
      " |      \n",
      " |      HOWEVER, do not use this method for such cases because the\n",
      " |      count() method is much for efficient.\n",
      " |  \n",
      " |  endswith(self, suffix, start=0, end=9223372036854775807)\n",
      " |      Return True if the Seq ends with the given suffix, False otherwise.\n",
      " |      \n",
      " |      This behaves like the python string method of the same name.\n",
      " |      \n",
      " |      Return True if the sequence ends with the specified suffix\n",
      " |      (a string or another Seq object), False otherwise.\n",
      " |      With optional start, test sequence beginning at that position.\n",
      " |      With optional end, stop comparing sequence at that position.\n",
      " |      suffix can also be a tuple of strings to try.  e.g.\n",
      " |      \n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> my_rna = Seq(\"GUCAUGGCCAUUGUAAUGGGCCGCUGAAAGGGUGCCCGAUAGUUG\")\n",
      " |      >>> my_rna.endswith(\"UUG\")\n",
      " |      True\n",
      " |      >>> my_rna.endswith(\"AUG\")\n",
      " |      False\n",
      " |      >>> my_rna.endswith(\"AUG\", 0, 18)\n",
      " |      True\n",
      " |      >>> my_rna.endswith((\"UCC\", \"UCA\", \"UUG\"))\n",
      " |      True\n",
      " |  \n",
      " |  find(self, sub, start=0, end=9223372036854775807)\n",
      " |      Find method, like that of a python string.\n",
      " |      \n",
      " |      This behaves like the python string method of the same name.\n",
      " |      \n",
      " |      Returns an integer, the index of the first occurrence of substring\n",
      " |      argument sub in the (sub)sequence given by [start:end].\n",
      " |      \n",
      " |      Arguments:\n",
      " |       - sub - a string or another Seq object to look for\n",
      " |       - start - optional integer, slice start\n",
      " |       - end - optional integer, slice end\n",
      " |      \n",
      " |      Returns -1 if the subsequence is NOT found.\n",
      " |      \n",
      " |      e.g. Locating the first typical start codon, AUG, in an RNA sequence:\n",
      " |      \n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> my_rna = Seq(\"GUCAUGGCCAUUGUAAUGGGCCGCUGAAAGGGUGCCCGAUAGUUG\")\n",
      " |      >>> my_rna.find(\"AUG\")\n",
      " |      3\n",
      " |  \n",
      " |  lower(self)\n",
      " |      Return a lower case copy of the sequence.\n",
      " |      \n",
      " |      This will adjust the alphabet if required. Note that the IUPAC\n",
      " |      alphabets are upper case only, and thus a generic alphabet must be\n",
      " |      substituted.\n",
      " |      \n",
      " |      >>> from Bio.Alphabet import Gapped, generic_dna\n",
      " |      >>> from Bio.Alphabet import IUPAC\n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> my_seq = Seq(\"CGGTACGCTTATGTCACGTAG*AAAAAA\",\n",
      " |      ...          Gapped(IUPAC.unambiguous_dna, \"*\"))\n",
      " |      >>> my_seq\n",
      " |      Seq('CGGTACGCTTATGTCACGTAG*AAAAAA', Gapped(IUPACUnambiguousDNA(), '*'))\n",
      " |      >>> my_seq.lower()\n",
      " |      Seq('cggtacgcttatgtcacgtag*aaaaaa', Gapped(DNAAlphabet(), '*'))\n",
      " |      \n",
      " |      See also the upper method.\n",
      " |  \n",
      " |  lstrip(self, chars=None)\n",
      " |      Return a new Seq object with leading (left) end stripped.\n",
      " |      \n",
      " |      This behaves like the python string method of the same name.\n",
      " |      \n",
      " |      Optional argument chars defines which characters to remove.  If\n",
      " |      omitted or None (default) then as for the python string method,\n",
      " |      this defaults to removing any white space.\n",
      " |      \n",
      " |      e.g. print(my_seq.lstrip(\"-\"))\n",
      " |      \n",
      " |      See also the strip and rstrip methods.\n",
      " |  \n",
      " |  reverse_complement(self)\n",
      " |      Return the reverse complement sequence by creating a new Seq object.\n",
      " |      \n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> from Bio.Alphabet import IUPAC\n",
      " |      >>> my_dna = Seq(\"CCCCCGATAGNR\", IUPAC.ambiguous_dna)\n",
      " |      >>> my_dna\n",
      " |      Seq('CCCCCGATAGNR', IUPACAmbiguousDNA())\n",
      " |      >>> my_dna.reverse_complement()\n",
      " |      Seq('YNCTATCGGGGG', IUPACAmbiguousDNA())\n",
      " |      \n",
      " |      Note in the above example, since R = G or A, its complement\n",
      " |      is Y (which denotes C or T).\n",
      " |      \n",
      " |      You can of course used mixed case sequences,\n",
      " |      \n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> from Bio.Alphabet import generic_dna\n",
      " |      >>> my_dna = Seq(\"CCCCCgatA-G\", generic_dna)\n",
      " |      >>> my_dna\n",
      " |      Seq('CCCCCgatA-G', DNAAlphabet())\n",
      " |      >>> my_dna.reverse_complement()\n",
      " |      Seq('C-TatcGGGGG', DNAAlphabet())\n",
      " |      \n",
      " |      Trying to complement a protein sequence raises an exception:\n",
      " |      \n",
      " |      >>> my_protein = Seq(\"MAIVMGR\", IUPAC.protein)\n",
      " |      >>> my_protein.reverse_complement()\n",
      " |      Traceback (most recent call last):\n",
      " |         ...\n",
      " |      ValueError: Proteins do not have complements!\n",
      " |  \n",
      " |  rfind(self, sub, start=0, end=9223372036854775807)\n",
      " |      Find from right method, like that of a python string.\n",
      " |      \n",
      " |      This behaves like the python string method of the same name.\n",
      " |      \n",
      " |      Returns an integer, the index of the last (right most) occurrence of\n",
      " |      substring argument sub in the (sub)sequence given by [start:end].\n",
      " |      \n",
      " |      Arguments:\n",
      " |       - sub - a string or another Seq object to look for\n",
      " |       - start - optional integer, slice start\n",
      " |       - end - optional integer, slice end\n",
      " |      \n",
      " |      Returns -1 if the subsequence is NOT found.\n",
      " |      \n",
      " |      e.g. Locating the last typical start codon, AUG, in an RNA sequence:\n",
      " |      \n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> my_rna = Seq(\"GUCAUGGCCAUUGUAAUGGGCCGCUGAAAGGGUGCCCGAUAGUUG\")\n",
      " |      >>> my_rna.rfind(\"AUG\")\n",
      " |      15\n",
      " |  \n",
      " |  rsplit(self, sep=None, maxsplit=-1)\n",
      " |      Do a right split method, like that of a python string.\n",
      " |      \n",
      " |      This behaves like the python string method of the same name.\n",
      " |      \n",
      " |      Return a list of the 'words' in the string (as Seq objects),\n",
      " |      using sep as the delimiter string.  If maxsplit is given, at\n",
      " |      most maxsplit splits are done COUNTING FROM THE RIGHT.\n",
      " |      If maxsplit is omitted, all splits are made.\n",
      " |      \n",
      " |      Following the python string method, sep will by default be any\n",
      " |      white space (tabs, spaces, newlines) but this is unlikely to\n",
      " |      apply to biological sequences.\n",
      " |      \n",
      " |      e.g. print(my_seq.rsplit(\"*\",1))\n",
      " |      \n",
      " |      See also the split method.\n",
      " |  \n",
      " |  rstrip(self, chars=None)\n",
      " |      Return a new Seq object with trailing (right) end stripped.\n",
      " |      \n",
      " |      This behaves like the python string method of the same name.\n",
      " |      \n",
      " |      Optional argument chars defines which characters to remove.  If\n",
      " |      omitted or None (default) then as for the python string method,\n",
      " |      this defaults to removing any white space.\n",
      " |      \n",
      " |      e.g. Removing a nucleotide sequence's polyadenylation (poly-A tail):\n",
      " |      \n",
      " |      >>> from Bio.Alphabet import IUPAC\n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> my_seq = Seq(\"CGGTACGCTTATGTCACGTAGAAAAAA\", IUPAC.unambiguous_dna)\n",
      " |      >>> my_seq\n",
      " |      Seq('CGGTACGCTTATGTCACGTAGAAAAAA', IUPACUnambiguousDNA())\n",
      " |      >>> my_seq.rstrip(\"A\")\n",
      " |      Seq('CGGTACGCTTATGTCACGTAG', IUPACUnambiguousDNA())\n",
      " |      \n",
      " |      See also the strip and lstrip methods.\n",
      " |  \n",
      " |  split(self, sep=None, maxsplit=-1)\n",
      " |      Split method, like that of a python string.\n",
      " |      \n",
      " |      This behaves like the python string method of the same name.\n",
      " |      \n",
      " |      Return a list of the 'words' in the string (as Seq objects),\n",
      " |      using sep as the delimiter string.  If maxsplit is given, at\n",
      " |      most maxsplit splits are done.  If maxsplit is omitted, all\n",
      " |      splits are made.\n",
      " |      \n",
      " |      Following the python string method, sep will by default be any\n",
      " |      white space (tabs, spaces, newlines) but this is unlikely to\n",
      " |      apply to biological sequences.\n",
      " |      \n",
      " |      e.g.\n",
      " |      \n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> my_rna = Seq(\"GUCAUGGCCAUUGUAAUGGGCCGCUGAAAGGGUGCCCGAUAGUUG\")\n",
      " |      >>> my_aa = my_rna.translate()\n",
      " |      >>> my_aa\n",
      " |      Seq('VMAIVMGR*KGAR*L', HasStopCodon(ExtendedIUPACProtein(), '*'))\n",
      " |      >>> for pep in my_aa.split(\"*\"):\n",
      " |      ...     pep\n",
      " |      Seq('VMAIVMGR', HasStopCodon(ExtendedIUPACProtein(), '*'))\n",
      " |      Seq('KGAR', HasStopCodon(ExtendedIUPACProtein(), '*'))\n",
      " |      Seq('L', HasStopCodon(ExtendedIUPACProtein(), '*'))\n",
      " |      >>> for pep in my_aa.split(\"*\", 1):\n",
      " |      ...     pep\n",
      " |      Seq('VMAIVMGR', HasStopCodon(ExtendedIUPACProtein(), '*'))\n",
      " |      Seq('KGAR*L', HasStopCodon(ExtendedIUPACProtein(), '*'))\n",
      " |      \n",
      " |      See also the rsplit method:\n",
      " |      \n",
      " |      >>> for pep in my_aa.rsplit(\"*\", 1):\n",
      " |      ...     pep\n",
      " |      Seq('VMAIVMGR*KGAR', HasStopCodon(ExtendedIUPACProtein(), '*'))\n",
      " |      Seq('L', HasStopCodon(ExtendedIUPACProtein(), '*'))\n",
      " |  \n",
      " |  startswith(self, prefix, start=0, end=9223372036854775807)\n",
      " |      Return True if the Seq starts with the given prefix, False otherwise.\n",
      " |      \n",
      " |      This behaves like the python string method of the same name.\n",
      " |      \n",
      " |      Return True if the sequence starts with the specified prefix\n",
      " |      (a string or another Seq object), False otherwise.\n",
      " |      With optional start, test sequence beginning at that position.\n",
      " |      With optional end, stop comparing sequence at that position.\n",
      " |      prefix can also be a tuple of strings to try.  e.g.\n",
      " |      \n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> my_rna = Seq(\"GUCAUGGCCAUUGUAAUGGGCCGCUGAAAGGGUGCCCGAUAGUUG\")\n",
      " |      >>> my_rna.startswith(\"GUC\")\n",
      " |      True\n",
      " |      >>> my_rna.startswith(\"AUG\")\n",
      " |      False\n",
      " |      >>> my_rna.startswith(\"AUG\", 3)\n",
      " |      True\n",
      " |      >>> my_rna.startswith((\"UCC\", \"UCA\", \"UCG\"), 1)\n",
      " |      True\n",
      " |  \n",
      " |  strip(self, chars=None)\n",
      " |      Return a new Seq object with leading and trailing ends stripped.\n",
      " |      \n",
      " |      This behaves like the python string method of the same name.\n",
      " |      \n",
      " |      Optional argument chars defines which characters to remove.  If\n",
      " |      omitted or None (default) then as for the python string method,\n",
      " |      this defaults to removing any white space.\n",
      " |      \n",
      " |      e.g. print(my_seq.strip(\"-\"))\n",
      " |      \n",
      " |      See also the lstrip and rstrip methods.\n",
      " |  \n",
      " |  tomutable(self)\n",
      " |      Return the full sequence as a MutableSeq object.\n",
      " |      \n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> from Bio.Alphabet import IUPAC\n",
      " |      >>> my_seq = Seq(\"MKQHKAMIVALIVICITAVVAAL\",\n",
      " |      ...              IUPAC.protein)\n",
      " |      >>> my_seq\n",
      " |      Seq('MKQHKAMIVALIVICITAVVAAL', IUPACProtein())\n",
      " |      >>> my_seq.tomutable()\n",
      " |      MutableSeq('MKQHKAMIVALIVICITAVVAAL', IUPACProtein())\n",
      " |      \n",
      " |      Note that the alphabet is preserved.\n",
      " |  \n",
      " |  transcribe(self)\n",
      " |      Return the RNA sequence from a DNA sequence by creating a new Seq object.\n",
      " |      \n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> from Bio.Alphabet import IUPAC\n",
      " |      >>> coding_dna = Seq(\"ATGGCCATTGTAATGGGCCGCTGAAAGGGTGCCCGATAG\",\n",
      " |      ...                  IUPAC.unambiguous_dna)\n",
      " |      >>> coding_dna\n",
      " |      Seq('ATGGCCATTGTAATGGGCCGCTGAAAGGGTGCCCGATAG', IUPACUnambiguousDNA())\n",
      " |      >>> coding_dna.transcribe()\n",
      " |      Seq('AUGGCCAUUGUAAUGGGCCGCUGAAAGGGUGCCCGAUAG', IUPACUnambiguousRNA())\n",
      " |      \n",
      " |      Trying to transcribe a protein or RNA sequence raises an exception:\n",
      " |      \n",
      " |      >>> my_protein = Seq(\"MAIVMGR\", IUPAC.protein)\n",
      " |      >>> my_protein.transcribe()\n",
      " |      Traceback (most recent call last):\n",
      " |         ...\n",
      " |      ValueError: Proteins cannot be transcribed!\n",
      " |  \n",
      " |  translate(self, table='Standard', stop_symbol='*', to_stop=False, cds=False, gap=None)\n",
      " |      Turn a nucleotide sequence into a protein sequence by creating a new Seq object.\n",
      " |      \n",
      " |      This method will translate DNA or RNA sequences, and those with a\n",
      " |      nucleotide or generic alphabet.  Trying to translate a protein\n",
      " |      sequence raises an exception.\n",
      " |      \n",
      " |      Arguments:\n",
      " |       - table - Which codon table to use?  This can be either a name\n",
      " |         (string), an NCBI identifier (integer), or a CodonTable\n",
      " |         object (useful for non-standard genetic codes).  This\n",
      " |         defaults to the \"Standard\" table.\n",
      " |       - stop_symbol - Single character string, what to use for\n",
      " |         terminators.  This defaults to the asterisk, \"*\".\n",
      " |       - to_stop - Boolean, defaults to False meaning do a full\n",
      " |         translation continuing on past any stop codons (translated as the\n",
      " |         specified stop_symbol).  If True, translation is terminated at\n",
      " |         the first in frame stop codon (and the stop_symbol is not\n",
      " |         appended to the returned protein sequence).\n",
      " |       - cds - Boolean, indicates this is a complete CDS.  If True,\n",
      " |         this checks the sequence starts with a valid alternative start\n",
      " |         codon (which will be translated as methionine, M), that the\n",
      " |         sequence length is a multiple of three, and that there is a\n",
      " |         single in frame stop codon at the end (this will be excluded\n",
      " |         from the protein sequence, regardless of the to_stop option).\n",
      " |         If these tests fail, an exception is raised.\n",
      " |       - gap - Single character string to denote symbol used for gaps.\n",
      " |         It will try to guess the gap character from the alphabet.\n",
      " |      \n",
      " |      e.g. Using the standard table:\n",
      " |      \n",
      " |      >>> coding_dna = Seq(\"GTGGCCATTGTAATGGGCCGCTGAAAGGGTGCCCGATAG\")\n",
      " |      >>> coding_dna.translate()\n",
      " |      Seq('VAIVMGR*KGAR*', HasStopCodon(ExtendedIUPACProtein(), '*'))\n",
      " |      >>> coding_dna.translate(stop_symbol=\"@\")\n",
      " |      Seq('VAIVMGR@KGAR@', HasStopCodon(ExtendedIUPACProtein(), '@'))\n",
      " |      >>> coding_dna.translate(to_stop=True)\n",
      " |      Seq('VAIVMGR', ExtendedIUPACProtein())\n",
      " |      \n",
      " |      Now using NCBI table 2, where TGA is not a stop codon:\n",
      " |      \n",
      " |      >>> coding_dna.translate(table=2)\n",
      " |      Seq('VAIVMGRWKGAR*', HasStopCodon(ExtendedIUPACProtein(), '*'))\n",
      " |      >>> coding_dna.translate(table=2, to_stop=True)\n",
      " |      Seq('VAIVMGRWKGAR', ExtendedIUPACProtein())\n",
      " |      \n",
      " |      In fact, GTG is an alternative start codon under NCBI table 2, meaning\n",
      " |      this sequence could be a complete CDS:\n",
      " |      \n",
      " |      >>> coding_dna.translate(table=2, cds=True)\n",
      " |      Seq('MAIVMGRWKGAR', ExtendedIUPACProtein())\n",
      " |      \n",
      " |      It isn't a valid CDS under NCBI table 1, due to both the start codon\n",
      " |      and also the in frame stop codons:\n",
      " |      \n",
      " |      >>> coding_dna.translate(table=1, cds=True)\n",
      " |      Traceback (most recent call last):\n",
      " |          ...\n",
      " |      Bio.Data.CodonTable.TranslationError: First codon 'GTG' is not a start codon\n",
      " |      \n",
      " |      If the sequence has no in-frame stop codon, then the to_stop argument\n",
      " |      has no effect:\n",
      " |      \n",
      " |      >>> coding_dna2 = Seq(\"TTGGCCATTGTAATGGGCCGC\")\n",
      " |      >>> coding_dna2.translate()\n",
      " |      Seq('LAIVMGR', ExtendedIUPACProtein())\n",
      " |      >>> coding_dna2.translate(to_stop=True)\n",
      " |      Seq('LAIVMGR', ExtendedIUPACProtein())\n",
      " |      \n",
      " |      When translating gapped sequences, the gap character is inferred from\n",
      " |      the alphabet:\n",
      " |      \n",
      " |      >>> from Bio.Alphabet import Gapped\n",
      " |      >>> coding_dna3 = Seq(\"GTG---GCCATT\", Gapped(IUPAC.unambiguous_dna))\n",
      " |      >>> coding_dna3.translate()\n",
      " |      Seq('V-AI', Gapped(ExtendedIUPACProtein(), '-'))\n",
      " |      \n",
      " |      It is possible to pass the gap character when the alphabet is missing:\n",
      " |      \n",
      " |      >>> coding_dna4 = Seq(\"GTG---GCCATT\")\n",
      " |      >>> coding_dna4.translate(gap='-')\n",
      " |      Seq('V-AI', Gapped(ExtendedIUPACProtein(), '-'))\n",
      " |      \n",
      " |      NOTE - Ambiguous codons like \"TAN\" or \"NNN\" could be an amino acid\n",
      " |      or a stop codon.  These are translated as \"X\".  Any invalid codon\n",
      " |      (e.g. \"TA?\" or \"T-A\") will throw a TranslationError.\n",
      " |      \n",
      " |      NOTE - This does NOT behave like the python string's translate\n",
      " |      method.  For that use str(my_seq).translate(...) instead.\n",
      " |  \n",
      " |  ungap(self, gap=None)\n",
      " |      Return a copy of the sequence without the gap character(s).\n",
      " |      \n",
      " |      The gap character can be specified in two ways - either as an explicit\n",
      " |      argument, or via the sequence's alphabet. For example:\n",
      " |      \n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> from Bio.Alphabet import generic_dna\n",
      " |      >>> my_dna = Seq(\"-ATA--TGAAAT-TTGAAAA\", generic_dna)\n",
      " |      >>> my_dna\n",
      " |      Seq('-ATA--TGAAAT-TTGAAAA', DNAAlphabet())\n",
      " |      >>> my_dna.ungap(\"-\")\n",
      " |      Seq('ATATGAAATTTGAAAA', DNAAlphabet())\n",
      " |      \n",
      " |      If the gap character is not given as an argument, it will be taken from\n",
      " |      the sequence's alphabet (if defined). Notice that the returned\n",
      " |      sequence's alphabet is adjusted since it no longer requires a gapped\n",
      " |      alphabet:\n",
      " |      \n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> from Bio.Alphabet import IUPAC, Gapped, HasStopCodon\n",
      " |      >>> my_pro = Seq(\"MVVLE=AD*\", HasStopCodon(Gapped(IUPAC.protein, \"=\")))\n",
      " |      >>> my_pro\n",
      " |      Seq('MVVLE=AD*', HasStopCodon(Gapped(IUPACProtein(), '='), '*'))\n",
      " |      >>> my_pro.ungap()\n",
      " |      Seq('MVVLEAD*', HasStopCodon(IUPACProtein(), '*'))\n",
      " |      \n",
      " |      Or, with a simpler gapped DNA example:\n",
      " |      \n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> from Bio.Alphabet import IUPAC, Gapped\n",
      " |      >>> my_seq = Seq(\"CGGGTAG=AAAAAA\", Gapped(IUPAC.unambiguous_dna, \"=\"))\n",
      " |      >>> my_seq\n",
      " |      Seq('CGGGTAG=AAAAAA', Gapped(IUPACUnambiguousDNA(), '='))\n",
      " |      >>> my_seq.ungap()\n",
      " |      Seq('CGGGTAGAAAAAA', IUPACUnambiguousDNA())\n",
      " |      \n",
      " |      As long as it is consistent with the alphabet, although it is\n",
      " |      redundant, you can still supply the gap character as an argument to\n",
      " |      this method:\n",
      " |      \n",
      " |      >>> my_seq\n",
      " |      Seq('CGGGTAG=AAAAAA', Gapped(IUPACUnambiguousDNA(), '='))\n",
      " |      >>> my_seq.ungap(\"=\")\n",
      " |      Seq('CGGGTAGAAAAAA', IUPACUnambiguousDNA())\n",
      " |      \n",
      " |      However, if the gap character given as the argument disagrees with that\n",
      " |      declared in the alphabet, an exception is raised:\n",
      " |      \n",
      " |      >>> my_seq\n",
      " |      Seq('CGGGTAG=AAAAAA', Gapped(IUPACUnambiguousDNA(), '='))\n",
      " |      >>> my_seq.ungap(\"-\")\n",
      " |      Traceback (most recent call last):\n",
      " |         ...\n",
      " |      ValueError: Gap '-' does not match '=' from alphabet\n",
      " |      \n",
      " |      Finally, if a gap character is not supplied, and the alphabet does not\n",
      " |      define one, an exception is raised:\n",
      " |      \n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> from Bio.Alphabet import generic_dna\n",
      " |      >>> my_dna = Seq(\"ATA--TGAAAT-TTGAAAA\", generic_dna)\n",
      " |      >>> my_dna\n",
      " |      Seq('ATA--TGAAAT-TTGAAAA', DNAAlphabet())\n",
      " |      >>> my_dna.ungap()\n",
      " |      Traceback (most recent call last):\n",
      " |         ...\n",
      " |      ValueError: Gap character not given and not defined in alphabet\n",
      " |  \n",
      " |  upper(self)\n",
      " |      Return an upper case copy of the sequence.\n",
      " |      \n",
      " |      >>> from Bio.Alphabet import HasStopCodon, generic_protein\n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> my_seq = Seq(\"VHLTPeeK*\", HasStopCodon(generic_protein))\n",
      " |      >>> my_seq\n",
      " |      Seq('VHLTPeeK*', HasStopCodon(ProteinAlphabet(), '*'))\n",
      " |      >>> my_seq.lower()\n",
      " |      Seq('vhltpeek*', HasStopCodon(ProteinAlphabet(), '*'))\n",
      " |      >>> my_seq.upper()\n",
      " |      Seq('VHLTPEEK*', HasStopCodon(ProteinAlphabet(), '*'))\n",
      " |      \n",
      " |      This will adjust the alphabet if required. See also the lower method.\n",
      " |  \n",
      " |  ----------------------------------------------------------------------\n",
      " |  Data descriptors defined here:\n",
      " |  \n",
      " |  __dict__\n",
      " |      dictionary for instance variables (if defined)\n",
      " |  \n",
      " |  __weakref__\n",
      " |      list of weak references to the object (if defined)\n",
      "\n"
     ]
    }
   ],
   "source": [
    "from Bio.Seq import Seq\n",
    "help(Seq)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "\"Seq('AGTACACTGGGGT')\""
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "my_seq = Seq(\"AGTACACTGGGGT\")\n",
    "repr(my_seq)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Alphabet()"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "my_seq.alphabet"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Podemos notar que essa sequência não possui um alfabeto. Como podemos saber se a sequência acima é de DNA ou aminoácidos?\n",
    "Para definirmos um alfabeto na criação da sequência precisamos primeiro importar a classe IUPAC da biblioteca Alphabet do Biopython. O IUPAC contém as nomenclaturas e terminologias usadas em química. No nosso caso podemos definir o alfabeto específico para DNA ou aminoácidos.\n",
    "\n",
    "| Código de nucleotídeos IUPAC | Base                |\n",
    "|------------------------------|---------------------|\n",
    "| A                            | Adenine             |\n",
    "| C                            | Cytosine            |\n",
    "| G                            | Guanine             |\n",
    "| T (or U)                     | Thymine (or Uracil) |\n",
    "| R                            | A or G              |\n",
    "| Y                            | C or T              |\n",
    "| S                            | G or C              |\n",
    "| W                            | A or T              |\n",
    "| K                            | G or T              |\n",
    "| M                            | A or C              |\n",
    "| B                            | C or G or T         |\n",
    "| D                            | A or G or T         |\n",
    "| H                            | A or C or T         |\n",
    "| V                            | A or C or G         |\n",
    "| N                            | any base            |\n",
    "| . or -                       | gap                 |\n",
    "\n",
    "\n",
    "| Código de aminoácidos IUPAC | Código de três letras | Aminoácido    |\n",
    "|-----------------------------|-----------------------|---------------|\n",
    "| A                           | Ala                   | Alanine       |\n",
    "| C                           | Cys                   | Cysteine      |\n",
    "| D                           | Asp                   | Aspartic Acid |\n",
    "| E                           | Glu                   | Glutamic Acid |\n",
    "| F                           | Phe                   | Phenylalanine |\n",
    "| G                           | Gly                   | Glycine       |\n",
    "| H                           | His                   | Histidine     |\n",
    "| I                           | Ile                   | Isoleucine    |\n",
    "| K                           | Lys                   | Lysine        |\n",
    "| L                           | Leu                   | Leucine       |\n",
    "| M                           | Met                   | Methionine    |\n",
    "| N                           | Asn                   | Asparagine    |\n",
    "| P                           | Pro                   | Proline       |\n",
    "| Q                           | Gln                   | Glutamine     |\n",
    "| R                           | Arg                   | Arginine      |\n",
    "| S                           | Ser                   | Serine        |\n",
    "| T                           | Thr                   | Threonine     |\n",
    "| V                           | Val                   | Valine        |\n",
    "| W                           | Trp                   | Tryptophan    |\n",
    "| Y                           | Tyr                   | Tyrosine      |"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Aminoácidos: ACDEFGHIKLMNPQRSTVWY\n",
      "DNA: GATC\n",
      "RNA: GAUC\n"
     ]
    }
   ],
   "source": [
    "from Bio.Alphabet import IUPAC\n",
    "print(\"Aminoácidos: \"+IUPAC.protein.letters)\n",
    "print(\"DNA: \"+IUPAC.unambiguous_dna.letters)\n",
    "print(\"RNA: \"+IUPAC.unambiguous_rna.letters)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "IUPACUnambiguousDNA()"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "my_seq = Seq(\"AGTACACTGGGGT\", IUPAC.unambiguous_dna)\n",
    "my_seq.alphabet"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "IUPACProtein()"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "my_prot = Seq(\"AGTACACGGGGT\", IUPAC.protein)\n",
    "my_prot.alphabet"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 1.2 Operações sobre sequências\n",
    "#### 1.2.1 Iterar sobre sequência\n",
    "A classe `Seq` também tem comportamento similar ao da classe string do Python. Sendo assim podemos iterar sobre elas:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0 A\n",
      "1 G\n",
      "2 T\n",
      "3 A\n",
      "4 C\n",
      "5 A\n",
      "6 C\n",
      "7 T\n",
      "8 G\n",
      "9 G\n",
      "10 G\n",
      "11 G\n",
      "12 T\n"
     ]
    }
   ],
   "source": [
    "my_seq = Seq(\"AGTACACTGGGGT\", IUPAC.unambiguous_dna)\n",
    "for index, letter in enumerate(my_seq):\n",
    "    print(\"%i %s\" % (index, letter))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Aqui é um ponto interessante para lembrar que os índices das sequências começam em Zero.\n",
    "\n",
    "#### 1.2.2 Retornar o tamanho da sequência\n",
    "Também podemos utilizar o método len() para verificar o tamanho da sequência."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "13\n"
     ]
    }
   ],
   "source": [
    "print(len(my_seq))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 1.2.3 Retornar determinada posição da sequência\n",
    "O acesso a qualquer caracter da sequência pode ser feito exatamente como em uma string, sendo que o índice -1 retorna o último caracter da sequência."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "A\n",
      "C\n",
      "T\n"
     ]
    }
   ],
   "source": [
    "print(my_seq[0])\n",
    "print(my_seq[4])\n",
    "print(my_seq[-1])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 1.2.4 Contar o número de vezes que uma subsequência aparece\n",
    "A classe sequência também implementa o método `count()`, onde é possível contar o número de ocorrência de uma string na sua sequência. Aqui é importante ressaltar que existem dois métodos para contar, um sem overlap e outro com overlap. Quando estamos contando apenas um caracter os dois métodos terão o mesmo resultado."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Sem overlap: 5\n",
      "Com overlap: 5\n"
     ]
    }
   ],
   "source": [
    "my_seq = Seq(\"AGTACACTGGGGT\", IUPAC.unambiguous_dna)\n",
    "\n",
    "sem_overlap = my_seq.count('G')\n",
    "com_overlap = my_seq.count_overlap('G')\n",
    "\n",
    "print(f\"Sem overlap: {sem_overlap}\")\n",
    "print(f\"Com overlap: {com_overlap}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Porém, quando contamos mais que um caracter podemos ter resultados distintos. No exemplo abaixo, enquanto o método ser overlap vai contando os 'GG' e dando um passo de tamanho 2, no método com overlap ele conta os 'GG' e sempre dá um passo de tamanho 1."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Sem overlap: 2\n",
      "Com overlap: 3\n"
     ]
    }
   ],
   "source": [
    "sem_overlap = my_seq.count('GG')\n",
    "com_overlap = my_seq.count_overlap('GG')\n",
    "\n",
    "print(f\"Sem overlap: {sem_overlap}\")\n",
    "print(f\"Com overlap: {com_overlap}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Uma aplicação simples do uso do método `count()` é para calcular o conteúdo GC de uma sequência. Simplesmente somando o número de Gs e Cs e dividindo pelo tamanho total da sequência, que como vimos antes pode ser obtido pelo método `len()`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "53.84615384615385"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "100 * int(my_seq.count(\"G\") + my_seq.count(\"C\")) / len(my_seq)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Mas o Biopython também tem um método já pronto para calcular o conteúdo GC de uma sequência dentro do pacote SeqUtils."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "53.84615384615385"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from Bio.SeqUtils import GC\n",
    "GC(my_seq)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 1.2.5 Slicing com sequências\n",
    "Assim como a classe string, a classe Seq também aceita slices. Aqui vale lembrar que a primeira posição é 0 e que o intervalo é fechado na posição inicial e aberto na posição final."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Seq('CAC', IUPACUnambiguousDNA())"
      ]
     },
     "execution_count": 17,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "my_seq = Seq(\"AGTACACTGGGGT\", IUPAC.unambiguous_dna)\n",
    "my_seq[4:7]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Outra aplicação interessante do slice é inverter a sequência de maneira simples."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "AGTACACTGGGGT\n",
      "TGGGGTCACATGA\n"
     ]
    }
   ],
   "source": [
    "my_seq_rev = my_seq[::-1]\n",
    "print(my_seq)\n",
    "print(my_seq_rev)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 1.2.6 Concatenar sequências\n",
    "A classe `Seq` também pode te salvar de cometer o erro ao concatenar duas sequências de alfabetos diferentes. Se estivessemos trabalhando com a classe string a operação abaixo seria permitada, porém a classe Seq lança um erro:\n",
    "\n",
    "`TypeError: Incompatible alphabets IUPACProtein() and IUPACUnambiguousDNA()`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "ename": "TypeError",
     "evalue": "Incompatible alphabets IUPACProtein() and IUPACUnambiguousDNA()",
     "output_type": "error",
     "traceback": [
      "\u001b[1;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[1;31mTypeError\u001b[0m                                 Traceback (most recent call last)",
      "\u001b[1;32m<ipython-input-18-493fc9b96e2a>\u001b[0m in \u001b[0;36m<module>\u001b[1;34m\u001b[0m\n\u001b[0;32m      1\u001b[0m \u001b[0mprotein_seq\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mSeq\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;34m\"EVRNAK\"\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mIUPAC\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mprotein\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m      2\u001b[0m \u001b[0mdna_seq\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mSeq\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;34m\"ACGT\"\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mIUPAC\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0munambiguous_dna\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m----> 3\u001b[1;33m \u001b[0mprotein_seq\u001b[0m \u001b[1;33m+\u001b[0m \u001b[0mdna_seq\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m",
      "\u001b[1;32m~\\Anaconda3\\lib\\site-packages\\Bio\\Seq.py\u001b[0m in \u001b[0;36m__add__\u001b[1;34m(self, other)\u001b[0m\n\u001b[0;32m    343\u001b[0m                 raise TypeError(\n\u001b[0;32m    344\u001b[0m                     \"Incompatible alphabets {0!r} and {1!r}\".format(\n\u001b[1;32m--> 345\u001b[1;33m                         self.alphabet, other.alphabet))\n\u001b[0m\u001b[0;32m    346\u001b[0m             \u001b[1;31m# They should be the same sequence type (or one of them is generic)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m    347\u001b[0m             \u001b[0ma\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mAlphabet\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0m_consensus_alphabet\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;33m[\u001b[0m\u001b[0mself\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0malphabet\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mother\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0malphabet\u001b[0m\u001b[1;33m]\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n",
      "\u001b[1;31mTypeError\u001b[0m: Incompatible alphabets IUPACProtein() and IUPACUnambiguousDNA()"
     ]
    }
   ],
   "source": [
    "protein_seq = Seq(\"EVRNAK\", IUPAC.protein)\n",
    "dna_seq = Seq(\"ACGT\", IUPAC.unambiguous_dna)\n",
    "protein_seq + dna_seq"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 1.2.7 Mascarar sequências utilizando os métodos upper() e lower()\n",
    "A tarefa de mudar a sequência para caixa alta ou baixa também pode ser facilmente feita utilizando os métodos upper() e lower(). Esses métodos podem ser úteis para mascarar sequências repetitivas ou de baixa complexidade, que normalmente são representadas por letras minúsculas."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Seq('AGTACACTGGGGT', IUPACUnambiguousDNA())"
      ]
     },
     "execution_count": 17,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "my_seq.upper()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Seq('agtacactggggt', DNAAlphabet())"
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "my_seq.lower()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 1.2.8 Reverso complementar\n",
    "O reverso complementar de uma sequência também pode ser retornado pelo método reverse_complement()."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Sequência original: AGTACACTGGGGT\n",
      "Sequência complementar: TCATGTGACCCCA\n",
      "Sequência reversa complementar: ACCCCAGTGTACT\n"
     ]
    }
   ],
   "source": [
    "print(\"Sequência original: \" + my_seq)\n",
    "\n",
    "print(\"Sequência complementar: \" + my_seq.complement())\n",
    "\n",
    "print(\"Sequência reversa complementar: \" + my_seq.reverse_complement())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "E aqui somos salvos mais uma vez de cometer um erro. Se tentamos gerar o reverso complementar de uma sequência de aminoácidos ele lança o erro:\n",
    "\n",
    "`ValueError: Proteins do not have complements!`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "ename": "ValueError",
     "evalue": "Proteins do not have complements!",
     "output_type": "error",
     "traceback": [
      "\u001b[1;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[1;31mValueError\u001b[0m                                Traceback (most recent call last)",
      "\u001b[1;32m<ipython-input-20-79572275fc13>\u001b[0m in \u001b[0;36m<module>\u001b[1;34m\u001b[0m\n\u001b[0;32m      1\u001b[0m \u001b[0mprotein_seq\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mSeq\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;34m\"EVRNAK\"\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mIUPAC\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mprotein\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m----> 2\u001b[1;33m \u001b[0mprotein_seq\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mreverse_complement\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m",
      "\u001b[1;32m~\\Anaconda3\\lib\\site-packages\\Bio\\Seq.py\u001b[0m in \u001b[0;36mreverse_complement\u001b[1;34m(self)\u001b[0m\n\u001b[0;32m    968\u001b[0m         \"\"\"\n\u001b[0;32m    969\u001b[0m         \u001b[1;31m# Use -1 stride/step to reverse the complement\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m--> 970\u001b[1;33m         \u001b[1;32mreturn\u001b[0m \u001b[0mself\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mcomplement\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m[\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m-\u001b[0m\u001b[1;36m1\u001b[0m\u001b[1;33m]\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m\u001b[0;32m    971\u001b[0m \u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m    972\u001b[0m     \u001b[1;32mdef\u001b[0m \u001b[0mtranscribe\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mself\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n",
      "\u001b[1;32m~\\Anaconda3\\lib\\site-packages\\Bio\\Seq.py\u001b[0m in \u001b[0;36mcomplement\u001b[1;34m(self)\u001b[0m\n\u001b[0;32m    918\u001b[0m         \u001b[0mbase\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mAlphabet\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0m_get_base_alphabet\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mself\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0malphabet\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m    919\u001b[0m         \u001b[1;32mif\u001b[0m \u001b[0misinstance\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mbase\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mAlphabet\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mProteinAlphabet\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m--> 920\u001b[1;33m             \u001b[1;32mraise\u001b[0m \u001b[0mValueError\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;34m\"Proteins do not have complements!\"\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m\u001b[0;32m    921\u001b[0m         \u001b[1;32mif\u001b[0m \u001b[0misinstance\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mbase\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mAlphabet\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mDNAAlphabet\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m    922\u001b[0m             \u001b[0mttable\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0m_dna_complement_table\u001b[0m\u001b[1;33m\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n",
      "\u001b[1;31mValueError\u001b[0m: Proteins do not have complements!"
     ]
    }
   ],
   "source": [
    "protein_seq = Seq(\"EVRNAK\", IUPAC.protein)\n",
    "protein_seq.reverse_complement()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 1.2.9 Transcrição\n",
    "\n",
    "A transcrição é obtida utilizando o método `transcribe()`. O que ele realmente faz é apenas trocar os Ts por Us."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Seq('AUGGCCAUUGUAAUGGGCCGCUGAAAGGGUGCCCGAUAG', IUPACUnambiguousRNA())"
      ]
     },
     "execution_count": 24,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "coding_dna = Seq(\"ATGGCCATTGTAATGGGCCGCTGAAAGGGTGCCCGATAG\", IUPAC.unambiguous_dna)\n",
    "messenger_rna = coding_dna.transcribe()\n",
    "messenger_rna"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 1.2.10 Tradução\n",
    "\n",
    "Utilizando a sequência transcrita acima podemos gerar a tradução dela utilizando o método `translate()`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Seq('MAIVMGR*KGAR*', HasStopCodon(IUPACProtein(), '*'))"
      ]
     },
     "execution_count": 25,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "protein = messenger_rna.translate()\n",
    "protein"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Mas a tradução não precisa ser feita a partir de uma sequência de RNA, ela pode ser obtida também diretamente da sequência de DNA."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Seq('MAIVMGR*KGAR*', HasStopCodon(IUPACProtein(), '*'))"
      ]
     },
     "execution_count": 26,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "coding_dna.translate()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Como visto no exemplo acima, a tradução é feita na sequência toda, então a tradução pode ter multiplos * que é o símbolo que representa o STOP CODON. Se quisermos fazer tradução somente até o primeiro STOP, podemos adicionar o parâmetro `to_stop=True` no método `translate()`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Seq('MAIVMGR', IUPACProtein())"
      ]
     },
     "execution_count": 27,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "coding_dna.translate(to_stop=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Outro parâmetro interessante do método translate() é o table. Com ele você pode especificar qual tabela de codons será utilizada. Esse parâmetro pode ser o número ou nome da tabela. Essa tabelas são baseadas nesta tabela do [NCBI](https://www.ncbi.nlm.nih.gov/Taxonomy/Utils/wprintgc.cgi)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Seq('MAIVMGRWKGAR*', HasStopCodon(IUPACProtein(), '*'))"
      ]
     },
     "execution_count": 25,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "coding_dna.translate(table=\"Vertebrate Mitochondrial\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Tradução: MAIVMGR*KGAR*\n",
      "Tradução até o STOP: MAIVMGR\n",
      "Tradução usando uma tabela diferente: MAIVMGRWKGAR*\n",
      "Tradução até o STOP usando uma tabela de diferente: MAIVMGRWKGAR\n"
     ]
    }
   ],
   "source": [
    "print(\"Tradução: \" + coding_dna.translate())\n",
    "print(\"Tradução até o STOP: \" + coding_dna.translate(to_stop=True))\n",
    "print(\"Tradução usando uma tabela diferente: \" + coding_dna.translate(table=2))\n",
    "print(\"Tradução até o STOP usando uma tabela de diferente: \" + coding_dna.translate(table=2, to_stop=True))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 1.2.11 Tabela de códons\n",
    "\n",
    "As tabelas de códons também podem ser importadas do pacote Data do Biopython. Assim é possível visualiza-la  com um layout bem amigável."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Table 1 Standard, SGC0\n",
      "\n",
      "  |  T      |  C      |  A      |  G      |\n",
      "--+---------+---------+---------+---------+--\n",
      "T | TTT F   | TCT S   | TAT Y   | TGT C   | T\n",
      "T | TTC F   | TCC S   | TAC Y   | TGC C   | C\n",
      "T | TTA L   | TCA S   | TAA Stop| TGA Stop| A\n",
      "T | TTG L(s)| TCG S   | TAG Stop| TGG W   | G\n",
      "--+---------+---------+---------+---------+--\n",
      "C | CTT L   | CCT P   | CAT H   | CGT R   | T\n",
      "C | CTC L   | CCC P   | CAC H   | CGC R   | C\n",
      "C | CTA L   | CCA P   | CAA Q   | CGA R   | A\n",
      "C | CTG L(s)| CCG P   | CAG Q   | CGG R   | G\n",
      "--+---------+---------+---------+---------+--\n",
      "A | ATT I   | ACT T   | AAT N   | AGT S   | T\n",
      "A | ATC I   | ACC T   | AAC N   | AGC S   | C\n",
      "A | ATA I   | ACA T   | AAA K   | AGA R   | A\n",
      "A | ATG M(s)| ACG T   | AAG K   | AGG R   | G\n",
      "--+---------+---------+---------+---------+--\n",
      "G | GTT V   | GCT A   | GAT D   | GGT G   | T\n",
      "G | GTC V   | GCC A   | GAC D   | GGC G   | C\n",
      "G | GTA V   | GCA A   | GAA E   | GGA G   | A\n",
      "G | GTG V   | GCG A   | GAG E   | GGG G   | G\n",
      "--+---------+---------+---------+---------+--\n",
      "Table 2 Vertebrate Mitochondrial, SGC1\n",
      "\n",
      "  |  T      |  C      |  A      |  G      |\n",
      "--+---------+---------+---------+---------+--\n",
      "T | TTT F   | TCT S   | TAT Y   | TGT C   | T\n",
      "T | TTC F   | TCC S   | TAC Y   | TGC C   | C\n",
      "T | TTA L   | TCA S   | TAA Stop| TGA W   | A\n",
      "T | TTG L   | TCG S   | TAG Stop| TGG W   | G\n",
      "--+---------+---------+---------+---------+--\n",
      "C | CTT L   | CCT P   | CAT H   | CGT R   | T\n",
      "C | CTC L   | CCC P   | CAC H   | CGC R   | C\n",
      "C | CTA L   | CCA P   | CAA Q   | CGA R   | A\n",
      "C | CTG L   | CCG P   | CAG Q   | CGG R   | G\n",
      "--+---------+---------+---------+---------+--\n",
      "A | ATT I(s)| ACT T   | AAT N   | AGT S   | T\n",
      "A | ATC I(s)| ACC T   | AAC N   | AGC S   | C\n",
      "A | ATA M(s)| ACA T   | AAA K   | AGA Stop| A\n",
      "A | ATG M(s)| ACG T   | AAG K   | AGG Stop| G\n",
      "--+---------+---------+---------+---------+--\n",
      "G | GTT V   | GCT A   | GAT D   | GGT G   | T\n",
      "G | GTC V   | GCC A   | GAC D   | GGC G   | C\n",
      "G | GTA V   | GCA A   | GAA E   | GGA G   | A\n",
      "G | GTG V(s)| GCG A   | GAG E   | GGG G   | G\n",
      "--+---------+---------+---------+---------+--\n"
     ]
    }
   ],
   "source": [
    "from Bio.Data import CodonTable\n",
    "standard_table = CodonTable.unambiguous_dna_by_name[\"Standard\"]\n",
    "mitochondrial_table = CodonTable.unambiguous_dna_by_name[\"Vertebrate Mitochondrial\"]\n",
    "print(standard_table)\n",
    "print(mitochondrial_table)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Também é possível acessar outros atríbutos da tabela que podem ser úteis em outras aplicações."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "['TAA', 'TAG', 'AGA', 'AGG']"
      ]
     },
     "execution_count": 30,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "mitochondrial_table.stop_codons"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "['ATT', 'ATC', 'ATA', 'ATG', 'GTG']"
      ]
     },
     "execution_count": 31,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "mitochondrial_table.start_codons"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Exercício 1\n",
    "Dada uma lista de sequências de DNA, calcule o conteúdo GC de cada uma e se ela tiver um conteúdo GC >= 45, traduza a sequência até o primeiro STOP Códon e imprima na tela."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "seqs = [Seq(\"ATGGCGTCGATCGATTTTACAT\"), Seq(\"CGATTCCGCGCGACGCGCAGCTGC\"), Seq(\"CCATACAAACATCTCAG\"), Seq(\"CGCAGCTCGACTCGCTCTTC\")]\n",
    "\n",
    "for seq in seqs:\n",
    "    gc = GC(seq)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2. SeqRecord\n",
    "\n",
    "Até aqui vimos como trabalhar com um tipo bem simples de dado que é uma sequência e a classe `Seq` é muito útil para isso, mas na prática somente a sequência é pouco informativa. Normalmente essa sequência precisa de id, um nome, uma descrição e outros dados. Para essa aplicações temos a classe `SeqRecord` que possui esses atributos e muitos mais.\n",
    "Então por que aprender a classe Seq? Porque um dos atributos da classe `SeqRecord` é justamente um objeto da classe `Seq`. No help abaixo vocês podem ver que os principais atributos são id e seq e esses são os atributos básicos que um arquivo no formato FASTA tem.\n",
    "\n",
    "Até aqui vimos como trabalhar com um tipo bem simples de dado que a sequência e classe Seq é muito útil para isso, mas na prática somente a sequência é pouco informativa. Normalmente essa sequência precisa de um id, nome, descrição, ... . Para essas aplicações temos a classe SeqRecord que possuí esses atributos e muitos mais.\n",
    "Então por que aprender a classe Seq? Por que um dos atributos da classe SeqRecord é justamente um objeto da classe Seq. No help abaixo vocês podem ver que os principais atributos são id e seq e esses são os atributos básico que um arquivo no formato FASTA tem."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Help on class SeqRecord in module Bio.SeqRecord:\n",
      "\n",
      "class SeqRecord(builtins.object)\n",
      " |  SeqRecord(seq, id='<unknown id>', name='<unknown name>', description='<unknown description>', dbxrefs=None, features=None, annotations=None, letter_annotations=None)\n",
      " |  \n",
      " |  A SeqRecord object holds a sequence and information about it.\n",
      " |  \n",
      " |  Main attributes:\n",
      " |   - id          - Identifier such as a locus tag (string)\n",
      " |   - seq         - The sequence itself (Seq object or similar)\n",
      " |  \n",
      " |  Additional attributes:\n",
      " |   - name        - Sequence name, e.g. gene name (string)\n",
      " |   - description - Additional text (string)\n",
      " |   - dbxrefs     - List of database cross references (list of strings)\n",
      " |   - features    - Any (sub)features defined (list of SeqFeature objects)\n",
      " |   - annotations - Further information about the whole sequence (dictionary).\n",
      " |     Most entries are strings, or lists of strings.\n",
      " |   - letter_annotations - Per letter/symbol annotation (restricted\n",
      " |     dictionary). This holds Python sequences (lists, strings\n",
      " |     or tuples) whose length matches that of the sequence.\n",
      " |     A typical use would be to hold a list of integers\n",
      " |     representing sequencing quality scores, or a string\n",
      " |     representing the secondary structure.\n",
      " |  \n",
      " |  You will typically use Bio.SeqIO to read in sequences from files as\n",
      " |  SeqRecord objects.  However, you may want to create your own SeqRecord\n",
      " |  objects directly (see the __init__ method for further details):\n",
      " |  \n",
      " |  >>> from Bio.Seq import Seq\n",
      " |  >>> from Bio.SeqRecord import SeqRecord\n",
      " |  >>> from Bio.Alphabet import IUPAC\n",
      " |  >>> record = SeqRecord(Seq(\"MKQHKAMIVALIVICITAVVAALVTRKDLCEVHIRTGQTEVAVF\",\n",
      " |  ...                         IUPAC.protein),\n",
      " |  ...                    id=\"YP_025292.1\", name=\"HokC\",\n",
      " |  ...                    description=\"toxic membrane protein\")\n",
      " |  >>> print(record)\n",
      " |  ID: YP_025292.1\n",
      " |  Name: HokC\n",
      " |  Description: toxic membrane protein\n",
      " |  Number of features: 0\n",
      " |  Seq('MKQHKAMIVALIVICITAVVAALVTRKDLCEVHIRTGQTEVAVF', IUPACProtein())\n",
      " |  \n",
      " |  If you want to save SeqRecord objects to a sequence file, use Bio.SeqIO\n",
      " |  for this.  For the special case where you want the SeqRecord turned into\n",
      " |  a string in a particular file format there is a format method which uses\n",
      " |  Bio.SeqIO internally:\n",
      " |  \n",
      " |  >>> print(record.format(\"fasta\"))\n",
      " |  >YP_025292.1 toxic membrane protein\n",
      " |  MKQHKAMIVALIVICITAVVAALVTRKDLCEVHIRTGQTEVAVF\n",
      " |  <BLANKLINE>\n",
      " |  \n",
      " |  You can also do things like slicing a SeqRecord, checking its length, etc\n",
      " |  \n",
      " |  >>> len(record)\n",
      " |  44\n",
      " |  >>> edited = record[:10] + record[11:]\n",
      " |  >>> print(edited.seq)\n",
      " |  MKQHKAMIVAIVICITAVVAALVTRKDLCEVHIRTGQTEVAVF\n",
      " |  >>> print(record.seq)\n",
      " |  MKQHKAMIVALIVICITAVVAALVTRKDLCEVHIRTGQTEVAVF\n",
      " |  \n",
      " |  Methods defined here:\n",
      " |  \n",
      " |  __add__(self, other)\n",
      " |      Add another sequence or string to this sequence.\n",
      " |      \n",
      " |      The other sequence can be a SeqRecord object, a Seq object (or\n",
      " |      similar, e.g. a MutableSeq) or a plain Python string. If you add\n",
      " |      a plain string or a Seq (like) object, the new SeqRecord will simply\n",
      " |      have this appended to the existing data. However, any per letter\n",
      " |      annotation will be lost:\n",
      " |      \n",
      " |      >>> from Bio import SeqIO\n",
      " |      >>> record = SeqIO.read(\"Quality/solexa_faked.fastq\", \"fastq-solexa\")\n",
      " |      >>> print(\"%s %s\" % (record.id, record.seq))\n",
      " |      slxa_0001_1_0001_01 ACGTACGTACGTACGTACGTACGTACGTACGTACGTACGTNNNNNN\n",
      " |      >>> print(list(record.letter_annotations))\n",
      " |      ['solexa_quality']\n",
      " |      \n",
      " |      >>> new = record + \"ACT\"\n",
      " |      >>> print(\"%s %s\" % (new.id, new.seq))\n",
      " |      slxa_0001_1_0001_01 ACGTACGTACGTACGTACGTACGTACGTACGTACGTACGTNNNNNNACT\n",
      " |      >>> print(list(new.letter_annotations))\n",
      " |      []\n",
      " |      \n",
      " |      The new record will attempt to combine the annotation, but for any\n",
      " |      ambiguities (e.g. different names) it defaults to omitting that\n",
      " |      annotation.\n",
      " |      \n",
      " |      >>> from Bio import SeqIO\n",
      " |      >>> with open(\"GenBank/pBAD30.gb\") as handle:\n",
      " |      ...     plasmid = SeqIO.read(handle, \"gb\")\n",
      " |      >>> print(\"%s %i\" % (plasmid.id, len(plasmid)))\n",
      " |      pBAD30 4923\n",
      " |      \n",
      " |      Now let's cut the plasmid into two pieces, and join them back up the\n",
      " |      other way round (i.e. shift the starting point on this plasmid, have\n",
      " |      a look at the annotated features in the original file to see why this\n",
      " |      particular split point might make sense):\n",
      " |      \n",
      " |      >>> left = plasmid[:3765]\n",
      " |      >>> right = plasmid[3765:]\n",
      " |      >>> new = right + left\n",
      " |      >>> print(\"%s %i\" % (new.id, len(new)))\n",
      " |      pBAD30 4923\n",
      " |      >>> str(new.seq) == str(right.seq + left.seq)\n",
      " |      True\n",
      " |      >>> len(new.features) == len(left.features) + len(right.features)\n",
      " |      True\n",
      " |      \n",
      " |      When we add the left and right SeqRecord objects, their annotation\n",
      " |      is all consistent, so it is all conserved in the new SeqRecord:\n",
      " |      \n",
      " |      >>> new.id == left.id == right.id == plasmid.id\n",
      " |      True\n",
      " |      >>> new.name == left.name == right.name == plasmid.name\n",
      " |      True\n",
      " |      >>> new.description == plasmid.description\n",
      " |      True\n",
      " |      >>> new.annotations == left.annotations == right.annotations\n",
      " |      True\n",
      " |      >>> new.letter_annotations == plasmid.letter_annotations\n",
      " |      True\n",
      " |      >>> new.dbxrefs == left.dbxrefs == right.dbxrefs\n",
      " |      True\n",
      " |      \n",
      " |      However, we should point out that when we sliced the SeqRecord,\n",
      " |      any annotations dictionary or dbxrefs list entries were lost.\n",
      " |      You can explicitly copy them like this:\n",
      " |      \n",
      " |      >>> new.annotations = plasmid.annotations.copy()\n",
      " |      >>> new.dbxrefs = plasmid.dbxrefs[:]\n",
      " |  \n",
      " |  __bool__(self)\n",
      " |      Boolean value of an instance of this class (True).\n",
      " |      \n",
      " |      This behaviour is for backwards compatibility, since until the\n",
      " |      __len__ method was added, a SeqRecord always evaluated as True.\n",
      " |      \n",
      " |      Note that in comparison, a Seq object will evaluate to False if it\n",
      " |      has a zero length sequence.\n",
      " |      \n",
      " |      WARNING: The SeqRecord may in future evaluate to False when its\n",
      " |      sequence is of zero length (in order to better match the Seq\n",
      " |      object behaviour)!\n",
      " |  \n",
      " |  __contains__(self, char)\n",
      " |      Implement the 'in' keyword, searches the sequence.\n",
      " |      \n",
      " |      e.g.\n",
      " |      \n",
      " |      >>> from Bio import SeqIO\n",
      " |      >>> record = SeqIO.read(\"Fasta/sweetpea.nu\", \"fasta\")\n",
      " |      >>> \"GAATTC\" in record\n",
      " |      False\n",
      " |      >>> \"AAA\" in record\n",
      " |      True\n",
      " |      \n",
      " |      This essentially acts as a proxy for using \"in\" on the sequence:\n",
      " |      \n",
      " |      >>> \"GAATTC\" in record.seq\n",
      " |      False\n",
      " |      >>> \"AAA\" in record.seq\n",
      " |      True\n",
      " |      \n",
      " |      Note that you can also use Seq objects as the query,\n",
      " |      \n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> from Bio.Alphabet import generic_dna\n",
      " |      >>> Seq(\"AAA\") in record\n",
      " |      True\n",
      " |      >>> Seq(\"AAA\", generic_dna) in record\n",
      " |      True\n",
      " |      \n",
      " |      See also the Seq object's __contains__ method.\n",
      " |  \n",
      " |  __eq__(self, other)\n",
      " |      Return self==value.\n",
      " |  \n",
      " |  __format__(self, format_spec)\n",
      " |      Return the record as a string in the specified file format.\n",
      " |      \n",
      " |      This method supports the python format() function added in\n",
      " |      Python 2.6/3.0.  The format_spec should be a lower case string\n",
      " |      supported by Bio.SeqIO as an output file format. See also the\n",
      " |      SeqRecord's format() method.\n",
      " |      \n",
      " |      Under Python 3 please note that for binary formats a bytes\n",
      " |      string is returned, otherwise a (unicode) string is returned.\n",
      " |  \n",
      " |  __ge__(self, other)\n",
      " |      Return self>=value.\n",
      " |  \n",
      " |  __getitem__(self, index)\n",
      " |      Return a sub-sequence or an individual letter.\n",
      " |      \n",
      " |      Slicing, e.g. my_record[5:10], returns a new SeqRecord for\n",
      " |      that sub-sequence with some annotation preserved as follows:\n",
      " |      \n",
      " |      * The name, id and description are kept as-is.\n",
      " |      * Any per-letter-annotations are sliced to match the requested\n",
      " |        sub-sequence.\n",
      " |      * Unless a stride is used, all those features which fall fully\n",
      " |        within the subsequence are included (with their locations\n",
      " |        adjusted accordingly). If you want to preserve any truncated\n",
      " |        features (e.g. GenBank/EMBL source features), you must\n",
      " |        explicitly add them to the new SeqRecord yourself.\n",
      " |      * The annotations dictionary and the dbxrefs list are not used\n",
      " |        for the new SeqRecord, as in general they may not apply to the\n",
      " |        subsequence. If you want to preserve them, you must explicitly\n",
      " |        copy them to the new SeqRecord yourself.\n",
      " |      \n",
      " |      Using an integer index, e.g. my_record[5] is shorthand for\n",
      " |      extracting that letter from the sequence, my_record.seq[5].\n",
      " |      \n",
      " |      For example, consider this short protein and its secondary\n",
      " |      structure as encoded by the PDB (e.g. H for alpha helices),\n",
      " |      plus a simple feature for its histidine self phosphorylation\n",
      " |      site:\n",
      " |      \n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> from Bio.SeqRecord import SeqRecord\n",
      " |      >>> from Bio.SeqFeature import SeqFeature, FeatureLocation\n",
      " |      >>> from Bio.Alphabet import IUPAC\n",
      " |      >>> rec = SeqRecord(Seq(\"MAAGVKQLADDRTLLMAGVSHDLRTPLTRIRLAT\"\n",
      " |      ...                     \"EMMSEQDGYLAESINKDIEECNAIIEQFIDYLR\",\n",
      " |      ...                     IUPAC.protein),\n",
      " |      ...                 id=\"1JOY\", name=\"EnvZ\",\n",
      " |      ...                 description=\"Homodimeric domain of EnvZ from E. coli\")\n",
      " |      >>> rec.letter_annotations[\"secondary_structure\"] = \"  S  SSSSSSHHHHHTTTHHHHHHHHHHHHHHHHHHHHHHTHHHHHHHHHHHHHHHHHHHHHTT  \"\n",
      " |      >>> rec.features.append(SeqFeature(FeatureLocation(20, 21),\n",
      " |      ...                     type = \"Site\"))\n",
      " |      \n",
      " |      Now let's have a quick look at the full record,\n",
      " |      \n",
      " |      >>> print(rec)\n",
      " |      ID: 1JOY\n",
      " |      Name: EnvZ\n",
      " |      Description: Homodimeric domain of EnvZ from E. coli\n",
      " |      Number of features: 1\n",
      " |      Per letter annotation for: secondary_structure\n",
      " |      Seq('MAAGVKQLADDRTLLMAGVSHDLRTPLTRIRLATEMMSEQDGYLAESINKDIEE...YLR', IUPACProtein())\n",
      " |      >>> rec.letter_annotations[\"secondary_structure\"]\n",
      " |      '  S  SSSSSSHHHHHTTTHHHHHHHHHHHHHHHHHHHHHHTHHHHHHHHHHHHHHHHHHHHHTT  '\n",
      " |      >>> print(rec.features[0].location)\n",
      " |      [20:21]\n",
      " |      \n",
      " |      Now let's take a sub sequence, here chosen as the first (fractured)\n",
      " |      alpha helix which includes the histidine phosphorylation site:\n",
      " |      \n",
      " |      >>> sub = rec[11:41]\n",
      " |      >>> print(sub)\n",
      " |      ID: 1JOY\n",
      " |      Name: EnvZ\n",
      " |      Description: Homodimeric domain of EnvZ from E. coli\n",
      " |      Number of features: 1\n",
      " |      Per letter annotation for: secondary_structure\n",
      " |      Seq('RTLLMAGVSHDLRTPLTRIRLATEMMSEQD', IUPACProtein())\n",
      " |      >>> sub.letter_annotations[\"secondary_structure\"]\n",
      " |      'HHHHHTTTHHHHHHHHHHHHHHHHHHHHHH'\n",
      " |      >>> print(sub.features[0].location)\n",
      " |      [9:10]\n",
      " |      \n",
      " |      You can also of course omit the start or end values, for\n",
      " |      example to get the first ten letters only:\n",
      " |      \n",
      " |      >>> print(rec[:10])\n",
      " |      ID: 1JOY\n",
      " |      Name: EnvZ\n",
      " |      Description: Homodimeric domain of EnvZ from E. coli\n",
      " |      Number of features: 0\n",
      " |      Per letter annotation for: secondary_structure\n",
      " |      Seq('MAAGVKQLAD', IUPACProtein())\n",
      " |      \n",
      " |      Or for the last ten letters:\n",
      " |      \n",
      " |      >>> print(rec[-10:])\n",
      " |      ID: 1JOY\n",
      " |      Name: EnvZ\n",
      " |      Description: Homodimeric domain of EnvZ from E. coli\n",
      " |      Number of features: 0\n",
      " |      Per letter annotation for: secondary_structure\n",
      " |      Seq('IIEQFIDYLR', IUPACProtein())\n",
      " |      \n",
      " |      If you omit both, then you get a copy of the original record (although\n",
      " |      lacking the annotations and dbxrefs):\n",
      " |      \n",
      " |      >>> print(rec[:])\n",
      " |      ID: 1JOY\n",
      " |      Name: EnvZ\n",
      " |      Description: Homodimeric domain of EnvZ from E. coli\n",
      " |      Number of features: 1\n",
      " |      Per letter annotation for: secondary_structure\n",
      " |      Seq('MAAGVKQLADDRTLLMAGVSHDLRTPLTRIRLATEMMSEQDGYLAESINKDIEE...YLR', IUPACProtein())\n",
      " |      \n",
      " |      Finally, indexing with a simple integer is shorthand for pulling out\n",
      " |      that letter from the sequence directly:\n",
      " |      \n",
      " |      >>> rec[5]\n",
      " |      'K'\n",
      " |      >>> rec.seq[5]\n",
      " |      'K'\n",
      " |  \n",
      " |  __gt__(self, other)\n",
      " |      Return self>value.\n",
      " |  \n",
      " |  __init__(self, seq, id='<unknown id>', name='<unknown name>', description='<unknown description>', dbxrefs=None, features=None, annotations=None, letter_annotations=None)\n",
      " |      Create a SeqRecord.\n",
      " |      \n",
      " |      Arguments:\n",
      " |       - seq         - Sequence, required (Seq, MutableSeq or UnknownSeq)\n",
      " |       - id          - Sequence identifier, recommended (string)\n",
      " |       - name        - Sequence name, optional (string)\n",
      " |       - description - Sequence description, optional (string)\n",
      " |       - dbxrefs     - Database cross references, optional (list of strings)\n",
      " |       - features    - Any (sub)features, optional (list of SeqFeature objects)\n",
      " |       - annotations - Dictionary of annotations for the whole sequence\n",
      " |       - letter_annotations - Dictionary of per-letter-annotations, values\n",
      " |         should be strings, list or tuples of the same length as the full\n",
      " |         sequence.\n",
      " |      \n",
      " |      You will typically use Bio.SeqIO to read in sequences from files as\n",
      " |      SeqRecord objects.  However, you may want to create your own SeqRecord\n",
      " |      objects directly.\n",
      " |      \n",
      " |      Note that while an id is optional, we strongly recommend you supply a\n",
      " |      unique id string for each record.  This is especially important\n",
      " |      if you wish to write your sequences to a file.\n",
      " |      \n",
      " |      If you don't have the actual sequence, but you do know its length,\n",
      " |      then using the UnknownSeq object from Bio.Seq is appropriate.\n",
      " |      \n",
      " |      You can create a 'blank' SeqRecord object, and then populate the\n",
      " |      attributes later.\n",
      " |  \n",
      " |  __iter__(self)\n",
      " |      Iterate over the letters in the sequence.\n",
      " |      \n",
      " |      For example, using Bio.SeqIO to read in a protein FASTA file:\n",
      " |      \n",
      " |      >>> from Bio import SeqIO\n",
      " |      >>> record = SeqIO.read(\"Fasta/loveliesbleeding.pro\", \"fasta\")\n",
      " |      >>> for amino in record:\n",
      " |      ...     print(amino)\n",
      " |      ...     if amino == \"L\": break\n",
      " |      X\n",
      " |      A\n",
      " |      G\n",
      " |      L\n",
      " |      >>> print(record.seq[3])\n",
      " |      L\n",
      " |      \n",
      " |      This is just a shortcut for iterating over the sequence directly:\n",
      " |      \n",
      " |      >>> for amino in record.seq:\n",
      " |      ...     print(amino)\n",
      " |      ...     if amino == \"L\": break\n",
      " |      X\n",
      " |      A\n",
      " |      G\n",
      " |      L\n",
      " |      >>> print(record.seq[3])\n",
      " |      L\n",
      " |      \n",
      " |      Note that this does not facilitate iteration together with any\n",
      " |      per-letter-annotation.  However, you can achieve that using the\n",
      " |      python zip function on the record (or its sequence) and the relevant\n",
      " |      per-letter-annotation:\n",
      " |      \n",
      " |      >>> from Bio import SeqIO\n",
      " |      >>> rec = SeqIO.read(\"Quality/solexa_faked.fastq\", \"fastq-solexa\")\n",
      " |      >>> print(\"%s %s\" % (rec.id, rec.seq))\n",
      " |      slxa_0001_1_0001_01 ACGTACGTACGTACGTACGTACGTACGTACGTACGTACGTNNNNNN\n",
      " |      >>> print(list(rec.letter_annotations))\n",
      " |      ['solexa_quality']\n",
      " |      >>> for nuc, qual in zip(rec, rec.letter_annotations[\"solexa_quality\"]):\n",
      " |      ...     if qual > 35:\n",
      " |      ...         print(\"%s %i\" % (nuc, qual))\n",
      " |      A 40\n",
      " |      C 39\n",
      " |      G 38\n",
      " |      T 37\n",
      " |      A 36\n",
      " |      \n",
      " |      You may agree that using zip(rec.seq, ...) is more explicit than using\n",
      " |      zip(rec, ...) as shown above.\n",
      " |  \n",
      " |  __le___(self, other)\n",
      " |  \n",
      " |  __len__(self)\n",
      " |      Return the length of the sequence.\n",
      " |      \n",
      " |      For example, using Bio.SeqIO to read in a FASTA nucleotide file:\n",
      " |      \n",
      " |      >>> from Bio import SeqIO\n",
      " |      >>> record = SeqIO.read(\"Fasta/sweetpea.nu\", \"fasta\")\n",
      " |      >>> len(record)\n",
      " |      309\n",
      " |      >>> len(record.seq)\n",
      " |      309\n",
      " |  \n",
      " |  __lt__(self, other)\n",
      " |      Return self<value.\n",
      " |  \n",
      " |  __ne__(self, other)\n",
      " |      Return self!=value.\n",
      " |  \n",
      " |  __nonzero__ = __bool__(self)\n",
      " |  \n",
      " |  __radd__(self, other)\n",
      " |      Add another sequence or string to this sequence (from the left).\n",
      " |      \n",
      " |      This method handles adding a Seq object (or similar, e.g. MutableSeq)\n",
      " |      or a plain Python string (on the left) to a SeqRecord (on the right).\n",
      " |      See the __add__ method for more details, but for example:\n",
      " |      \n",
      " |      >>> from Bio import SeqIO\n",
      " |      >>> record = SeqIO.read(\"Quality/solexa_faked.fastq\", \"fastq-solexa\")\n",
      " |      >>> print(\"%s %s\" % (record.id, record.seq))\n",
      " |      slxa_0001_1_0001_01 ACGTACGTACGTACGTACGTACGTACGTACGTACGTACGTNNNNNN\n",
      " |      >>> print(list(record.letter_annotations))\n",
      " |      ['solexa_quality']\n",
      " |      \n",
      " |      >>> new = \"ACT\" + record\n",
      " |      >>> print(\"%s %s\" % (new.id, new.seq))\n",
      " |      slxa_0001_1_0001_01 ACTACGTACGTACGTACGTACGTACGTACGTACGTACGTACGTNNNNNN\n",
      " |      >>> print(list(new.letter_annotations))\n",
      " |      []\n",
      " |  \n",
      " |  __repr__(self)\n",
      " |      Return a concise summary of the record for debugging (string).\n",
      " |      \n",
      " |      The python built in function repr works by calling the object's ___repr__\n",
      " |      method.  e.g.\n",
      " |      \n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> from Bio.SeqRecord import SeqRecord\n",
      " |      >>> from Bio.Alphabet import generic_protein\n",
      " |      >>> rec = SeqRecord(Seq(\"MASRGVNKVILVGNLGQDPEVRYMPNGGAVANITLATSESWRDKAT\"\n",
      " |      ...                    +\"GEMKEQTEWHRVVLFGKLAEVASEYLRKGSQVYIEGQLRTRKWTDQ\"\n",
      " |      ...                    +\"SGQDRYTTEVVVNVGGTMQMLGGRQGGGAPAGGNIGGGQPQGGWGQ\"\n",
      " |      ...                    +\"PQQPQGGNQFSGGAQSRPQQSAPAAPSNEPPMDFDDDIPF\",\n",
      " |      ...                    generic_protein),\n",
      " |      ...                 id=\"NP_418483.1\", name=\"b4059\",\n",
      " |      ...                 description=\"ssDNA-binding protein\",\n",
      " |      ...                 dbxrefs=[\"ASAP:13298\", \"GI:16131885\", \"GeneID:948570\"])\n",
      " |      >>> print(repr(rec))\n",
      " |      SeqRecord(seq=Seq('MASRGVNKVILVGNLGQDPEVRYMPNGGAVANITLATSESWRDKATGEMKEQTE...IPF', ProteinAlphabet()), id='NP_418483.1', name='b4059', description='ssDNA-binding protein', dbxrefs=['ASAP:13298', 'GI:16131885', 'GeneID:948570'])\n",
      " |      \n",
      " |      At the python prompt you can also use this shorthand:\n",
      " |      \n",
      " |      >>> rec\n",
      " |      SeqRecord(seq=Seq('MASRGVNKVILVGNLGQDPEVRYMPNGGAVANITLATSESWRDKATGEMKEQTE...IPF', ProteinAlphabet()), id='NP_418483.1', name='b4059', description='ssDNA-binding protein', dbxrefs=['ASAP:13298', 'GI:16131885', 'GeneID:948570'])\n",
      " |      \n",
      " |      Note that long sequences are shown truncated. Also note that any\n",
      " |      annotations, letter_annotations and features are not shown (as they\n",
      " |      would lead to a very long string).\n",
      " |  \n",
      " |  __str__(self)\n",
      " |      Return a human readable summary of the record and its annotation (string).\n",
      " |      \n",
      " |      The python built in function str works by calling the object's ___str__\n",
      " |      method.  e.g.\n",
      " |      \n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> from Bio.SeqRecord import SeqRecord\n",
      " |      >>> from Bio.Alphabet import IUPAC\n",
      " |      >>> record = SeqRecord(Seq(\"MKQHKAMIVALIVICITAVVAALVTRKDLCEVHIRTGQTEVAVF\",\n",
      " |      ...                         IUPAC.protein),\n",
      " |      ...                    id=\"YP_025292.1\", name=\"HokC\",\n",
      " |      ...                    description=\"toxic membrane protein, small\")\n",
      " |      >>> print(str(record))\n",
      " |      ID: YP_025292.1\n",
      " |      Name: HokC\n",
      " |      Description: toxic membrane protein, small\n",
      " |      Number of features: 0\n",
      " |      Seq('MKQHKAMIVALIVICITAVVAALVTRKDLCEVHIRTGQTEVAVF', IUPACProtein())\n",
      " |      \n",
      " |      In this example you don't actually need to call str explicity, as the\n",
      " |      print command does this automatically:\n",
      " |      \n",
      " |      >>> print(record)\n",
      " |      ID: YP_025292.1\n",
      " |      Name: HokC\n",
      " |      Description: toxic membrane protein, small\n",
      " |      Number of features: 0\n",
      " |      Seq('MKQHKAMIVALIVICITAVVAALVTRKDLCEVHIRTGQTEVAVF', IUPACProtein())\n",
      " |      \n",
      " |      Note that long sequences are shown truncated.\n",
      " |  \n",
      " |  format(self, format)\n",
      " |      Return the record as a string in the specified file format.\n",
      " |      \n",
      " |      The format should be a lower case string supported as an output\n",
      " |      format by Bio.SeqIO, which is used to turn the SeqRecord into a\n",
      " |      string.  e.g.\n",
      " |      \n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> from Bio.SeqRecord import SeqRecord\n",
      " |      >>> from Bio.Alphabet import IUPAC\n",
      " |      >>> record = SeqRecord(Seq(\"MKQHKAMIVALIVICITAVVAALVTRKDLCEVHIRTGQTEVAVF\",\n",
      " |      ...                         IUPAC.protein),\n",
      " |      ...                    id=\"YP_025292.1\", name=\"HokC\",\n",
      " |      ...                    description=\"toxic membrane protein\")\n",
      " |      >>> record.format(\"fasta\")\n",
      " |      '>YP_025292.1 toxic membrane protein\\nMKQHKAMIVALIVICITAVVAALVTRKDLCEVHIRTGQTEVAVF\\n'\n",
      " |      >>> print(record.format(\"fasta\"))\n",
      " |      >YP_025292.1 toxic membrane protein\n",
      " |      MKQHKAMIVALIVICITAVVAALVTRKDLCEVHIRTGQTEVAVF\n",
      " |      <BLANKLINE>\n",
      " |      \n",
      " |      The python print command automatically appends a new line, meaning\n",
      " |      in this example a blank line is shown.  If you look at the string\n",
      " |      representation you can see there is a trailing new line (shown as\n",
      " |      slash n) which is important when writing to a file or if\n",
      " |      concatenating multiple sequence strings together.\n",
      " |      \n",
      " |      Note that this method will NOT work on every possible file format\n",
      " |      supported by Bio.SeqIO (e.g. some are for multiple sequences only).\n",
      " |  \n",
      " |  lower(self)\n",
      " |      Return a copy of the record with a lower case sequence.\n",
      " |      \n",
      " |      All the annotation is preserved unchanged. e.g.\n",
      " |      \n",
      " |      >>> from Bio import SeqIO\n",
      " |      >>> record = SeqIO.read(\"Fasta/aster.pro\", \"fasta\")\n",
      " |      >>> print(record.format(\"fasta\"))\n",
      " |      >gi|3298468|dbj|BAA31520.1| SAMIPF\n",
      " |      GGHVNPAVTFGAFVGGNITLLRGIVYIIAQLLGSTVACLLLKFVTNDMAVGVFSLSAGVG\n",
      " |      VTNALVFEIVMTFGLVYTVYATAIDPKKGSLGTIAPIAIGFIVGANI\n",
      " |      <BLANKLINE>\n",
      " |      >>> print(record.lower().format(\"fasta\"))\n",
      " |      >gi|3298468|dbj|BAA31520.1| SAMIPF\n",
      " |      gghvnpavtfgafvggnitllrgivyiiaqllgstvaclllkfvtndmavgvfslsagvg\n",
      " |      vtnalvfeivmtfglvytvyataidpkkgslgtiapiaigfivgani\n",
      " |      <BLANKLINE>\n",
      " |      \n",
      " |      To take a more annotation rich example,\n",
      " |      \n",
      " |      >>> from Bio import SeqIO\n",
      " |      >>> old = SeqIO.read(\"EMBL/TRBG361.embl\", \"embl\")\n",
      " |      >>> len(old.features)\n",
      " |      3\n",
      " |      >>> new = old.lower()\n",
      " |      >>> len(old.features) == len(new.features)\n",
      " |      True\n",
      " |      >>> old.annotations[\"organism\"] == new.annotations[\"organism\"]\n",
      " |      True\n",
      " |      >>> old.dbxrefs == new.dbxrefs\n",
      " |      True\n",
      " |  \n",
      " |  reverse_complement(self, id=False, name=False, description=False, features=True, annotations=False, letter_annotations=True, dbxrefs=False)\n",
      " |      Return new SeqRecord with reverse complement sequence.\n",
      " |      \n",
      " |      By default the new record does NOT preserve the sequence identifier,\n",
      " |      name, description, general annotation or database cross-references -\n",
      " |      these are unlikely to apply to the reversed sequence.\n",
      " |      \n",
      " |      You can specify the returned record's id, name and description as\n",
      " |      strings, or True to keep that of the parent, or False for a default.\n",
      " |      \n",
      " |      You can specify the returned record's features with a list of\n",
      " |      SeqFeature objects, or True to keep that of the parent, or False to\n",
      " |      omit them. The default is to keep the original features (with the\n",
      " |      strand and locations adjusted).\n",
      " |      \n",
      " |      You can also specify both the returned record's annotations and\n",
      " |      letter_annotations as dictionaries, True to keep that of the parent,\n",
      " |      or False to omit them. The default is to keep the original\n",
      " |      annotations (with the letter annotations reversed).\n",
      " |      \n",
      " |      To show what happens to the pre-letter annotations, consider an\n",
      " |      example Solexa variant FASTQ file with a single entry, which we'll\n",
      " |      read in as a SeqRecord:\n",
      " |      \n",
      " |      >>> from Bio import SeqIO\n",
      " |      >>> record = SeqIO.read(\"Quality/solexa_faked.fastq\", \"fastq-solexa\")\n",
      " |      >>> print(\"%s %s\" % (record.id, record.seq))\n",
      " |      slxa_0001_1_0001_01 ACGTACGTACGTACGTACGTACGTACGTACGTACGTACGTNNNNNN\n",
      " |      >>> print(list(record.letter_annotations))\n",
      " |      ['solexa_quality']\n",
      " |      >>> print(record.letter_annotations[\"solexa_quality\"])\n",
      " |      [40, 39, 38, 37, 36, 35, 34, 33, 32, 31, 30, 29, 28, 27, 26, 25, 24, 23, 22, 21, 20, 19, 18, 17, 16, 15, 14, 13, 12, 11, 10, 9, 8, 7, 6, 5, 4, 3, 2, 1, 0, -1, -2, -3, -4, -5]\n",
      " |      \n",
      " |      Now take the reverse complement, here we explicitly give a new\n",
      " |      identifier (the old identifier with a suffix):\n",
      " |      \n",
      " |      >>> rc_record = record.reverse_complement(id=record.id + \"_rc\")\n",
      " |      >>> print(\"%s %s\" % (rc_record.id, rc_record.seq))\n",
      " |      slxa_0001_1_0001_01_rc NNNNNNACGTACGTACGTACGTACGTACGTACGTACGTACGTACGT\n",
      " |      \n",
      " |      Notice that the per-letter-annotations have also been reversed,\n",
      " |      although this may not be appropriate for all cases.\n",
      " |      \n",
      " |      >>> print(rc_record.letter_annotations[\"solexa_quality\"])\n",
      " |      [-5, -4, -3, -2, -1, 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40]\n",
      " |      \n",
      " |      Now for the features, we need a different example. Parsing a GenBank\n",
      " |      file is probably the easiest way to get an nice example with features\n",
      " |      in it...\n",
      " |      \n",
      " |      >>> from Bio import SeqIO\n",
      " |      >>> with open(\"GenBank/pBAD30.gb\") as handle:\n",
      " |      ...     plasmid = SeqIO.read(handle, \"gb\")\n",
      " |      >>> print(\"%s %i\" % (plasmid.id, len(plasmid)))\n",
      " |      pBAD30 4923\n",
      " |      >>> plasmid.seq\n",
      " |      Seq('GCTAGCGGAGTGTATACTGGCTTACTATGTTGGCACTGATGAGGGTGTCAGTGA...ATG', IUPACAmbiguousDNA())\n",
      " |      >>> len(plasmid.features)\n",
      " |      13\n",
      " |      \n",
      " |      Now, let's take the reverse complement of this whole plasmid:\n",
      " |      \n",
      " |      >>> rc_plasmid = plasmid.reverse_complement(id=plasmid.id+\"_rc\")\n",
      " |      >>> print(\"%s %i\" % (rc_plasmid.id, len(rc_plasmid)))\n",
      " |      pBAD30_rc 4923\n",
      " |      >>> rc_plasmid.seq\n",
      " |      Seq('CATGGGCAAATATTATACGCAAGGCGACAAGGTGCTGATGCCGCTGGCGATTCA...AGC', IUPACAmbiguousDNA())\n",
      " |      >>> len(rc_plasmid.features)\n",
      " |      13\n",
      " |      \n",
      " |      Let's compare the first CDS feature - it has gone from being the\n",
      " |      second feature (index 1) to the second last feature (index -2), its\n",
      " |      strand has changed, and the location switched round.\n",
      " |      \n",
      " |      >>> print(plasmid.features[1])\n",
      " |      type: CDS\n",
      " |      location: [1081:1960](-)\n",
      " |      qualifiers:\n",
      " |          Key: label, Value: ['araC']\n",
      " |          Key: note, Value: ['araC regulator of the arabinose BAD promoter']\n",
      " |          Key: vntifkey, Value: ['4']\n",
      " |      <BLANKLINE>\n",
      " |      >>> print(rc_plasmid.features[-2])\n",
      " |      type: CDS\n",
      " |      location: [2963:3842](+)\n",
      " |      qualifiers:\n",
      " |          Key: label, Value: ['araC']\n",
      " |          Key: note, Value: ['araC regulator of the arabinose BAD promoter']\n",
      " |          Key: vntifkey, Value: ['4']\n",
      " |      <BLANKLINE>\n",
      " |      \n",
      " |      You can check this new location, based on the length of the plasmid:\n",
      " |      \n",
      " |      >>> len(plasmid) - 1081\n",
      " |      3842\n",
      " |      >>> len(plasmid) - 1960\n",
      " |      2963\n",
      " |      \n",
      " |      Note that if the SeqFeature annotation includes any strand specific\n",
      " |      information (e.g. base changes for a SNP), this information is not\n",
      " |      amended, and would need correction after the reverse complement.\n",
      " |      \n",
      " |      Note trying to reverse complement a protein SeqRecord raises an\n",
      " |      exception:\n",
      " |      \n",
      " |      >>> from Bio.SeqRecord import SeqRecord\n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> from Bio.Alphabet import IUPAC\n",
      " |      >>> protein_rec = SeqRecord(Seq(\"MAIVMGR\", IUPAC.protein), id=\"Test\")\n",
      " |      >>> protein_rec.reverse_complement()\n",
      " |      Traceback (most recent call last):\n",
      " |         ...\n",
      " |      ValueError: Proteins do not have complements!\n",
      " |      \n",
      " |      Also note you can reverse complement a SeqRecord using a MutableSeq:\n",
      " |      \n",
      " |      >>> from Bio.SeqRecord import SeqRecord\n",
      " |      >>> from Bio.Seq import MutableSeq\n",
      " |      >>> from Bio.Alphabet import generic_dna\n",
      " |      >>> rec = SeqRecord(MutableSeq(\"ACGT\", generic_dna), id=\"Test\")\n",
      " |      >>> rec.seq[0] = \"T\"\n",
      " |      >>> print(\"%s %s\" % (rec.id, rec.seq))\n",
      " |      Test TCGT\n",
      " |      >>> rc = rec.reverse_complement(id=True)\n",
      " |      >>> print(\"%s %s\" % (rc.id, rc.seq))\n",
      " |      Test ACGA\n",
      " |  \n",
      " |  translate(self, table='Standard', stop_symbol='*', to_stop=False, cds=False, gap=None, id=False, name=False, description=False, features=False, annotations=False, letter_annotations=False, dbxrefs=False)\n",
      " |      Return new SeqRecord with translated sequence.\n",
      " |      \n",
      " |      This calls the record's .seq.translate() method (which describes\n",
      " |      the translation related arguments, like table for the genetic code),\n",
      " |      \n",
      " |      By default the new record does NOT preserve the sequence identifier,\n",
      " |      name, description, general annotation or database cross-references -\n",
      " |      these are unlikely to apply to the translated sequence.\n",
      " |      \n",
      " |      You can specify the returned record's id, name and description as\n",
      " |      strings, or True to keep that of the parent, or False for a default.\n",
      " |      \n",
      " |      You can specify the returned record's features with a list of\n",
      " |      SeqFeature objects, or False (default) to omit them.\n",
      " |      \n",
      " |      You can also specify both the returned record's annotations and\n",
      " |      letter_annotations as dictionaries, True to keep that of the parent\n",
      " |      (annotations only), or False (default) to omit them.\n",
      " |      \n",
      " |      e.g. Loading a FASTA gene and translating it,\n",
      " |      \n",
      " |      >>> from Bio import SeqIO\n",
      " |      >>> gene_record = SeqIO.read(\"Fasta/sweetpea.nu\", \"fasta\")\n",
      " |      >>> print(gene_record.format(\"fasta\"))\n",
      " |      >gi|3176602|gb|U78617.1|LOU78617 Lathyrus odoratus phytochrome A (PHYA) gene, partial cds\n",
      " |      CAGGCTGCGCGGTTTCTATTTATGAAGAACAAGGTCCGTATGATAGTTGATTGTCATGCA\n",
      " |      AAACATGTGAAGGTTCTTCAAGACGAAAAACTCCCATTTGATTTGACTCTGTGCGGTTCG\n",
      " |      ACCTTAAGAGCTCCACATAGTTGCCATTTGCAGTACATGGCTAACATGGATTCAATTGCT\n",
      " |      TCATTGGTTATGGCAGTGGTCGTCAATGACAGCGATGAAGATGGAGATAGCCGTGACGCA\n",
      " |      GTTCTACCACAAAAGAAAAAGAGACTTTGGGGTTTGGTAGTTTGTCATAACACTACTCCG\n",
      " |      AGGTTTGTT\n",
      " |      <BLANKLINE>\n",
      " |      \n",
      " |      And now translating the record, specifying the new ID and description:\n",
      " |      \n",
      " |      >>> protein_record = gene_record.translate(table=11,\n",
      " |      ...                                        id=\"phya\",\n",
      " |      ...                                        description=\"translation\")\n",
      " |      >>> print(protein_record.format(\"fasta\"))\n",
      " |      >phya translation\n",
      " |      QAARFLFMKNKVRMIVDCHAKHVKVLQDEKLPFDLTLCGSTLRAPHSCHLQYMANMDSIA\n",
      " |      SLVMAVVVNDSDEDGDSRDAVLPQKKKRLWGLVVCHNTTPRFV\n",
      " |      <BLANKLINE>\n",
      " |  \n",
      " |  upper(self)\n",
      " |      Return a copy of the record with an upper case sequence.\n",
      " |      \n",
      " |      All the annotation is preserved unchanged. e.g.\n",
      " |      \n",
      " |      >>> from Bio.Alphabet import generic_dna\n",
      " |      >>> from Bio.Seq import Seq\n",
      " |      >>> from Bio.SeqRecord import SeqRecord\n",
      " |      >>> record = SeqRecord(Seq(\"acgtACGT\", generic_dna), id=\"Test\",\n",
      " |      ...                    description = \"Made up for this example\")\n",
      " |      >>> record.letter_annotations[\"phred_quality\"] = [1, 2, 3, 4, 5, 6, 7, 8]\n",
      " |      >>> print(record.upper().format(\"fastq\"))\n",
      " |      @Test Made up for this example\n",
      " |      ACGTACGT\n",
      " |      +\n",
      " |      \"#$%&'()\n",
      " |      <BLANKLINE>\n",
      " |      \n",
      " |      Naturally, there is a matching lower method:\n",
      " |      \n",
      " |      >>> print(record.lower().format(\"fastq\"))\n",
      " |      @Test Made up for this example\n",
      " |      acgtacgt\n",
      " |      +\n",
      " |      \"#$%&'()\n",
      " |      <BLANKLINE>\n",
      " |  \n",
      " |  ----------------------------------------------------------------------\n",
      " |  Data descriptors defined here:\n",
      " |  \n",
      " |  __dict__\n",
      " |      dictionary for instance variables (if defined)\n",
      " |  \n",
      " |  __weakref__\n",
      " |      list of weak references to the object (if defined)\n",
      " |  \n",
      " |  letter_annotations\n",
      " |      Dictionary of per-letter-annotation for the sequence.\n",
      " |      \n",
      " |      For example, this can hold quality scores used in FASTQ or QUAL files.\n",
      " |      Consider this example using Bio.SeqIO to read in an example Solexa\n",
      " |      variant FASTQ file as a SeqRecord:\n",
      " |      \n",
      " |      >>> from Bio import SeqIO\n",
      " |      >>> record = SeqIO.read(\"Quality/solexa_faked.fastq\", \"fastq-solexa\")\n",
      " |      >>> print(\"%s %s\" % (record.id, record.seq))\n",
      " |      slxa_0001_1_0001_01 ACGTACGTACGTACGTACGTACGTACGTACGTACGTACGTNNNNNN\n",
      " |      >>> print(list(record.letter_annotations))\n",
      " |      ['solexa_quality']\n",
      " |      >>> print(record.letter_annotations[\"solexa_quality\"])\n",
      " |      [40, 39, 38, 37, 36, 35, 34, 33, 32, 31, 30, 29, 28, 27, 26, 25, 24, 23, 22, 21, 20, 19, 18, 17, 16, 15, 14, 13, 12, 11, 10, 9, 8, 7, 6, 5, 4, 3, 2, 1, 0, -1, -2, -3, -4, -5]\n",
      " |      \n",
      " |      The letter_annotations get sliced automatically if you slice the\n",
      " |      parent SeqRecord, for example taking the last ten bases:\n",
      " |      \n",
      " |      >>> sub_record = record[-10:]\n",
      " |      >>> print(\"%s %s\" % (sub_record.id, sub_record.seq))\n",
      " |      slxa_0001_1_0001_01 ACGTNNNNNN\n",
      " |      >>> print(sub_record.letter_annotations[\"solexa_quality\"])\n",
      " |      [4, 3, 2, 1, 0, -1, -2, -3, -4, -5]\n",
      " |      \n",
      " |      Any python sequence (i.e. list, tuple or string) can be recorded in\n",
      " |      the SeqRecord's letter_annotations dictionary as long as the length\n",
      " |      matches that of the SeqRecord's sequence.  e.g.\n",
      " |      \n",
      " |      >>> len(sub_record.letter_annotations)\n",
      " |      1\n",
      " |      >>> sub_record.letter_annotations[\"dummy\"] = \"abcdefghij\"\n",
      " |      >>> len(sub_record.letter_annotations)\n",
      " |      2\n",
      " |      \n",
      " |      You can delete entries from the letter_annotations dictionary as usual:\n",
      " |      \n",
      " |      >>> del sub_record.letter_annotations[\"solexa_quality\"]\n",
      " |      >>> sub_record.letter_annotations\n",
      " |      {'dummy': 'abcdefghij'}\n",
      " |      \n",
      " |      You can completely clear the dictionary easily as follows:\n",
      " |      \n",
      " |      >>> sub_record.letter_annotations = {}\n",
      " |      >>> sub_record.letter_annotations\n",
      " |      {}\n",
      " |      \n",
      " |      Note that if replacing the record's sequence with a sequence of a\n",
      " |      different length you must first clear the letter_annotations dict.\n",
      " |  \n",
      " |  seq\n",
      " |      The sequence itself, as a Seq or MutableSeq object.\n",
      " |  \n",
      " |  ----------------------------------------------------------------------\n",
      " |  Data and other attributes defined here:\n",
      " |  \n",
      " |  __hash__ = None\n",
      "\n"
     ]
    }
   ],
   "source": [
    "from Bio.SeqRecord import SeqRecord\n",
    "help(SeqRecord)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2.1 Criando um objeto SeqRecord\n",
    "Vamos colocar a mão na massa e criar um objeto SeqRecord. Primeiro criamos um objeto Seq e em seguida criamos um objeto SeqRecord passando como parâmetro o objeto Seq."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "SeqRecord(seq=Seq('GATC'), id='<unknown id>', name='<unknown name>', description='<unknown description>', dbxrefs=[])"
      ]
     },
     "execution_count": 35,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "simple_seq = Seq(\"GATC\")\n",
    "\n",
    "simple_seq_r = SeqRecord(simple_seq)\n",
    "simple_seq_r"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Percebam que esse `SeqRecord` tem apenas a sequência, mas podemos ir adicionando os outros atributos como id e descrição."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "SeqRecord(seq=Seq('GATC'), id='AC12345', name='<unknown name>', description='Made up sequence I wish I could write a paper about', dbxrefs=[])"
      ]
     },
     "execution_count": 36,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "simple_seq_r.id = \"AC12345\"\n",
    "simple_seq_r.description = \"Made up sequence I wish I could write a paper about\"\n",
    "simple_seq_r"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 2.1.1 Adicionando anotações\n",
    "Também podemos adicionar uma anotação para esse `SeqRecord`. O atributo anotação é um dicionário, assim qualquer tipo de anotação pode ser adicionada, sendo composto por uma chave e um valor. No exemplo abaixo adicionamos uma anotação onde a chave é \"evidence\" e o valor é \"None. I just made it up.\""
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "{'evidence': 'None. I just made it up.'}\n",
      "None. I just made it up.\n"
     ]
    }
   ],
   "source": [
    "simple_seq_r.annotations[\"evidence\"] = \"None. I just made it up.\"\n",
    "print(simple_seq_r.annotations)\n",
    "\n",
    "print(simple_seq_r.annotations[\"evidence\"])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 2.1.2 Adicionando anotações por letras\n",
    "Outro tipo de anotação que pode ser bem interessante é a anotação das letras da sequência que também é um dicionário. Um exemplo prático é adicionar a qualidade de cada uma das bases. Assim adicionamos a chave \"phred_quality\" e passamos como valor uma lista com a qualidade de cada base."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "{'phred_quality': [40, 40, 38, 40]}\n",
      "[40, 40, 38, 40]\n"
     ]
    }
   ],
   "source": [
    "simple_seq_r.letter_annotations[\"phred_quality\"] = [40, 40, 38, 40]\n",
    "print(simple_seq_r.letter_annotations)\n",
    "\n",
    "print(simple_seq_r.letter_annotations[\"phred_quality\"])\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2.2 Carregando um SeqRecord completo\n",
    "Vamos carregar um exemplo com mais atributos. Aqui vamos utilizar a Biblioteca `SeqIO` para carregar as informações diretamente de um arquivo com anotação do Genbank. Mais para frente vamos falar sobre o `SeqIO`, mas por enquanto vamos apenas entender que o comando abaixo faz um parser do arquivo Genbank [nc_005816](https://www.ncbi.nlm.nih.gov/nuccore/NC_005816.1?report=gbwithparts&log$=seqview) e cria um objeto SeqRecord com as informações."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "SeqRecord(seq=Seq('TGTAACGAACGGTGCAATAGTGATCCACACCCAACGCCTGAAATCAGATCCAGG...CTG', IUPACAmbiguousDNA()), id='NC_005816.1', name='NC_005816', description='Yersinia pestis biovar Microtus str. 91001 plasmid pPCP1, complete sequence', dbxrefs=['BioProject:PRJNA224116', 'BioSample:SAMN02602970', 'Assembly:GCF_000007885.1'])"
      ]
     },
     "execution_count": 31,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from Bio import SeqIO\n",
    "record = SeqIO.read(\"sequence.gb\", \"genbank\")\n",
    "record"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "{'molecule_type': 'DNA',\n",
       " 'topology': 'circular',\n",
       " 'data_file_division': 'CON',\n",
       " 'date': '11-JAN-2018',\n",
       " 'accessions': ['NC_005816'],\n",
       " 'sequence_version': 1,\n",
       " 'gi': '45478711',\n",
       " 'keywords': ['RefSeq'],\n",
       " 'source': 'Yersinia pestis biovar Microtus str. 91001',\n",
       " 'organism': 'Yersinia pestis biovar Microtus str. 91001',\n",
       " 'taxonomy': ['Bacteria',\n",
       "  'Proteobacteria',\n",
       "  'Gammaproteobacteria',\n",
       "  'Enterobacterales',\n",
       "  'Yersiniaceae',\n",
       "  'Yersinia'],\n",
       " 'references': [Reference(title='Genetics of metabolic variations between Yersinia pestis biovars and the proposal of a new biovar, microtus', ...),\n",
       "  Reference(title='Complete genome sequence of Yersinia pestis strain 91001, an isolate avirulent to humans', ...),\n",
       "  Reference(title='Direct Submission', ...)],\n",
       " 'comment': 'REFSEQ INFORMATION: The reference sequence was derived from\\nAE017046.\\nAnnotation was added by the NCBI Prokaryotic Genome Annotation\\nPipeline (released 2013). Information about the Pipeline can be\\nfound here: https://www.ncbi.nlm.nih.gov/genome/annotation_prok/\\n                                     Annotation Pipeline\\n                                     set; GeneMarkS+\\n                                     repeat_region\\nCOMPLETENESS: full length.',\n",
       " 'structured_comment': OrderedDict([('Genome-Annotation-Data',\n",
       "               OrderedDict([('Annotation Provider', 'NCBI'),\n",
       "                            ('Annotation Date', '01/10/2018 21:05:30'),\n",
       "                            ('Annotation Pipeline', 'NCBI Prokaryotic Genome'),\n",
       "                            ('Annotation Method',\n",
       "                             'Best-placed reference protein'),\n",
       "                            ('Annotation Software revision', '4.4'),\n",
       "                            ('Features Annotated',\n",
       "                             'Gene; CDS; rRNA; tRNA; ncRNA;'),\n",
       "                            ('Genes (total)', '4,592'),\n",
       "                            ('CDS (total)', '4,489'),\n",
       "                            ('Genes (coding)', '4,280'),\n",
       "                            ('CDS (coding)', '4,280'),\n",
       "                            ('Genes (RNA)', '103'),\n",
       "                            ('rRNAs', '8, 7, 7 (5S, 16S, 23S)'),\n",
       "                            ('complete rRNAs', '8, 7, 7 (5S, 16S, 23S)'),\n",
       "                            ('tRNAs', '71'),\n",
       "                            ('ncRNAs', '10'),\n",
       "                            ('Pseudo Genes (total)', '209'),\n",
       "                            ('Pseudo Genes (ambiguous residues)', '0 of 209'),\n",
       "                            ('Pseudo Genes (frameshifted)', '91 of 209'),\n",
       "                            ('Pseudo Genes (incomplete)', '125 of 209'),\n",
       "                            ('Pseudo Genes (internal stop)', '31 of 209'),\n",
       "                            ('Pseudo Genes (multiple problems)', '35 of 209'),\n",
       "                            ('CRISPR Arrays', '1')]))])}"
      ]
     },
     "execution_count": 32,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "record.annotations"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "['Bacteria',\n",
       " 'Proteobacteria',\n",
       " 'Gammaproteobacteria',\n",
       " 'Enterobacterales',\n",
       " 'Yersiniaceae',\n",
       " 'Yersinia']"
      ]
     },
     "execution_count": 33,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "record.annotations[\"taxonomy\"]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "['BioProject:PRJNA224116',\n",
       " 'BioSample:SAMN02602970',\n",
       " 'Assembly:GCF_000007885.1']"
      ]
     },
     "execution_count": 34,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "record.dbxrefs"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2.3 Acessando as SeqFeature do objeto SeqRecord\n",
    "\n",
    "A classe `SeqRecord` tem o atributo features. Este atributo recebe um objeto da classe `SeqFeature`, sendo que ele basicamente marca regiões de interesse no sua sequência, por exemplo, a região de um gene."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "type: source\n",
      "location: [0:9609](+)\n",
      "qualifiers:\n",
      "    Key: biovar, Value: ['Microtus']\n",
      "    Key: db_xref, Value: ['taxon:229193']\n",
      "    Key: mol_type, Value: ['genomic DNA']\n",
      "    Key: organism, Value: ['Yersinia pestis biovar Microtus str. 91001']\n",
      "    Key: plasmid, Value: ['pPCP1']\n",
      "    Key: strain, Value: ['91001']\n",
      "    Key: sub_species, Value: ['microtus']\n",
      "\n",
      "type: gene\n",
      "location: [86:1109](+)\n",
      "qualifiers:\n",
      "    Key: locus_tag, Value: ['YP_RS22210']\n",
      "    Key: old_locus_tag, Value: ['pPCP01', 'YP_pPCP01']\n",
      "\n",
      "type: CDS\n",
      "location: [86:1109](+)\n",
      "qualifiers:\n",
      "    Key: codon_start, Value: ['1']\n",
      "    Key: db_xref, Value: ['GI:446178089']\n",
      "    Key: inference, Value: ['COORDINATES: similar to AA sequence:RefSeq:WP_000627495.1']\n",
      "    Key: locus_tag, Value: ['YP_RS22210']\n",
      "    Key: note, Value: ['Derived by automated computational analysis using gene prediction method: Protein Homology.']\n",
      "    Key: old_locus_tag, Value: ['pPCP01', 'YP_pPCP01']\n",
      "    Key: product, Value: ['IS21-like element IS100 family transposase']\n",
      "    Key: protein_id, Value: ['WP_000255944.1']\n",
      "    Key: transl_table, Value: ['11']\n",
      "    Key: translation, Value: ['MVTFETVMEIKILHKQGMSSRAIARELGISRNTVKRYLQAKSEPPKYTPRPAVASLLDEYRDYIRQRIADAHPYKIPATVIAREIRDQGYRGGMTILRAFIRSLSVPQEQEPAVRFETEPGRQMQVDWGTMRNGRSPLHVFVAVLGYSRMLYIEFTDNMRYDTLETCHRNAFRFFGGVPREVLYDNMKTVVLQRDAYQTGQHRFHPSLWQFGKEMGFSPRLCRPFRAQTKGKVERMVQYTRNSFYIPLMTRLRPMGITVDVETANRHGLRWLHDVANQRKHETIQARPCDRWLEEQQSMLALPPEKKEYDVHLDENLVNFDKHPLHHPLSIYDSFCRGVA']\n",
      "\n",
      "type: gene\n",
      "location: [1105:1888](+)\n",
      "qualifiers:\n",
      "    Key: locus_tag, Value: ['YP_RS22215']\n",
      "    Key: old_locus_tag, Value: ['pPCP02', 'YP_pPCP02']\n",
      "\n",
      "type: CDS\n",
      "location: [1105:1888](+)\n",
      "qualifiers:\n",
      "    Key: codon_start, Value: ['1']\n",
      "    Key: db_xref, Value: ['GI:485682829']\n",
      "    Key: inference, Value: ['COORDINATES: similar to AA sequence:RefSeq:NP_395228.1']\n",
      "    Key: locus_tag, Value: ['YP_RS22215']\n",
      "    Key: note, Value: ['Derived by automated computational analysis using gene prediction method: Protein Homology.']\n",
      "    Key: old_locus_tag, Value: ['pPCP02', 'YP_pPCP02']\n",
      "    Key: product, Value: ['AAA family ATPase']\n",
      "    Key: protein_id, Value: ['WP_001317493.1']\n",
      "    Key: transl_table, Value: ['11']\n",
      "    Key: translation, Value: ['MMMELQHQRLMALAGQLQLESLISAAPALSQQAVDQEWSYMDFLEHLLHEEKLARHQRKQAMYTRMAAFPAVKTFEEYDFTFATGAPQKQLQSLRSLSFIERNENIVLLGPSGVGKTHLAIAMGYEAVRAGIKVRFTTAADLLLQLSTAQRQGRYKTTLQRGVMAPRLLIIDEIGYLPFSQEEAKLFFQVIAKRYEKSAMILTSNLPFGQWDQTFAGDAALTSAMLDRILHHSHVVQIKGESYRLRQKRKAGVIAEANPE']\n",
      "\n",
      "type: gene\n",
      "location: [2924:3119](+)\n",
      "qualifiers:\n",
      "    Key: locus_tag, Value: ['YP_RS22220']\n",
      "    Key: old_locus_tag, Value: ['pPCP03', 'YP_pPCP03']\n",
      "\n",
      "type: CDS\n",
      "location: [2924:3119](+)\n",
      "qualifiers:\n",
      "    Key: codon_start, Value: ['1']\n",
      "    Key: db_xref, Value: ['GI:488142225']\n",
      "    Key: inference, Value: ['COORDINATES: similar to AA sequence:RefSeq:NP_395229.1']\n",
      "    Key: locus_tag, Value: ['YP_RS22220']\n",
      "    Key: note, Value: ['Derived by automated computational analysis using gene prediction method: Protein Homology.']\n",
      "    Key: old_locus_tag, Value: ['pPCP03', 'YP_pPCP03']\n",
      "    Key: product, Value: ['Rop family plasmid primer RNA-binding protein']\n",
      "    Key: protein_id, Value: ['WP_002213433.1']\n",
      "    Key: transl_table, Value: ['11']\n",
      "    Key: translation, Value: ['MNKQQQTALNMARFIRSQSLILLEKLDALDADEQAAMCERLHELAEELQNSIQARFEAESETGT']\n",
      "\n",
      "type: gene\n",
      "location: [<3358:3442](+)\n",
      "qualifiers:\n",
      "    Key: locus_tag, Value: ['YP_RS24050']\n",
      "    Key: pseudo, Value: ['']\n",
      "\n",
      "type: CDS\n",
      "location: [<3358:3442](+)\n",
      "qualifiers:\n",
      "    Key: codon_start, Value: ['1']\n",
      "    Key: inference, Value: ['COORDINATES: similar to AA sequence:RefSeq:WP_001985485.1']\n",
      "    Key: locus_tag, Value: ['YP_RS24050']\n",
      "    Key: note, Value: ['incomplete; partial in the middle of a contig; missing start; Derived by automated computational analysis using gene prediction method: Protein Homology.']\n",
      "    Key: product, Value: ['SocE']\n",
      "    Key: pseudo, Value: ['']\n",
      "    Key: transl_table, Value: ['11']\n",
      "\n",
      "type: gene\n",
      "location: [4354:4780](+)\n",
      "qualifiers:\n",
      "    Key: locus_tag, Value: ['YP_RS22225']\n",
      "    Key: old_locus_tag, Value: ['pPCP05', 'YP_pPCP05']\n",
      "\n",
      "type: CDS\n",
      "location: [4354:4780](+)\n",
      "qualifiers:\n",
      "    Key: codon_start, Value: ['1']\n",
      "    Key: db_xref, Value: ['GI:488147304']\n",
      "    Key: inference, Value: ['COORDINATES: similar to AA sequence:RefSeq:NP_395230.1']\n",
      "    Key: locus_tag, Value: ['YP_RS22225']\n",
      "    Key: note, Value: ['Derived by automated computational analysis using gene prediction method: Protein Homology.']\n",
      "    Key: old_locus_tag, Value: ['pPCP05', 'YP_pPCP05']\n",
      "    Key: product, Value: ['pesticin immunity protein']\n",
      "    Key: protein_id, Value: ['WP_002218512.1']\n",
      "    Key: transl_table, Value: ['11']\n",
      "    Key: translation, Value: ['MISKLFCLALIFLSSSGLAEKNTYTAKDILQNLELNTFGNSLSHGIYGKQTTFKQTEFTNIKSNTKKHIALINKDNSWMISLKILGIKRDEYTVCFEDFSLIRPPTYVAIHPLLIKKVKSGNFIVVKEIKKSIPGCTVYYH']\n",
      "\n",
      "type: gene\n",
      "location: [4814:5888](-)\n",
      "qualifiers:\n",
      "    Key: locus_tag, Value: ['YP_RS22230']\n",
      "    Key: old_locus_tag, Value: ['pPCP06', 'YP_pPCP06']\n",
      "\n",
      "type: CDS\n",
      "location: [4814:5888](-)\n",
      "qualifiers:\n",
      "    Key: codon_start, Value: ['1']\n",
      "    Key: db_xref, Value: ['GI:488147301']\n",
      "    Key: inference, Value: ['COORDINATES: similar to AA sequence:RefSeq:NP_395231.1']\n",
      "    Key: locus_tag, Value: ['YP_RS22230']\n",
      "    Key: note, Value: ['Derived by automated computational analysis using gene prediction method: Protein Homology.']\n",
      "    Key: old_locus_tag, Value: ['pPCP06', 'YP_pPCP06']\n",
      "    Key: product, Value: ['pesticin']\n",
      "    Key: protein_id, Value: ['WP_002218509.1']\n",
      "    Key: transl_table, Value: ['11']\n",
      "    Key: translation, Value: ['MSDTMVVNGSGGVPAFLFSGSTLSSYRPNFEANSITIALPHYVDLPGRSNFKLMYIMGFPIDTEMEKDSEYSNKIRQESKISKTEGTVSYEQKITVETGQEKDGVKVYRVMVLEGTIAESIEHLDKKENEDILNNNRNRIVLADNTVINFDNISQLKEFLRRSVNIVDHDIFSSNGFEGFNPTSHFPSNPSSDYFNSTGVTFGSGVDLGQRSKQDLLNDGVPQYIADRLDGYYMLRGKEAYDKVRTAPLTLSDNEAHLLSNIYIDKFSHKIEGLFNDANIGLRFSDLPLRTRTALVSIGYQKGFKLSRTAPTVWNKVIAKDWNGLVNAFNNIVDGMSDRRKREGALVQKDIDSGLLK']\n",
      "\n",
      "type: gene\n",
      "location: [6004:6421](+)\n",
      "qualifiers:\n",
      "    Key: locus_tag, Value: ['YP_RS22235']\n",
      "    Key: old_locus_tag, Value: ['pPCP07']\n",
      "\n",
      "type: CDS\n",
      "location: [6004:6421](+)\n",
      "qualifiers:\n",
      "    Key: codon_start, Value: ['1']\n",
      "    Key: db_xref, Value: ['GI:488142243']\n",
      "    Key: inference, Value: ['COORDINATES: similar to AA sequence:RefSeq:NP_395232.1']\n",
      "    Key: locus_tag, Value: ['YP_RS22235']\n",
      "    Key: note, Value: ['Derived by automated computational analysis using gene prediction method: Protein Homology.']\n",
      "    Key: old_locus_tag, Value: ['pPCP07']\n",
      "    Key: product, Value: ['hypothetical protein']\n",
      "    Key: protein_id, Value: ['WP_002213451.1']\n",
      "    Key: transl_table, Value: ['11']\n",
      "    Key: translation, Value: ['MKFHFCDLNHSYKNQEGKIRSRKTAPGNIRKKQKGDNVSKTKSGRHRLSKTDKRLLAALVVAGYEERTARDLIQKHVYTLTQADLRHLVSEISNGVGQSQAYDAIYQARRIRLARKYLSGKKPEGVEPREGQEREDLP']\n",
      "\n",
      "type: gene\n",
      "location: [6663:7602](+)\n",
      "qualifiers:\n",
      "    Key: locus_tag, Value: ['YP_RS22240']\n",
      "    Key: old_locus_tag, Value: ['pPCP08', 'YP_pPCP08']\n",
      "\n",
      "type: CDS\n",
      "location: [6663:7602](+)\n",
      "qualifiers:\n",
      "    Key: codon_start, Value: ['1']\n",
      "    Key: db_xref, Value: ['GI:488159966']\n",
      "    Key: inference, Value: ['COORDINATES: similar to AA sequence:RefSeq:NP_395233.1']\n",
      "    Key: locus_tag, Value: ['YP_RS22240']\n",
      "    Key: note, Value: ['Derived by automated computational analysis using gene prediction method: Protein Homology.']\n",
      "    Key: old_locus_tag, Value: ['pPCP08', 'YP_pPCP08']\n",
      "    Key: product, Value: ['omptin family plasminogen activator Pla']\n",
      "    Key: protein_id, Value: ['WP_002231174.1']\n",
      "    Key: transl_table, Value: ['11']\n",
      "    Key: translation, Value: ['MKKSSIVATIITILSGSANAASSQLIPNISPDSFTVAASTGMLSGKSHEMLYDAETGRKISQLDWKIKNVAILKGDISWDPYSFLTLNARGWTSLASGSGNMDDYDWMNENQSEWTDHSSHPATNVNHANEYDLNVKGWLLQDENYKAGITAGYQETRFSWTATGGSYSYNNGAYTGNFPKGVRVIGYNQRFSMPYIGLAGQYRINDFELNALFKFSDWVRAHDNDEHYMRDLTFREKTSGSRYYGTVINAGYYVTPNAKVFAEFTYSKYDEGKGGTQIIDKNSGDSVSIGGDAAGISNKNYTVTAGLQYRF']\n",
      "\n",
      "type: gene\n",
      "location: [7788:8088](-)\n",
      "qualifiers:\n",
      "    Key: locus_tag, Value: ['YP_RS22245']\n",
      "    Key: old_locus_tag, Value: ['pPCP09', 'YP_pPCP09']\n",
      "\n",
      "type: CDS\n",
      "location: [7788:8088](-)\n",
      "qualifiers:\n",
      "    Key: codon_start, Value: ['1']\n",
      "    Key: db_xref, Value: ['GI:488142238']\n",
      "    Key: inference, Value: ['COORDINATES: similar to AA sequence:RefSeq:NP_395234.1']\n",
      "    Key: locus_tag, Value: ['YP_RS22245']\n",
      "    Key: note, Value: ['Derived by automated computational analysis using gene prediction method: Protein Homology.']\n",
      "    Key: old_locus_tag, Value: ['pPCP09', 'YP_pPCP09']\n",
      "    Key: product, Value: ['helix-turn-helix domain-containing protein']\n",
      "    Key: protein_id, Value: ['WP_002213446.1']\n",
      "    Key: transl_table, Value: ['11']\n",
      "    Key: translation, Value: ['MRTLDEVIASRSPESQTRIKEMADEMILEVGLQMMREELQLSQKQVAEAMGISQPAVTKLEQRGNDLKLATLKRYVEAMGGKLSLDVELPTGRRVAFHV']\n",
      "\n",
      "type: gene\n",
      "location: [8087:8435](-)\n",
      "qualifiers:\n",
      "    Key: locus_tag, Value: ['YP_RS22250']\n",
      "    Key: old_locus_tag, Value: ['pPCP10']\n",
      "\n",
      "type: CDS\n",
      "location: [8087:8435](-)\n",
      "qualifiers:\n",
      "    Key: codon_start, Value: ['1']\n",
      "    Key: db_xref, Value: ['GI:488151648']\n",
      "    Key: inference, Value: ['COORDINATES: similar to AA sequence:RefSeq:NP_395235.1']\n",
      "    Key: locus_tag, Value: ['YP_RS22250']\n",
      "    Key: note, Value: ['Derived by automated computational analysis using gene prediction method: Protein Homology.']\n",
      "    Key: old_locus_tag, Value: ['pPCP10']\n",
      "    Key: product, Value: ['type II toxin-antitoxin system RelE/ParE family toxin']\n",
      "    Key: protein_id, Value: ['WP_002222856.1']\n",
      "    Key: transl_table, Value: ['11']\n",
      "    Key: translation, Value: ['MWMVLFSQRFDDWLNEQEDALQEKVLADLKKLQVYGPELPRPYADTVKGSRYKNMKELRVQFSGRPIRAFYAFDPIRRAIVLCAGDKSNDKRFYEKLVRIAEDEFTAHLNTLESK']\n",
      "\n",
      "type: gene\n",
      "location: [8507:8639](-)\n",
      "qualifiers:\n",
      "    Key: locus_tag, Value: ['YP_RS23685']\n",
      "\n",
      "type: CDS\n",
      "location: [8507:8639](-)\n",
      "qualifiers:\n",
      "    Key: codon_start, Value: ['1']\n",
      "    Key: db_xref, Value: ['GI:1101975813']\n",
      "    Key: inference, Value: ['COORDINATES: similar to AA sequence:RefSeq:WP_001503282.1']\n",
      "    Key: locus_tag, Value: ['YP_RS23685']\n",
      "    Key: note, Value: ['Derived by automated computational analysis using gene prediction method: Protein Homology.']\n",
      "    Key: product, Value: ['hypothetical protein']\n",
      "    Key: protein_id, Value: ['WP_071525773.1']\n",
      "    Key: transl_table, Value: ['11']\n",
      "    Key: translation, Value: ['MRLKFIRRDRGGFSGGLLLFLPVCCRDRAERISVVRTNCRYVK']\n",
      "\n",
      "type: gene\n",
      "location: [9264:9558](-)\n",
      "qualifiers:\n",
      "    Key: locus_tag, Value: ['YP_RS23690']\n",
      "\n",
      "type: CDS\n",
      "location: [9264:9558](-)\n",
      "qualifiers:\n",
      "    Key: codon_start, Value: ['1']\n",
      "    Key: db_xref, Value: ['GI:503338034']\n",
      "    Key: inference, Value: ['COORDINATES: ab initio prediction:GeneMarkS+']\n",
      "    Key: locus_tag, Value: ['YP_RS23690']\n",
      "    Key: note, Value: ['Derived by automated computational analysis using gene prediction method: GeneMarkS+.']\n",
      "    Key: product, Value: ['hypothetical protein']\n",
      "    Key: protein_id, Value: ['WP_013572695.1']\n",
      "    Key: transl_table, Value: ['11']\n",
      "    Key: translation, Value: ['MSHFYFLVLHAVKQGFIFGDAKVRTRHINHRFSPYHTPSQNTSCEGLFYVVAKKKPLCLKEGLYVSKKIPGVIHLRKFQTNHKMEHRQCPEVKECSD']\n",
      "\n"
     ]
    }
   ],
   "source": [
    "for feature in record.features:\n",
    "    print(feature)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "No exemplo abaixo verificamos se um SNP em uma determinada posição está dentro de um intervalo do CDS de um gene."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 53,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "<class 'Bio.SeqFeature.SeqFeature'> CDS ['YP_RS22225']\n"
     ]
    }
   ],
   "source": [
    "my_snp = 4360\n",
    "for feature in record.features:\n",
    "    if my_snp in feature and feature.type == 'CDS':\n",
    "        print(f\"{type(feature)} {feature.type} {feature.qualifiers.get('locus_tag')}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "A classe SeqFeature tem o atributo location. Este atributo recece um objeto da classe FeatureLocation, sendo que ele basicamente marca o inicio e fim de uma feature, como no exemplo abaixo."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 39,
   "metadata": {},
   "outputs": [],
   "source": [
    "from Bio.SeqFeature import SeqFeature, FeatureLocation\n",
    "example_feature = SeqFeature(FeatureLocation(1, 18), type=\"gene\", strand=-1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Exercício 2\n",
    "Dado o arquivo gb do exemplo anterior imprima o nome dos genes e quantos genes estão anotados. O nome do gene está representado no campo `feature.qualifiers.get('locus_tag')`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from Bio import SeqIO\n",
    "record = SeqIO.read(\"sequence.gb\", \"genbank\")\n",
    "contador_genes = 0\n",
    "for feature in record.features:\n",
    "    if feature.type == \"gene\":"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 2.4 Extraindo a sequência de um gene\n",
    "Agora imagine que você tem uma sequência de um genoma e queira pegar somente a sequência da feature `example_feature` do exemplo anterior . Uma das formas é pegar o início e fim da feature e usar o slice como já vimos. No final gera o reverso complementar, pois o gene está na fita oposta (strand=-1)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 40,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "AGCCTTTGCCGTCTCGG\n"
     ]
    }
   ],
   "source": [
    "example_parent = Seq(\"ACCGAGACGGCAAAGGCTAGCATAGGTATGAGACTTCCTTCCTGCCAGTGCTGAGGAACTGGGAGCCTAC\")\n",
    "feature_seq = example_parent[example_feature.location.start:example_feature.location.end].reverse_complement()\n",
    "print(feature_seq)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Mas a classe `SeqFeature` já implementa um método chamado `extract()` onde é possível passar uma sequência e ele retorna somente a sequência do intervalo da feature."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 41,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "AGCCTTTGCCGTCTCGG\n"
     ]
    }
   ],
   "source": [
    "feature_seq = example_feature.extract(example_parent)\n",
    "print(feature_seq)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Utilizar o método extract() começa a fazer mais sentido quando nossa feature tem mais de um intervalo. Vamos ao exemplo abaixo onde temos duas FeatureLocation que podem ser os intervalos de dois exons de um gene. Podemos combina-los utilizando o operador '+' e dessa forma criamos um novo objeto que chamamos de join. Esse objeto guarda esses dois FeatureLocation em um objeto da classe CompoundLocation. Sem seguida podemos atribuir essas localizações ao SeqFeture."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 42,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "CompoundLocation([FeatureLocation(ExactPosition(0), ExactPosition(10)), FeatureLocation(ExactPosition(20), ExactPosition(30))], 'join')"
      ]
     },
     "execution_count": 42,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "f1 = FeatureLocation(0,10)\n",
    "f2 = FeatureLocation(20,30)\n",
    "join = f1 + f2\n",
    "join"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 43,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "SeqFeature(CompoundLocation([FeatureLocation(ExactPosition(0), ExactPosition(10), strand=1), FeatureLocation(ExactPosition(20), ExactPosition(30), strand=1)], 'join'), type='gene', location_operator='join')"
      ]
     },
     "execution_count": 43,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "example_feature = SeqFeature(join, type=\"gene\", strand=1)\n",
    "example_feature"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Em seguida podemos utilizar o método `extract()` para pegar as sequências dos dois exons."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 44,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "ACCGAGACGGCATAGGTATG\n"
     ]
    }
   ],
   "source": [
    "feature_seq = example_feature.extract(example_parent)\n",
    "print(feature_seq)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 45,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "20\n",
      "20\n"
     ]
    }
   ],
   "source": [
    "print(len(join))\n",
    "print(len(feature_seq))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Exercício 3\n",
    "Dado o arquivo gb do exemplo anterior e para cada feature do tipo gene extraia e imprima a sequência.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from Bio import SeqIO\n",
    "record = SeqIO.read(\"sequence.gb\", \"genbank\")\n",
    "for feature in record.features:\n",
    "    if feature.type == \"gene\":\n",
    "        "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Exercício 4\n",
    "Abra o arquivo \"sequence.gb\" do exemplo anterior e para cada feature 'gene' extraia a sequência de DNA. Em seguida calcule o conteúdo GC dessa sequência e se ele for maior ou igual a 45 então traduza ela utilizando a tabela de códons de número 11 até o primeiro STOP códon. No final retorne a sequência da proteína no formato FASTA, sendo que o identificar será o atríbuto 'locus_tag' dos qualifiers da classe Seqfeature e a descrição terá o conteúdo GC.\n",
    "Ex:\n",
    "\n",
    "\n",
    "\\>YP_RS22210 GC=52.78592375366569\n",
    "MVTFETVMEIKILHKQGMSSRAIARELGISRNTVKRYLQAKSEPPKYTPRPAVASLLDEYRDYIRQRIADAHPYKIPATVIAREIRDQGYRGGMTILRAFIRSLSVPQEQEPAVRFETEPGRQMQVDWGTMRNGRSPLHVFVAVLGYSRMLYIEFTDNMRYDTLETCHRNAFRFFGGVPREVLYDNMKTVVLQRDAYQTGQHRFHPSLWQFGKEMGFSPRLCRPFRAQTKGKVERMVQYTRNSFYIPLMTRLRPMGITVDVETANRHGLRWLHDVANQRKHETIQARPCDRWLEEQQSMLALPPEKKEYDVHLDENLVNFDKHPLHHPLSIYDSFCRGVA'"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 49,
   "metadata": {},
   "outputs": [],
   "source": [
    "from Bio import SeqIO\n",
    "from Bio.SeqUtils import GC\n",
    "record = SeqIO.read(\"sequence.gb\", \"genbank\")"
   ]
  }
 ],
 "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.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}